Skip to content

Commit

Permalink
add private members, initialize them if needed
Browse files Browse the repository at this point in the history
  • Loading branch information
juliannz committed Dec 19, 2023
1 parent ce20634 commit 459ea13
Show file tree
Hide file tree
Showing 3 changed files with 148 additions and 75 deletions.
152 changes: 104 additions & 48 deletions mt-kahypar/partition/refinement/spectral/solvers/slepc_gevp.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -29,21 +29,43 @@
namespace mt_kahypar {
namespace spectral {

SLEPcGEVPSolver::SLEPcGEVPSolver() {
void SLEPcGEVPSolver::start_slepc() {
PetscFunctionBeginUser;

/* TODO parse args */
/* TODO parse cl args */

PetscBool slepc_running_but_with_pets_datatype;
SlepcInitialized(&slepc_running_but_with_pets_datatype);

if (slepc_running_but_with_pets_datatype == PETSC_FALSE) {
CallPetsc(SlepcInitialize(0, nullptr, (char*)0, nullptr)); /* TODO pass cl args? */
} else {
/* TODO
if(EPSIsCreated(&eps)) {
CallPetsc(EPSDestroy(&eps));
}*/
}

CallPetsc(SlepcInitialize(0, nullptr, (char*)0, nullptr)); /* TODO cl args */
CallPetsc(EPSCreate(GLOBAL_COMMUNICATOR, &eps));

SLEPcGEVPSolver::getN = GET_N_NAIVE;
SLEPcGEVPSolver::getProblemType = GET_PBT_NAIVE;
getN = GET_N_NAIVE;
getProblemType = GET_PBT_NAIVE;
getEpsType = GET_TYPE_NAIVE;

slepc_running = true;

PetscFunctionReturnVoid();
}

void SLEPcGEVPSolver::initialize(Operator& a, Operator& b) {
void SLEPcGEVPSolver::setProblem(Operator& a, Operator& b) {
PetscFunctionBeginUser;

// initialize

if (!slepc_running) {
start_slepc();
}

// validation

/*
Expand All @@ -55,32 +77,37 @@ void SLEPcGEVPSolver::initialize(Operator& a, Operator& b) {

op_a = &a;
op_b = &b;

solved = false;

size_t n = op_a->dimension();

// create matrices

CallPetsc(MatCreateShell(GLOBAL_COMMUNICATOR, SLEPcGEVPSolver::getN(n), SLEPcGEVPSolver::getN(n), SLEPcGEVPSolver::getN(n), SLEPcGEVPSolver::getN(n), &n, &mat_A));
CallPetsc(MatShellSetContext(mat_A, (void *) op_a));
reset_matrices();

// A
CallPetsc(MatCreateShell(GLOBAL_COMMUNICATOR, getN(n), getN(n), getN(n), getN(n), (void *) op_a, &mat_A));
CallPetsc(MatShellSetOperation(mat_A, MATOP_MULT, (void(*)(void)) &matMult));
CallPetsc(MatShellSetOperation(mat_A, MATOP_MULT_TRANSPOSE,(void(*)(void)) &matMultTranspose));
// CallPetsc(MatShellSetOperation(mat_A,MATOP_GET_DIAGONAL,(void(*)(void))MatGetDiagonal_Laplacian2D));
CallPetsc(MatShellSetOperation(mat_A,MATOP_GET_DIAGONAL,(void(*)(void)) &getDiagonal));

CallPetsc(MatCreateShell(GLOBAL_COMMUNICATOR, SLEPcGEVPSolver::getN(n), SLEPcGEVPSolver::getN(n), SLEPcGEVPSolver::getN(n), SLEPcGEVPSolver::getN(n), &n, &mat_B));
CallPetsc(MatShellSetContext(mat_B, (void *) op_b));
// B
CallPetsc(MatCreateShell(GLOBAL_COMMUNICATOR, getN(n), getN(n), getN(n), getN(n), (void *) op_b, &mat_B));
CallPetsc(MatShellSetOperation(mat_B, MATOP_MULT, (void(*)(void)) &matMult));
CallPetsc(MatShellSetOperation(mat_B, MATOP_MULT_TRANSPOSE,(void(*)(void)) &matMultTranspose));
// CallPetsc(MatShellSetOperation(mat_B,MATOP_GET_DIAGONAL,(void(*)(void))MatGetDiagonal_Laplacian2D));
CallPetsc(MatShellSetOperation(mat_B,MATOP_GET_DIAGONAL,(void(*)(void)) &getDiagonal));

// set options

CallPetsc(EPSSetType(eps, SLEPcGEVPSolver::getEpsType(*op_a, *op_b)));
CallPetsc(EPSSetProblemType(eps, SLEPcGEVPSolver::getProblemType(*op_a, *op_b)));
CallPetsc(EPSSetType(eps, getEpsType(*op_a, *op_b)));
CallPetsc(EPSSetProblemType(eps, getProblemType(*op_a, *op_b)));

CallPetsc(EPSSetOperators(eps, mat_A, mat_B));

CallPetsc(EPSSetFromOptions(eps)); /* TODO really needed? */
CallPetsc(EPSSetFromOptions(eps)); /* TODO really needed? */

PetscFunctionReturnVoid();
}

bool SLEPcGEVPSolver::nextEigenpair(Skalar& eval, Vector& evec) {
Expand All @@ -91,82 +118,111 @@ bool SLEPcGEVPSolver::nextEigenpair(Skalar& eval, Vector& evec) {
}

// ...
return true;
PetscFunctionReturn(true);
}

SLEPcGEVPSolver::~SLEPcGEVPSolver() {
void SLEPcGEVPSolver::reset_matrices() {
PetscFunctionBeginUser;

if (mat_A != nullptr) {
MatDestroy(&mat_A);
}
if (mat_B != nullptr) {
MatDestroy(&mat_B);
}

PetscFunctionReturnVoid();
}

void SLEPcGEVPSolver::end_slepc() {
PetscFunctionBeginUser;

reset_matrices();

CallPetsc(EPSDestroy(&eps));
CallPetsc(SlepcFinalize());

//CallPetsc(SlepcFinalize());

slepc_running = false;

PetscFunctionReturnVoid();
}

SLEPcGEVPSolver::~SLEPcGEVPSolver() {
if (slepc_running) {
end_slepc();
}
}

void SLEPcGEVPSolver::solve() {
PetscFunctionBeginUser;

// solve

/* CallPetsc(EPSSolve(eps)); */
// CallPetsc(EPSSolve(eps));
solved = true;

PetscFunctionReturnVoid();
}

void SLEPcGEVPSolver::vector2Vec(Vector& vector, Vec& vec) {
void SLEPcGEVPSolver::vector2Vec(Vector& vector, Vec vec) {
/* TODO */
}

void SLEPcGEVPSolver::vec2vector(Vec& vec, Vector& vector) {
void SLEPcGEVPSolver::vec2vector(Vec vec, Vector& vector) {
/* TODO */
}

int (*SLEPcGEVPSolver::matMult) (Mat&, Vec&, Vec&) = [](Mat& M, Vec& x, Vec& y) {
PetscErrorCode SLEPcGEVPSolver::matMult(Mat M, Vec x, Vec y) {
PetscFunctionBeginUser;

Operator *op_m;
getMatContext(M, op_m);
getMatContext(M, &op_m);

Vector sp_x(op_m->dimension());
Vector sp_y(op_m->dimension());
SLEPcGEVPSolver::vec2vector(x, sp_x);
vec2vector(x, sp_x);
op_m->apply(sp_x, sp_y);
SLEPcGEVPSolver::vector2Vec(sp_y, y);
return 0; /* TODO error code */
};
vector2Vec(sp_y, y);

PetscFunctionReturn(0); /* TODO error code */
}

int (*SLEPcGEVPSolver::matMultTranspose) (Mat&, Vec&, Vec&) = [](Mat& M, Vec& x, Vec& y) {
PetscErrorCode SLEPcGEVPSolver::matMultTranspose(Mat M, Vec x, Vec y) {
PetscFunctionBeginUser;

Operator *op_m;
getMatContext(M, op_m);
getMatContext(M, &op_m);

if (op_m->isSymmetric()) {
return SLEPcGEVPSolver::matMult(M, x, y);
PetscFunctionReturn(matMult(M, x, y));
} else {
return -1; /* TODO */
PetscFunctionReturn(-1); /* TODO */
}
};
}

void SLEPcGEVPSolver::getMatContext(Mat& M, Operator *op) {
PetscErrorCode SLEPcGEVPSolver::getDiagonal(Mat M, Vec diag) {
PetscFunctionBeginUser;

Operator *op_m;
getMatContext(M, &op_m);

void *ctx;
CallPetsc(MatShellGetContext(M, &ctx));
*op = *((Operator *)ctx);
}

Vector d(op_m->dimension());
op_m->getDiagonal(d);
vector2Vec(d, diag);

size_t (*SLEPcGEVPSolver::GET_N_NAIVE) (size_t) = [](size_t n) { return n*n; };
PetscFunctionReturn(0); /* TODO error code */
}

EPSProblemType (*SLEPcGEVPSolver::GET_PBT_NAIVE) (Operator&, Operator&) = [](Operator& a, Operator& b) {
return a.isSymmetric() && b.isSymmetric() ? EPS_GHEP : EPS_GNHEP;
};
void SLEPcGEVPSolver::getMatContext(Mat M, Operator **op) {
PetscFunctionBeginUser;

EPSType (*SLEPcGEVPSolver::GET_TYPE_NAIVE) (Operator&, Operator&) = [](Operator& a, Operator& b) {
return EPSLOBPCG;
};
void *ctx;
CallPetsc(MatShellGetContext(M, &ctx));
*op = (Operator *) ctx;

size_t (*SLEPcGEVPSolver::getN) (size_t) = GET_N_NAIVE;
EPSProblemType (*SLEPcGEVPSolver::getProblemType) (Operator&, Operator&) = GET_PBT_NAIVE;
EPSType (*SLEPcGEVPSolver::getEpsType) (Operator&, Operator&) = GET_TYPE_NAIVE;
PetscFunctionReturnVoid();
}

}
}
69 changes: 43 additions & 26 deletions mt-kahypar/partition/refinement/spectral/solvers/slepc_gevp.h
Original file line number Diff line number Diff line change
Expand Up @@ -43,47 +43,64 @@ class SLEPcGEVPSolver : public GEVPSolver {
private:
// attributes

Operator *op_a;
Operator *op_b;
bool slepc_running = false;
bool solved = false;

// custom matrices
Operator *op_a = nullptr;
Operator *op_b = nullptr;

Mat mat_A;
Mat mat_B;
// petsc matrices
Mat mat_A = nullptr;
Mat mat_B = nullptr;

// petsc eigen problem solver
EPS eps;
bool solved = false;


// methods

void start_slepc();
void solve();
void reset_matrices();
void end_slepc();

// conversion
/* TODO maybe not needed? */
static void vector2Vec(Vector& vector, Vec vec);
static void vec2vector(Vec vec, Vector& vector);

static void getMatContext(Mat M, Operator **op);

// algebraic operations
static PetscErrorCode matMult(Mat M, Vec x, Vec y);
static PetscErrorCode matMultTranspose(Mat M, Vec x, Vec y);
static PetscErrorCode getDiagonal(Mat M, Vec diag);

// lambdas
size_t (*getN) (size_t) = [](size_t n) { throw "unassigned"; return n; };
EPSProblemType (*getProblemType) (Operator&, Operator&) = [](Operator& a, Operator& b) { throw "unassigned"; return EPS_HEP; };
EPSType (*getEpsType) (Operator&, Operator&) = [](Operator& a, Operator& b) { throw "unassigned"; return EPSLOBPCG; };


// constants

static constexpr MPI_Comm& GLOBAL_COMMUNICATOR = PETSC_COMM_WORLD;

// matrix storage size calculation strategys
static size_t (*GET_N_NAIVE) (size_t);
static constexpr size_t (*GET_N_NAIVE) (size_t) = [](size_t n) { return n*n; };
/* ... */

// problem type strategys
static EPSProblemType (*GET_PBT_NAIVE) (Operator&, Operator&);
static constexpr EPSProblemType (*GET_PBT_NAIVE) (Operator&, Operator&) = [](Operator& a, Operator& b) {
return a.isSymmetric() && b.isSymmetric() ? EPS_GHEP : EPS_GNHEP;
};
/* ... */

// eps type strategys
static EPSType (*GET_TYPE_NAIVE) (Operator&, Operator&);

// methods

void solve();

static void vector2Vec(Vector& vector, Vec& vec);
static void vec2vector(Vec& vec, Vector& vector);

static void getMatContext(Mat& M, Operator *op);


// lambdas
static size_t (*getN) (size_t);
static EPSProblemType (*getProblemType) (Operator&, Operator&);
static EPSType (*getEpsType) (Operator&, Operator&);

static int (*matMult) (Mat& M, Vec& x, Vec& y);
static int (*matMultTranspose) (Mat& M, Vec& x, Vec& y);
static constexpr EPSType (*GET_TYPE_NAIVE) (Operator&, Operator&) = [](Operator& a, Operator& b) {
return EPSLOBPCG;
};
};

}
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -128,7 +128,7 @@ namespace mt_kahypar {
generateHintGraphLaplacian(hintSolution, hintGraphLaplacian);

spectral::SLEPcGEVPSolver solver; /* TODO get gevp variant otherwise */
solver.initialize(graphLaplacian, /* baseBalance + */ hintGraphLaplacian);
solver.setProblem(graphLaplacian, /* baseBalance + */ hintGraphLaplacian);

/* TODO */
// to see something:
Expand Down

0 comments on commit 459ea13

Please sign in to comment.