Skip to content

Code Improvements #228

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Closed
julianlitz wants to merge 2 commits into from
Closed

Code Improvements #228

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
77777cb
fix
julianlitz Apr 22, 2026
6da8785
fix
julianlitz Apr 22, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
12 changes: 8 additions & 4 deletions include/DirectSolver/DirectSolver-COO-MUMPS-Give/buildSolverMatrix.inl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -42,6 +42,7 @@ void DirectSolver_COO_MUMPS_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int
const double coeff2 = 0.5 * (k1 + k2) / h2;
const double coeff3 = 0.5 * (h1 + h2) / k1;
const double coeff4 = 0.5 * (h1 + h2) / k2;
const double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

const int center_nz_index = getSolverMatrixIndex(i_r, i_theta);
const int left_nz_index = getSolverMatrixIndex(i_r - 1, i_theta);
Expand All @@ -63,7 +64,7 @@ void DirectSolver_COO_MUMPS_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int

offset = CenterStencil[StencilPosition::Center];
col = center_index;
val = 0.25 * (h1 + h2) * (k1 + k2) * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
val = coeff5 * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
updateMatrixElement(solver_matrix, ptr, offset, row, col, val);

offset = CenterStencil[StencilPosition::Left];
Expand Down Expand Up @@ -269,6 +270,7 @@ void DirectSolver_COO_MUMPS_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int
double coeff2 = 0.5 * (k1 + k2) / h2;
double coeff3 = 0.5 * (h1 + h2) / k1;
double coeff4 = 0.5 * (h1 + h2) / k2;
double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

const int center_nz_index = getSolverMatrixIndex(i_r, i_theta);
const int left_nz_index = getSolverMatrixIndex(i_r, i_theta_AcrossOrigin);
Expand All @@ -289,7 +291,7 @@ void DirectSolver_COO_MUMPS_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int

offset = CenterStencil[StencilPosition::Center];
col = center_index;
val = 0.25 * (h1 + h2) * (k1 + k2) * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
val = coeff5 * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
updateMatrixElement(solver_matrix, ptr, offset, row, col, val);

offset = CenterStencil[StencilPosition::Left];
Expand Down Expand Up @@ -424,6 +426,7 @@ void DirectSolver_COO_MUMPS_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int
const double coeff2 = 0.5 * (k1 + k2) / h2;
const double coeff3 = 0.5 * (h1 + h2) / k1;
const double coeff4 = 0.5 * (h1 + h2) / k2;
const double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

const int i_theta_M1 = grid.wrapThetaIndex(i_theta - 1);
const int i_theta_P1 = grid.wrapThetaIndex(i_theta + 1);
Expand All @@ -448,7 +451,7 @@ void DirectSolver_COO_MUMPS_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int

offset = CenterStencil[StencilPosition::Center];
col = center_index;
val = 0.25 * (h1 + h2) * (k1 + k2) * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
val = coeff5 * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
updateMatrixElement(solver_matrix, ptr, offset, row, col, val);

/* REMOVED: Moved to the right hand side to make the matrix symmetric */
Expand Down Expand Up @@ -604,6 +607,7 @@ void DirectSolver_COO_MUMPS_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int
const double coeff2 = 0.5 * (k1 + k2) / h2;
const double coeff3 = 0.5 * (h1 + h2) / k1;
const double coeff4 = 0.5 * (h1 + h2) / k2;
const double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

const int i_theta_M1 = grid.wrapThetaIndex(i_theta - 1);
const int i_theta_P1 = grid.wrapThetaIndex(i_theta + 1);
Expand All @@ -628,7 +632,7 @@ void DirectSolver_COO_MUMPS_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int

offset = CenterStencil[StencilPosition::Center];
col = center_index;
val = 0.25 * (h1 + h2) * (k1 + k2) * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
val = coeff5 * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
updateMatrixElement(solver_matrix, ptr, offset, row, col, val);

offset = CenterStencil[StencilPosition::Left];
Expand Down
56 changes: 28 additions & 28 deletions include/DirectSolver/DirectSolver-COO-MUMPS-Take/buildSolverMatrix.inl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -41,6 +41,7 @@ void DirectSolver_COO_MUMPS_Take<LevelCacheType>::nodeBuildSolverMatrixTake(
double coeff2 = 0.5 * (k1 + k2) / h2;
double coeff3 = 0.5 * (h1 + h2) / k1;
double coeff4 = 0.5 * (h1 + h2) / k2;
double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

int center_nz_index = getSolverMatrixIndex(i_r, i_theta);

Expand All @@ -59,13 +60,12 @@ void DirectSolver_COO_MUMPS_Take<LevelCacheType>::nodeBuildSolverMatrixTake(
double bottom_value = -coeff3 * (att(center_index) + att(bottom_index)); /* Bottom */
double top_value = -coeff4 * (att(center_index) + att(top_index)); /* Top */

double center_value =
(+0.25 * (h1 + h2) * (k1 + k2) * coeff_beta[center_index] * fabs(detDF(center_index)) /* beta_{i,j} */
- left_value /* Center: (Left) */
- right_value /* Center: (Right) */
- bottom_value /* Center: (Bottom) */
- top_value /* Center: (Top) */
);
double center_value = (coeff5 * coeff_beta[center_index] * fabs(detDF(center_index)) /* beta_{i,j} */
- left_value /* Center: (Left) */
- right_value /* Center: (Right) */
- bottom_value /* Center: (Bottom) */
- top_value /* Center: (Top) */
);

double bottom_left_value = -0.25 * (art(left_index) + art(bottom_index)); /* Bottom Left */
double bottom_right_value = +0.25 * (art(right_index) + art(bottom_index)); /* Bottom Right */
Expand Down Expand Up @@ -168,6 +168,7 @@ void DirectSolver_COO_MUMPS_Take<LevelCacheType>::nodeBuildSolverMatrixTake(
double coeff2 = 0.5 * (k1 + k2) / h2;
double coeff3 = 0.5 * (h1 + h2) / k1;
double coeff4 = 0.5 * (h1 + h2) / k2;
double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

int center_nz_index = getSolverMatrixIndex(i_r, i_theta);

Expand All @@ -184,13 +185,12 @@ void DirectSolver_COO_MUMPS_Take<LevelCacheType>::nodeBuildSolverMatrixTake(
double bottom_value = -coeff3 * (att(center_index) + att(bottom_index)); /* Bottom */
double top_value = -coeff4 * (att(center_index) + att(top_index)); /* Top */

double center_value =
(+0.25 * (h1 + h2) * (k1 + k2) * coeff_beta[center_index] * fabs(detDF(center_index)) /* beta_{i,j} */
- left_value /* Center: (Left) */
- right_value /* Center: (Right) */
- bottom_value /* Center: (Bottom) */
- top_value /* Center: (Top) */
);
double center_value = (coeff5 * coeff_beta[center_index] * fabs(detDF(center_index)) /* beta_{i,j} */
- left_value /* Center: (Left) */
- right_value /* Center: (Right) */
- bottom_value /* Center: (Bottom) */
- top_value /* Center: (Top) */
);

double bottom_right_value = +0.25 * (art(right_index) + art(bottom_index)); /* Bottom Right */
double top_right_value = -0.25 * (art(right_index) + art(top_index)); /* Top Right */
Expand Down Expand Up @@ -257,6 +257,7 @@ void DirectSolver_COO_MUMPS_Take<LevelCacheType>::nodeBuildSolverMatrixTake(
double coeff2 = 0.5 * (k1 + k2) / h2;
double coeff3 = 0.5 * (h1 + h2) / k1;
double coeff4 = 0.5 * (h1 + h2) / k2;
double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

int center_nz_index = getSolverMatrixIndex(i_r, i_theta);

Expand All @@ -275,13 +276,12 @@ void DirectSolver_COO_MUMPS_Take<LevelCacheType>::nodeBuildSolverMatrixTake(
double bottom_value = -coeff3 * (att(center_index) + att(bottom_index)); /* Bottom */
double top_value = -coeff4 * (att(center_index) + att(top_index)); /* Top */

double center_value =
(+0.25 * (h1 + h2) * (k1 + k2) * coeff_beta[center_index] * fabs(detDF(center_index)) /* beta_{i,j} */
- left_value /* Center: (Left) */
- right_value /* Center: (Right) */
- bottom_value /* Center: (Bottom) */
- top_value /* Center: (Top) */
);
double center_value = (coeff5 * coeff_beta[center_index] * fabs(detDF(center_index)) /* beta_{i,j} */
- left_value /* Center: (Left) */
- right_value /* Center: (Right) */
- bottom_value /* Center: (Bottom) */
- top_value /* Center: (Top) */
);

double bottom_left_value = -0.25 * (art(left_index) + art(bottom_index)); /* Bottom Left */
double bottom_right_value = +0.25 * (art(right_index) + art(bottom_index)); /* Bottom Right */
Expand Down Expand Up @@ -364,6 +364,7 @@ void DirectSolver_COO_MUMPS_Take<LevelCacheType>::nodeBuildSolverMatrixTake(
double coeff2 = 0.5 * (k1 + k2) / h2;
double coeff3 = 0.5 * (h1 + h2) / k1;
double coeff4 = 0.5 * (h1 + h2) / k2;
double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

int center_nz_index = getSolverMatrixIndex(i_r, i_theta);

Expand All @@ -382,13 +383,12 @@ void DirectSolver_COO_MUMPS_Take<LevelCacheType>::nodeBuildSolverMatrixTake(
double bottom_value = -coeff3 * (att(center_index) + att(bottom_index)); /* Bottom */
double top_value = -coeff4 * (att(center_index) + att(top_index)); /* Top */

double center_value =
(+0.25 * (h1 + h2) * (k1 + k2) * coeff_beta[center_index] * fabs(detDF(center_index)) /* beta_{i,j} */
- left_value /* Center: (Left) */
- right_value /* Center: (Right) */
- bottom_value /* Center: (Bottom) */
- top_value /* Center: (Top) */
);
double center_value = (coeff5 * coeff_beta[center_index] * fabs(detDF(center_index)) /* beta_{i,j} */
- left_value /* Center: (Left) */
- right_value /* Center: (Right) */
- bottom_value /* Center: (Bottom) */
- top_value /* Center: (Top) */
);

double bottom_left_value = -0.25 * (art(left_index) + art(bottom_index)); /* Bottom Left */
double bottom_right_value = +0.25 * (art(right_index) + art(bottom_index)); /* Bottom Right */
Expand Down
21 changes: 13 additions & 8 deletions include/DirectSolver/DirectSolver-CSR-LU-Give/buildSolverMatrix.inl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -38,6 +38,7 @@ void DirectSolver_CSR_LU_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int i_r
const double coeff2 = 0.5 * (k1 + k2) / h2;
const double coeff3 = 0.5 * (h1 + h2) / k1;
const double coeff4 = 0.5 * (h1 + h2) / k2;
const double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

const int center_index = grid.index(i_r, i_theta);
const int left_index = grid.index(i_r - 1, i_theta);
Expand All @@ -52,7 +53,7 @@ void DirectSolver_CSR_LU_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int i_r

offset = CenterStencil[StencilPosition::Center];
col = center_index;
val = 0.25 * (h1 + h2) * (k1 + k2) * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
val = coeff5 * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
updateMatrixElement(solver_matrix, offset, row, col, val);

offset = CenterStencil[StencilPosition::Left];
Expand Down Expand Up @@ -250,14 +251,16 @@ void DirectSolver_CSR_LU_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int i_r
assert(grid.ntheta() % 2 == 0);
const int i_theta_AcrossOrigin = grid.wrapThetaIndex(i_theta + grid.ntheta() / 2);

double h1 = 2.0 * grid.radius(0);
double h2 = grid.radialSpacing(i_r);
double k1 = grid.angularSpacing(i_theta_M1);
double k2 = grid.angularSpacing(i_theta);
double h1 = 2.0 * grid.radius(0);
double h2 = grid.radialSpacing(i_r);
double k1 = grid.angularSpacing(i_theta_M1);
double k2 = grid.angularSpacing(i_theta);

double coeff1 = 0.5 * (k1 + k2) / h1;
double coeff2 = 0.5 * (k1 + k2) / h2;
double coeff3 = 0.5 * (h1 + h2) / k1;
double coeff4 = 0.5 * (h1 + h2) / k2;
double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

const int center_index = grid.index(i_r, i_theta);
const int left_index = grid.index(i_r, i_theta_AcrossOrigin);
Expand All @@ -272,7 +275,7 @@ void DirectSolver_CSR_LU_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int i_r

offset = CenterStencil[StencilPosition::Center];
col = center_index;
val = 0.25 * (h1 + h2) * (k1 + k2) * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
val = coeff5 * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
updateMatrixElement(solver_matrix, offset, row, col, val);

offset = CenterStencil[StencilPosition::Left];
Expand Down Expand Up @@ -404,6 +407,7 @@ void DirectSolver_CSR_LU_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int i_r
const double coeff2 = 0.5 * (k1 + k2) / h2;
const double coeff3 = 0.5 * (h1 + h2) / k1;
const double coeff4 = 0.5 * (h1 + h2) / k2;
const double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

const int i_theta_M1 = grid.wrapThetaIndex(i_theta - 1);
const int i_theta_P1 = grid.wrapThetaIndex(i_theta + 1);
Expand All @@ -421,7 +425,7 @@ void DirectSolver_CSR_LU_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int i_r

offset = CenterStencil[StencilPosition::Center];
col = center_index;
val = 0.25 * (h1 + h2) * (k1 + k2) * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
val = coeff5 * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
updateMatrixElement(solver_matrix, offset, row, col, val);

offset = CenterStencil[StencilPosition::Left];
Expand Down Expand Up @@ -564,6 +568,7 @@ void DirectSolver_CSR_LU_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int i_r
const double coeff2 = 0.5 * (k1 + k2) / h2;
const double coeff3 = 0.5 * (h1 + h2) / k1;
const double coeff4 = 0.5 * (h1 + h2) / k2;
const double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

const int i_theta_M1 = grid.wrapThetaIndex(i_theta - 1);
const int i_theta_P1 = grid.wrapThetaIndex(i_theta + 1);
Expand All @@ -581,7 +586,7 @@ void DirectSolver_CSR_LU_Give<LevelCacheType>::nodeBuildSolverMatrixGive(int i_r

offset = CenterStencil[StencilPosition::Center];
col = center_index;
val = 0.25 * (h1 + h2) * (k1 + k2) * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
val = coeff5 * coeff_beta * std::fabs(detDF); /* beta_{i,j} */
updateMatrixElement(solver_matrix, offset, row, col, val);

offset = CenterStencil[StencilPosition::Left];
Expand Down
28 changes: 14 additions & 14 deletions include/DirectSolver/DirectSolver-CSR-LU-Take/buildSolverMatrix.inl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -37,6 +37,7 @@ void DirectSolver_CSR_LU_Take<LevelCacheType>::nodeBuildSolverMatrixTake(
double coeff2 = 0.5 * (k1 + k2) / h2;
double coeff3 = 0.5 * (h1 + h2) / k1;
double coeff4 = 0.5 * (h1 + h2) / k2;
double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

int center_index = grid.index(i_r, i_theta);
int left_index = grid.index(i_r - 1, i_theta);
Expand All @@ -53,13 +54,12 @@ void DirectSolver_CSR_LU_Take<LevelCacheType>::nodeBuildSolverMatrixTake(
double bottom_value = -coeff3 * (att[center_index] + att[bottom_index]); /* Bottom */
double top_value = -coeff4 * (att[center_index] + att[top_index]); /* Top */

double center_value =
(+0.25 * (h1 + h2) * (k1 + k2) * coeff_beta[center_index] * std::fabs(detDF[center_index]) /* beta_{i,j} */
- left_value /* Center: (Left) */
- right_value /* Center: (Right) */
- bottom_value /* Center: (Bottom) */
- top_value /* Center: (Top) */
);
double center_value = (coeff5 * coeff_beta[center_index] * std::fabs(detDF[center_index]) /* beta_{i,j} */
- left_value /* Center: (Left) */
- right_value /* Center: (Right) */
- bottom_value /* Center: (Bottom) */
- top_value /* Center: (Top) */
);

double bottom_left_value = -0.25 * (art[left_index] + art[bottom_index]); /* Bottom Left */
double bottom_right_value = +0.25 * (art[right_index] + art[bottom_index]); /* Bottom Right */
Expand Down Expand Up @@ -158,6 +158,7 @@ void DirectSolver_CSR_LU_Take<LevelCacheType>::nodeBuildSolverMatrixTake(
double coeff2 = 0.5 * (k1 + k2) / h2;
double coeff3 = 0.5 * (h1 + h2) / k1;
double coeff4 = 0.5 * (h1 + h2) / k2;
double coeff5 = 0.25 * (h1 + h2) * (k1 + k2);

int center_index = grid.index(i_r, i_theta);
int left_index = grid.index(i_r, i_theta_AcrossOrigin);
Expand All @@ -172,13 +173,12 @@ void DirectSolver_CSR_LU_Take<LevelCacheType>::nodeBuildSolverMatrixTake(
double bottom_value = -coeff3 * (att[center_index] + att[bottom_index]); /* Bottom */
double top_value = -coeff4 * (att[center_index] + att[top_index]); /* Top */

double center_value =
(+0.25 * (h1 + h2) * (k1 + k2) * coeff_beta[center_index] * fabs(detDF[center_index]) /* beta_{i,j} */
- left_value /* Center: (Left) */
- right_value /* Center: (Right) */
- bottom_value /* Center: (Bottom) */
- top_value /* Center: (Top) */
);
double center_value = (coeff5 * coeff_beta[center_index] * fabs(detDF[center_index]) /* beta_{i,j} */
- left_value /* Center: (Left) */
- right_value /* Center: (Right) */
- bottom_value /* Center: (Bottom) */
- top_value /* Center: (Top) */
);

double bottom_right_value = +0.25 * (art[right_index] + art[bottom_index]); /* Bottom Right */
double top_right_value = -0.25 * (art[right_index] + art[top_index]); /* Top Right */
Expand Down
Loading
Loading