Skip to content

Commit

Permalink
Merge branch 'memory-ownership' into 'develop'
Browse files Browse the repository at this point in the history
Memory owneship

Closes #2

See merge request ecpcitest/exasgd/resolve!40
  • Loading branch information
pelesh committed Oct 7, 2023
2 parents 0bbfbbb + 7ee611e commit 4d6e99d
Show file tree
Hide file tree
Showing 28 changed files with 293 additions and 201 deletions.
8 changes: 4 additions & 4 deletions examples/r_KLU_GLU.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,8 @@
#include <resolve/LinSolverDirectCuSolverGLU.hpp>
#include <resolve/LinAlgWorkspace.hpp>

using namespace ReSolve::constants;

int main(int argc, char *argv[])
{
// Use the same data types as those you specified in ReSolve build.
Expand All @@ -19,6 +21,7 @@ int main(int argc, char *argv[])
using vector_type = ReSolve::vector::Vector;
using matrix_type = ReSolve::matrix::Sparse;

(void) argc; // TODO: Check if the number of input parameters is correct.
std::string matrixFileName = argv[1];
std::string rhsFileName = argv[2];

Expand All @@ -44,9 +47,6 @@ int main(int argc, char *argv[])
vector_type* vec_x;
vector_type* vec_r;

real_type one = 1.0;
real_type minusone = -1.0;

ReSolve::LinSolverDirectKLU* KLU = new ReSolve::LinSolverDirectKLU;
ReSolve::LinSolverDirectCuSolverGLU* GLU = new ReSolve::LinSolverDirectCuSolverGLU(workspace_CUDA);

Expand Down Expand Up @@ -146,7 +146,7 @@ int main(int argc, char *argv[])


matrix_handler->setValuesChanged(true);
matrix_handler->matvec(A, vec_x, vec_r, &one, &minusone,"csr", "cuda");
matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr", "cuda");

printf("\t 2-Norm of the residual: %16.16e\n", sqrt(vector_handler->dot(vec_r, vec_r, "cuda")));

Expand Down
8 changes: 4 additions & 4 deletions examples/r_KLU_GLU_matrix_values_update.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,8 @@

// this updates the matrix values to simulate what CFD/optimization software does.

using namespace ReSolve::constants;

int main(int argc, char *argv[])
{
// Use the same data types as those you specified in ReSolve build.
Expand All @@ -21,6 +23,7 @@ int main(int argc, char *argv[])
using vector_type = ReSolve::vector::Vector;
using matrix_type = ReSolve::matrix::Sparse;

(void) argc; // TODO: Check if the number of input parameters is correct.
std::string matrixFileName = argv[1];
std::string rhsFileName = argv[2];

Expand All @@ -47,9 +50,6 @@ int main(int argc, char *argv[])
vector_type* vec_x;
vector_type* vec_r;

real_type one = 1.0;
real_type minusone = -1.0;

ReSolve::LinSolverDirectKLU* KLU = new ReSolve::LinSolverDirectKLU;
ReSolve::LinSolverDirectCuSolverGLU* GLU = new ReSolve::LinSolverDirectCuSolverGLU(workspace_CUDA);

Expand Down Expand Up @@ -156,7 +156,7 @@ int main(int argc, char *argv[])


matrix_handler->setValuesChanged(true);
matrix_handler->matvec(A, vec_x, vec_r, &one, &minusone,"csr", "cuda");
matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr", "cuda");

printf("\t 2-Norm of the residual: %16.16e\n", sqrt(vector_handler->dot(vec_r, vec_r, "cuda")));

Expand Down
8 changes: 4 additions & 4 deletions examples/r_KLU_KLU.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,7 @@
#include <resolve/LinSolverDirectKLU.hpp>
#include <resolve/LinAlgWorkspace.hpp>

using namespace ReSolve::constants;

int main(int argc, char *argv[])
{
Expand All @@ -19,6 +20,7 @@ int main(int argc, char *argv[])
using real_type = ReSolve::real_type;
using vector_type = ReSolve::vector::Vector;

(void) argc; // TODO: Check if the number of input parameters is correct.
std::string matrixFileName = argv[1];
std::string rhsFileName = argv[2];

Expand All @@ -44,9 +46,6 @@ int main(int argc, char *argv[])
vector_type* vec_x;
vector_type* vec_r;

real_type one = 1.0;
real_type minusone = -1.0;

ReSolve::LinSolverDirectKLU* KLU = new ReSolve::LinSolverDirectKLU;

for (int i = 0; i < numSystems; ++i)
Expand Down Expand Up @@ -134,8 +133,9 @@ int main(int argc, char *argv[])

matrix_handler->setValuesChanged(true);

matrix_handler->matvec(A, vec_x, vec_r, &one, &minusone, "csr", "cuda");
matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE, "csr", "cuda");
real_type* test = vec_r->getData("cpu");
(void) test; // TODO: Do we need `test` variable in this example?

printf("\t 2-Norm of the residual: %16.16e\n", sqrt(vector_handler->dot(vec_r, vec_r, "cuda")));

Expand Down
121 changes: 60 additions & 61 deletions examples/r_KLU_KLU_standalone.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -11,14 +11,15 @@
#include <resolve/LinSolverDirectKLU.hpp>
#include <resolve/LinAlgWorkspace.hpp>

using namespace ReSolve::constants;

int main(int argc, char *argv[])
{
// Use the same data types as those you specified in ReSolve build.
using index_type = ReSolve::index_type;
using real_type = ReSolve::real_type;
using vector_type = ReSolve::vector::Vector;

(void) argc; // TODO: Check if the number of input parameters is correct.
std::string matrixFileName = argv[1];
std::string rhsFileName = argv[2];

Expand All @@ -41,70 +42,68 @@ int main(int argc, char *argv[])
vector_type* vec_x;
vector_type* vec_r;

real_type one = 1.0;
real_type minusone = -1.0;

ReSolve::LinSolverDirectKLU* KLU = new ReSolve::LinSolverDirectKLU;



// Read matrix first
std::cout << "========================================================================================================================"<<std::endl;
std::cout << "Reading: " << matrixFileName << std::endl;
std::cout << "========================================================================================================================"<<std::endl;
std::cout << std::endl;
// Read first matrix
std::ifstream mat_file(matrixFileName);
if(!mat_file.is_open())
{
std::cout << "Failed to open file " << matrixFileName << "\n";
return -1;
}
std::ifstream rhs_file(rhsFileName);
if(!rhs_file.is_open())
{
std::cout << "Failed to open file " << rhsFileName << "\n";
return -1;
}
A_coo = ReSolve::io::readMatrixFromFile(mat_file);
A = new ReSolve::matrix::Csr(A_coo->getNumRows(),
A_coo->getNumColumns(),
A_coo->getNnz(),
A_coo->symmetric(),
A_coo->expanded());

rhs = ReSolve::io::readRhsFromFile(rhs_file);
x = new real_type[A->getNumRows()];
vec_rhs = new vector_type(A->getNumRows());
vec_x = new vector_type(A->getNumRows());
vec_r = new vector_type(A->getNumRows());
std::cout<<"Finished reading the matrix and rhs, size: "<<A->getNumRows()<<" x "<<A->getNumColumns()<< ", nnz: "<< A->getNnz()<< ", symmetric? "<<A->symmetric()<< ", Expanded? "<<A->expanded()<<std::endl;
mat_file.close();
rhs_file.close();

//Now convert to CSR.
matrix_handler->coo2csr(A_coo, A, "cpu");
vec_rhs->update(rhs, "cpu", "cpu");
vec_rhs->setDataUpdated("cpu");
std::cout << "COO to CSR completed. Expanded NNZ: " << A->getNnzExpanded() << std::endl;
//Now call direct solver
KLU->setupParameters(1, 0.1, false);
int status;
KLU->setup(A);
status = KLU->analyze();
std::cout<<"KLU analysis status: "<<status<<std::endl;
status = KLU->factorize();
std::cout << "KLU factorization status: " << status << std::endl;
status = KLU->solve(vec_rhs, vec_x);
std::cout << "KLU solve status: " << status << std::endl;
vec_r->update(rhs, "cpu", "cuda");

matrix_handler->setValuesChanged(true);

matrix_handler->matvec(A, vec_x, vec_r, &one, &minusone, "csr", "cuda");
real_type* test = vec_r->getData("cpu");

printf("\t 2-Norm of the residual: %16.16e\n", sqrt(vector_handler->dot(vec_r, vec_r, "cuda")));
// Read matrix first
std::cout << "========================================================================================================================"<<std::endl;
std::cout << "Reading: " << matrixFileName << std::endl;
std::cout << "========================================================================================================================"<<std::endl;
std::cout << std::endl;
// Read first matrix
std::ifstream mat_file(matrixFileName);
if(!mat_file.is_open())
{
std::cout << "Failed to open file " << matrixFileName << "\n";
return -1;
}
std::ifstream rhs_file(rhsFileName);
if(!rhs_file.is_open())
{
std::cout << "Failed to open file " << rhsFileName << "\n";
return -1;
}
A_coo = ReSolve::io::readMatrixFromFile(mat_file);
A = new ReSolve::matrix::Csr(A_coo->getNumRows(),
A_coo->getNumColumns(),
A_coo->getNnz(),
A_coo->symmetric(),
A_coo->expanded());

rhs = ReSolve::io::readRhsFromFile(rhs_file);
x = new real_type[A->getNumRows()];
vec_rhs = new vector_type(A->getNumRows());
vec_x = new vector_type(A->getNumRows());
vec_r = new vector_type(A->getNumRows());
std::cout<<"Finished reading the matrix and rhs, size: "<<A->getNumRows()<<" x "<<A->getNumColumns()<< ", nnz: "<< A->getNnz()<< ", symmetric? "<<A->symmetric()<< ", Expanded? "<<A->expanded()<<std::endl;
mat_file.close();
rhs_file.close();

//Now convert to CSR.
matrix_handler->coo2csr(A_coo, A, "cpu");
vec_rhs->update(rhs, "cpu", "cpu");
vec_rhs->setDataUpdated("cpu");
std::cout << "COO to CSR completed. Expanded NNZ: " << A->getNnzExpanded() << std::endl;
//Now call direct solver
KLU->setupParameters(1, 0.1, false);
int status;
KLU->setup(A);
status = KLU->analyze();
std::cout<<"KLU analysis status: "<<status<<std::endl;
status = KLU->factorize();
std::cout << "KLU factorization status: " << status << std::endl;
status = KLU->solve(vec_rhs, vec_x);
std::cout << "KLU solve status: " << status << std::endl;
vec_r->update(rhs, "cpu", "cuda");

matrix_handler->setValuesChanged(true);

matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE, "csr", "cuda");
real_type* test = vec_r->getData("cpu");
(void) test; // TODO: Do we need `test` variable in this example?

printf("\t 2-Norm of the residual: %16.16e\n", sqrt(vector_handler->dot(vec_r, vec_r, "cuda")));



Expand Down
9 changes: 4 additions & 5 deletions examples/r_KLU_rf.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -12,13 +12,16 @@
#include <resolve/LinSolverDirectCuSolverRf.hpp>
#include <resolve/LinAlgWorkspace.hpp>

using namespace ReSolve::constants;

int main(int argc, char *argv[] )
{
// Use the same data types as those you specified in ReSolve build.
using index_type = ReSolve::index_type;
using real_type = ReSolve::real_type;
using vector_type = ReSolve::vector::Vector;

(void) argc; // TODO: Check if the number of input parameters is correct.
std::string matrixFileName = argv[1];
std::string rhsFileName = argv[2];

Expand All @@ -45,9 +48,6 @@ int main(int argc, char *argv[] )
vector_type* vec_x;
vector_type* vec_r;

real_type one = 1.0;
real_type minusone = -1.0;

ReSolve::LinSolverDirectKLU* KLU = new ReSolve::LinSolverDirectKLU;
ReSolve::LinSolverDirectCuSolverRf* Rf = new ReSolve::LinSolverDirectCuSolverRf;

Expand Down Expand Up @@ -160,8 +160,7 @@ int main(int argc, char *argv[] )

matrix_handler->setValuesChanged(true);

matrix_handler->matvec(A, vec_x, vec_r, &one, &minusone,"csr", "cuda");
real_type* test = vec_r->getData("cpu");
matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr", "cuda");

printf("\t 2-Norm of the residual: %16.16e\n", sqrt(vector_handler->dot(vec_r, vec_r, "cuda")));

Expand Down
14 changes: 7 additions & 7 deletions examples/r_KLU_rf_FGMRES.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -13,13 +13,16 @@
#include <resolve/LinSolverIterativeFGMRES.hpp>
#include <resolve/LinAlgWorkspace.hpp>

using namespace ReSolve::constants;

int main(int argc, char *argv[])
{
// Use the same data types as those you specified in ReSolve build.
using index_type = ReSolve::index_type;
using real_type = ReSolve::real_type;
using vector_type = ReSolve::vector::Vector;

(void) argc; // TODO: Check if the number of input parameters is correct.
std::string matrixFileName = argv[1];
std::string rhsFileName = argv[2];

Expand All @@ -38,16 +41,13 @@ int main(int argc, char *argv[])
workspace_CUDA->initializeHandles();
ReSolve::MatrixHandler* matrix_handler = new ReSolve::MatrixHandler(workspace_CUDA);
ReSolve::VectorHandler* vector_handler = new ReSolve::VectorHandler(workspace_CUDA);
real_type* rhs;
real_type* x;
real_type* rhs = nullptr;
real_type* x = nullptr;

vector_type* vec_rhs;
vector_type* vec_x;
vector_type* vec_r;

real_type one = 1.0;
real_type minusone = -1.0;

ReSolve::GramSchmidt* GS = new ReSolve::GramSchmidt(vector_handler, ReSolve::GramSchmidt::cgs2);
ReSolve::LinSolverDirectKLU* KLU = new ReSolve::LinSolverDirectKLU;
ReSolve::LinSolverDirectCuSolverRf* Rf = new ReSolve::LinSolverDirectCuSolverRf;
Expand Down Expand Up @@ -136,7 +136,7 @@ int main(int argc, char *argv[])
norm_b = vector_handler->dot(vec_r, vec_r, "cuda");
norm_b = sqrt(norm_b);
matrix_handler->setValuesChanged(true);
matrix_handler->matvec(A, vec_x, vec_r, &one, &minusone,"csr", "cuda");
matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr", "cuda");
printf("\t 2-Norm of the residual : %16.16e\n", sqrt(vector_handler->dot(vec_r, vec_r, "cuda"))/norm_b);
if (i == 1) {
ReSolve::matrix::Csc* L_csc = (ReSolve::matrix::Csc*) KLU->getLFactor();
Expand Down Expand Up @@ -169,7 +169,7 @@ int main(int argc, char *argv[])
FGMRES->resetMatrix(A);
FGMRES->setupPreconditioner("CuSolverRf", Rf);

matrix_handler->matvec(A, vec_x, vec_r, &one, &minusone,"csr", "cuda");
matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr", "cuda");

printf("\t 2-Norm of the residual (before IR): %16.16e\n", sqrt(vector_handler->dot(vec_r, vec_r, "cuda"))/norm_b);

Expand Down
14 changes: 7 additions & 7 deletions examples/r_KLU_rf_FGMRES_reuse_factorization.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -13,13 +13,16 @@
#include <resolve/LinSolverIterativeFGMRES.hpp>
#include <resolve/LinAlgWorkspace.hpp>

using namespace ReSolve::constants;

int main(int argc, char *argv[])
{
// Use the same data types as those you specified in ReSolve build.
using index_type = ReSolve::index_type;
using real_type = ReSolve::real_type;
using vector_type = ReSolve::vector::Vector;

(void) argc; // TODO: Check if the number of input parameters is correct.
std::string matrixFileName = argv[1];
std::string rhsFileName = argv[2];

Expand All @@ -39,16 +42,13 @@ int main(int argc, char *argv[])
workspace_CUDA->initializeHandles();
ReSolve::MatrixHandler* matrix_handler = new ReSolve::MatrixHandler(workspace_CUDA);
ReSolve::VectorHandler* vector_handler = new ReSolve::VectorHandler(workspace_CUDA);
real_type* rhs;
real_type* x;
real_type* rhs = nullptr;
real_type* x = nullptr;

vector_type* vec_rhs;
vector_type* vec_x;
vector_type* vec_r;

real_type one = 1.0;
real_type minusone = -1.0;

ReSolve::GramSchmidt* GS = new ReSolve::GramSchmidt(vector_handler, ReSolve::GramSchmidt::cgs2);

ReSolve::LinSolverDirectKLU* KLU = new ReSolve::LinSolverDirectKLU;
Expand Down Expand Up @@ -137,7 +137,7 @@ int main(int argc, char *argv[])
norm_b = vector_handler->dot(vec_r, vec_r, "cuda");
norm_b = sqrt(norm_b);
matrix_handler->setValuesChanged(true);
matrix_handler->matvec(A, vec_x, vec_r, &one, &minusone,"csr", "cuda");
matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr", "cuda");
printf("\t 2-Norm of the residual : %16.16e\n", sqrt(vector_handler->dot(vec_r, vec_r, "cuda"))/norm_b);
if (i == 1) {
ReSolve::matrix::Csc* L_csc = (ReSolve::matrix::Csc*) KLU->getLFactor();
Expand Down Expand Up @@ -181,7 +181,7 @@ int main(int argc, char *argv[])
matrix_handler->setValuesChanged(true);
FGMRES->resetMatrix(A);

matrix_handler->matvec(A, vec_x, vec_r, &one, &minusone,"csr", "cuda");
matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr", "cuda");

printf("\t 2-Norm of the residual (before IR): %16.16e\n", sqrt(vector_handler->dot(vec_r, vec_r, "cuda"))/norm_b);
printf("\t 2-Norm of the RIGHT HAND SIDE: %16.16e\n", norm_b);
Expand Down
Loading

0 comments on commit 4d6e99d

Please sign in to comment.