From 5ad47189ec88a1a9f3fccfd91537338834869714 Mon Sep 17 00:00:00 2001 From: David Bold Date: Wed, 28 Oct 2026 13:54:28 +0100 Subject: [PATCH 1/6] Check more PETSc return values --- src/solver/impls/imex-bdf2/imex-bdf2.cxx | 78 ++++++++++++------------ 1 file changed, 40 insertions(+), 38 deletions(-) diff --git a/src/solver/impls/imex-bdf2/imex-bdf2.cxx b/src/solver/impls/imex-bdf2/imex-bdf2.cxx index 26721da95b..e2f76572bb 100644 --- a/src/solver/impls/imex-bdf2/imex-bdf2.cxx +++ b/src/solver/impls/imex-bdf2/imex-bdf2.cxx @@ -492,15 +492,15 @@ void IMEXBDF2::constructSNES(SNES* snesIn) { } // Pre-allocate - MatMPIAIJSetPreallocation(Jmf, 0, d_nnz, 0, o_nnz); - MatSetUp(Jmf); - MatSetOption(Jmf, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE); - PetscFree(d_nnz); - PetscFree(o_nnz); + BOUT_DO_PETSC(MatMPIAIJSetPreallocation(Jmf, 0, d_nnz, 0, o_nnz)); + BOUT_DO_PETSC(MatSetUp(Jmf)); + BOUT_DO_PETSC(MatSetOption(Jmf, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE)); + BOUT_DO_PETSC(PetscFree(d_nnz)); + BOUT_DO_PETSC(PetscFree(o_nnz)); // Determine which row/columns of the matrix are locally owned int Istart, Iend; - MatGetOwnershipRange(Jmf, &Istart, &Iend); + BOUT_DO_PETSC(MatGetOwnershipRange(Jmf, &Istart, &Iend)); // Convert local into global indices index += Istart; @@ -547,7 +547,7 @@ void IMEXBDF2::constructSNES(SNES* snesIn) { PetscInt col = ind2 + j; // output.write("SETTING 1: {:d}, {:d}\n", row, col); - MatSetValues(Jmf, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(Jmf, 1, &row, 1, &col, &val, INSERT_VALUES)); } } } @@ -567,7 +567,7 @@ void IMEXBDF2::constructSNES(SNES* snesIn) { for (int j = 0; j < n2d; j++) { PetscInt col = ind0 + j; // output.write("SETTING 2: {:d}, {:d}\n", row, col); - MatSetValues(Jmf, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(Jmf, 1, &row, 1, &col, &val, INSERT_VALUES)); } // 5 point star pattern @@ -593,7 +593,7 @@ void IMEXBDF2::constructSNES(SNES* snesIn) { for (int j = 0; j < n3d; j++) { PetscInt col = ind2 + j; // output.write("SETTING 3: {:d}, {:d}\n", row, col); - MatSetValues(Jmf, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(Jmf, 1, &row, 1, &col, &val, INSERT_VALUES)); } } @@ -610,7 +610,7 @@ void IMEXBDF2::constructSNES(SNES* snesIn) { for (int j = 0; j < n3d; j++) { PetscInt col = ind2 + j; // output.write("SETTING 4: {:d}, {:d}\n", row, col); - MatSetValues(Jmf, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(Jmf, 1, &row, 1, &col, &val, INSERT_VALUES)); } int zm = (z - 1 + nz) % nz; @@ -621,7 +621,7 @@ void IMEXBDF2::constructSNES(SNES* snesIn) { for (int j = 0; j < n3d; j++) { PetscInt col = ind2 + j; // output.write("SETTING 5: {:d}, {:d}\n", row, col); - MatSetValues(Jmf, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(Jmf, 1, &row, 1, &col, &val, INSERT_VALUES)); } } } @@ -631,31 +631,32 @@ void IMEXBDF2::constructSNES(SNES* snesIn) { // Finished marking non-zero entries // Assemble Matrix - MatAssemblyBegin(Jmf, MAT_FINAL_ASSEMBLY); - MatAssemblyEnd(Jmf, MAT_FINAL_ASSEMBLY); + BOUT_DO_PETSC(MatAssemblyBegin(Jmf, MAT_FINAL_ASSEMBLY)); + BOUT_DO_PETSC(MatAssemblyEnd(Jmf, MAT_FINAL_ASSEMBLY)); ISColoring iscoloring; MatColoring coloring; // This new in PETSc 3.5 - MatColoringCreate(Jmf, &coloring); - MatColoringSetType(coloring, MATCOLORINGSL); - MatColoringSetFromOptions(coloring); + BOUT_DO_PETSC(MatColoringCreate(Jmf, &coloring)); + BOUT_DO_PETSC(MatColoringSetType(coloring, MATCOLORINGSL)); + BOUT_DO_PETSC(MatColoringSetFromOptions(coloring)); // Calculate index sets - MatColoringApply(coloring, &iscoloring); - MatColoringDestroy(&coloring); + BOUT_DO_PETSC(MatColoringApply(coloring, &iscoloring)); + BOUT_DO_PETSC(MatColoringDestroy(&coloring)); // Create data structure for SNESComputeJacobianDefaultColor - MatFDColoringCreate(Jmf, iscoloring, &fdcoloring); + BOUT_DO_PETSC(MatFDColoringCreate(Jmf, iscoloring, &fdcoloring)); // Set the function to difference - MatFDColoringSetFunction(fdcoloring, - bout::cast_MatFDColoringFn(FormFunctionForColoring), this); - MatFDColoringSetFromOptions(fdcoloring); - MatFDColoringSetUp(Jmf, iscoloring, fdcoloring); - ISColoringDestroy(&iscoloring); + BOUT_DO_PETSC(MatFDColoringSetFunction( + fdcoloring, bout::cast_MatFDColoringFn(FormFunctionForColoring), this)); + BOUT_DO_PETSC(MatFDColoringSetFromOptions(fdcoloring)); + BOUT_DO_PETSC(MatFDColoringSetUp(Jmf, iscoloring, fdcoloring)); + BOUT_DO_PETSC(ISColoringDestroy(&iscoloring)); - SNESSetJacobian(*snesIn, Jmf, Jmf, SNESComputeJacobianDefaultColor, fdcoloring); + BOUT_DO_PETSC(SNESSetJacobian(*snesIn, Jmf, Jmf, SNESComputeJacobianDefaultColor, + fdcoloring)); // Re-use Jacobian - SNESSetLagJacobian(*snesIn, lag_jacobian); + BOUT_DO_PETSC(SNESSetLagJacobian(*snesIn, lag_jacobian)); // MatView(Jmf, PETSC_VIEWER_DRAW_WORLD); // MatView(Jmf, PETSC_VIEWER_STDOUT_WORLD); @@ -663,40 +664,41 @@ void IMEXBDF2::constructSNES(SNES* snesIn) { // Brute force calculation // NOTE: Slow! - MatCreateAIJ( + BOUT_DO_PETSC(MatCreateAIJ( BoutComm::get(), nlocal, nlocal, // Local sizes PETSC_DETERMINE, PETSC_DETERMINE, // Global sizes 3, // Number of nonzero entries in diagonal portion of local submatrix nullptr, 0, // Number of nonzeros per row in off-diagonal portion of local submatrix - nullptr, &Jmf); + nullptr, &Jmf)); - SNESSetJacobian(*snesIn, Jmf, Jmf, SNESComputeJacobianDefault, this); + BOUT_DO_PETSC(SNESSetJacobian(*snesIn, Jmf, Jmf, SNESComputeJacobianDefault, this)); - MatSetOption(Jmf, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE); + BOUT_DO_PETSC(MatSetOption(Jmf, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE)); } } ///////////////////////////////////////////////////// // Set tolerances - SNESSetTolerances(*snesIn, atol, rtol, PETSC_DEFAULT, max_nonlinear_it, PETSC_DEFAULT); + BOUT_DO_PETSC(SNESSetTolerances(*snesIn, atol, rtol, PETSC_DEFAULT, max_nonlinear_it, + PETSC_DEFAULT)); ///////////////////////////////////////////////////// // Preconditioner // Get KSP context from SNES KSP ksp; - SNESGetKSP(*snesIn, &ksp); + BOUT_DO_PETSC(SNESGetKSP(*snesIn, &ksp)); if (kspsetinitialguessnonzero) { - KSPSetInitialGuessNonzero(ksp, PETSC_TRUE); + BOUT_DO_PETSC(KSPSetInitialGuessNonzero(ksp, PETSC_TRUE)); } - KSPSetTolerances(ksp, - PETSC_DEFAULT, // rtol - PETSC_DEFAULT, // abstol - PETSC_DEFAULT, // dtol (divergence tolerance) - maxl); // Maximum number of iterations + BOUT_DO_PETSC(KSPSetTolerances(ksp, + PETSC_DEFAULT, // rtol + PETSC_DEFAULT, // abstol + PETSC_DEFAULT, // dtol (divergence tolerance) + maxl)); // Maximum number of iterations // Get PC context from KSP PC pc; From 725e7f5f5067d8ef33a6aefb076c25cbeb392a77 Mon Sep 17 00:00:00 2001 From: David Bold Date: Tue, 3 Nov 2026 10:49:35 +0100 Subject: [PATCH 2/6] Check more return values --- .../laplacexy/impls/petsc/laplacexy-petsc.cxx | 4 +- .../laplacexz/impls/petsc/laplacexz-petsc.cxx | 109 +++++++++--------- 2 files changed, 57 insertions(+), 56 deletions(-) diff --git a/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx b/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx index 881a0b7fc0..540647d01d 100644 --- a/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx +++ b/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx @@ -190,8 +190,8 @@ LaplaceXYpetsc::LaplaceXYpetsc(Mesh* m, Options* opt, const CELL_LOC loc) MatMPIAIJSetPreallocation(MatA, 0, d_nnz, 0, o_nnz); MatSetUp(MatA); - PetscFree(d_nnz); - PetscFree(o_nnz); + BOUT_DO_PETSC(PetscFree(d_nnz)); + BOUT_DO_PETSC(PetscFree(o_nnz)); // Determine which row/columns of the matrix are locally owned int Istart = 0; diff --git a/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx b/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx index 3e55262b6d..7538c9d3f2 100644 --- a/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx +++ b/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx @@ -178,9 +178,10 @@ LaplaceXZpetsc::LaplaceXZpetsc(Mesh* m, Options* opt, const CELL_LOC loc) data.yindex = y; // Set size of Matrix on each processor to localN x localN - MatCreate(comm, &data.MatA); - MatSetSizes(data.MatA, localN, localN, PETSC_DETERMINE, PETSC_DETERMINE); - MatSetFromOptions(data.MatA); + BOUT_DO_PETSC(MatCreate(comm, &data.MatA)); + BOUT_DO_PETSC( + MatSetSizes(data.MatA, localN, localN, PETSC_DETERMINE, PETSC_DETERMINE)); + BOUT_DO_PETSC(MatSetFromOptions(data.MatA)); ////////////////////////////////////////////////// // Pre-allocate space for matrix elements @@ -223,28 +224,28 @@ LaplaceXZpetsc::LaplaceXZpetsc(Mesh* m, Options* opt, const CELL_LOC loc) } } - MatMPIAIJSetPreallocation(data.MatA, 0, d_nnz, 0, o_nnz); - MatSetUp(data.MatA); - PetscFree(d_nnz); - PetscFree(o_nnz); + BOUT_DO_PETSC(MatMPIAIJSetPreallocation(data.MatA, 0, d_nnz, 0, o_nnz)); + BOUT_DO_PETSC(MatSetUp(data.MatA)); + BOUT_DO_PETSC(PetscFree(d_nnz)); + BOUT_DO_PETSC(PetscFree(o_nnz)); ////////////////////////////////////////////////// // Declare KSP Context - KSPCreate(comm, &data.ksp); + BOUT_DO_PETSC(KSPCreate(comm, &data.ksp)); // Set KSP type - KSPSetType(data.ksp, ksptype.c_str()); + BOUT_DO_PETSC(KSPSetType(data.ksp, ksptype.c_str())); // Set KSP tolerances - KSPSetTolerances(data.ksp, rtol, atol, dtol, maxits); + BOUT_DO_PETSC(KSPSetTolerances(data.ksp, rtol, atol, dtol, maxits)); // Set KSP preconditioner PC pc; - KSPGetPC(data.ksp, &pc); - PCSetType(pc, pctype.c_str()); + BOUT_DO_PETSC(KSPGetPC(data.ksp, &pc)); + BOUT_DO_PETSC(PCSetType(pc, pctype.c_str())); #if PETSC_VERSION_GE(3, 9, 0) - PCFactorSetMatSolverType(pc, factor_package.c_str()); + BOUT_DO_PETSC(PCFactorSetMatSolverType(pc, factor_package.c_str())); #else - PCFactorSetMatSolverPackage(pc, factor_package.c_str()); + BOUT_DO_PETSC(PCFactorSetMatSolverPackage(pc, factor_package.c_str())); #endif lib.setOptionsFromInputFile(data.ksp); @@ -257,22 +258,22 @@ LaplaceXZpetsc::LaplaceXZpetsc(Mesh* m, Options* opt, const CELL_LOC loc) LaplaceXZpetsc::~LaplaceXZpetsc() { PetscBool petsc_is_finalised; - PetscFinalized(&petsc_is_finalised); + BOUT_DO_PETSC(PetscFinalized(&petsc_is_finalised)); for (auto& it : slice) { - MatDestroy(&it.MatA); + BOUT_DO_PETSC(MatDestroy(&it.MatA)); if (coefs_set) { - MatDestroy(&it.MatP); + BOUT_DO_PETSC(MatDestroy(&it.MatP)); } if (!petsc_is_finalised) { // PetscFinalize may already have destroyed this object - KSPDestroy(&it.ksp); + BOUT_DO_PETSC(KSPDestroy(&it.ksp)); } } - VecDestroy(&bs); - VecDestroy(&xs); + BOUT_DO_PETSC(VecDestroy(&bs)); + BOUT_DO_PETSC(VecDestroy(&xs)); } void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { @@ -316,7 +317,7 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { int y = it.yindex; int Istart, Iend; - MatGetOwnershipRange(it.MatA, &Istart, &Iend); + BOUT_DO_PETSC(MatGetOwnershipRange(it.MatA, &Istart, &Iend)); //////////////////////////////////////////////// // Inner X boundary (see note about BC in LaplaceXZ constructor) @@ -329,11 +330,11 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { */ for (int z = localmesh->zstart; z <= localmesh->zend; z++) { PetscScalar val = 1.0; - MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = row + (localmesh->LocalNz); // +1 in X val = -1.0; - MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); row++; } @@ -341,11 +342,11 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { // Setting BC from x0 for (int z = localmesh->zstart; z <= localmesh->zend; z++) { PetscScalar val = 1.0; - MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = row + (localmesh->LocalNz); // +1 in X val = 0.0; - MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); row++; } @@ -353,11 +354,11 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { // Setting BC from b for (int z = localmesh->zstart; z <= localmesh->zend; z++) { PetscScalar val = 1.0; - MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = row + (localmesh->LocalNz); // +1 in X val = 0.0; - MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); row++; } @@ -368,10 +369,10 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { */ for (int z = localmesh->zstart; z <= localmesh->zend; z++) { PetscScalar val = 0.5; - MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = row + (localmesh->LocalNz); // +1 in X - MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); row++; } @@ -597,15 +598,15 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { // // Set the centre (diagonal) - MatSetValues(it.MatA, 1, &row, 1, &row, &c, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &row, &c, INSERT_VALUES)); // X + 1 int col = row + (localmesh->LocalNz); - MatSetValues(it.MatA, 1, &row, 1, &col, &xp, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &col, &xp, INSERT_VALUES)); // X - 1 col = row - (localmesh->LocalNz); - MatSetValues(it.MatA, 1, &row, 1, &col, &xm, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &col, &xm, INSERT_VALUES)); // Z + 1 col = row + 1; @@ -613,30 +614,30 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { col -= localmesh->LocalNz; // Wrap around } - MatSetValues(it.MatA, 1, &row, 1, &col, &zp, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &col, &zp, INSERT_VALUES)); // X + 1, Z + 1 const int xpzp_col = col + localmesh->LocalNz; - MatSetValues(it.MatA, 1, &row, 1, &xpzp_col, &xpzp, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &xpzp_col, &xpzp, INSERT_VALUES)); // X - 1, Z + 1 const int xmzp_col = col - localmesh->LocalNz; - MatSetValues(it.MatA, 1, &row, 1, &xmzp_col, &xmzp, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &xmzp_col, &xmzp, INSERT_VALUES)); // Z - 1 col = row - 1; if (z == 0) { col += localmesh->LocalNz; // Wrap around } - MatSetValues(it.MatA, 1, &row, 1, &col, &zm, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &col, &zm, INSERT_VALUES)); // X + 1, Z - 1 const int xpzm_col = col + localmesh->LocalNz; - MatSetValues(it.MatA, 1, &row, 1, &xpzm_col, &xpzm, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &xpzm_col, &xpzm, INSERT_VALUES)); // X - 1, Z - 1 const int xmzm_col = col - localmesh->LocalNz; - MatSetValues(it.MatA, 1, &row, 1, &xmzm_col, &xmzm, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &xmzm_col, &xmzm, INSERT_VALUES)); row++; } @@ -649,11 +650,11 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { // Neumann 0 for (int z = localmesh->zstart; z <= localmesh->zend; z++) { PetscScalar val = 1.0; - MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = row - (localmesh->LocalNz); // -1 in X val = -1.0; - MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); row++; } @@ -661,11 +662,11 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { // Setting BC from x0 for (int z = localmesh->zstart; z <= localmesh->zend; z++) { PetscScalar val = 1.0; - MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = row - (localmesh->LocalNz); // -1 in X val = 0.0; - MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); row++; } @@ -673,11 +674,11 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { // Setting BC from b for (int z = localmesh->zstart; z <= localmesh->zend; z++) { PetscScalar val = 1.0; - MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = row - (localmesh->LocalNz); // -1 in X val = 0.0; - MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); row++; } @@ -686,10 +687,10 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { PetscScalar val = 0.5; for (int z = localmesh->zstart; z <= localmesh->zend; z++) { - MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = row - (localmesh->LocalNz); // -1 in X - MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(it.MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); row++; } @@ -699,8 +700,8 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { ASSERT1(row == Iend); // Check that row is currently on the last row // Assemble Matrix - MatAssemblyBegin(it.MatA, MAT_FINAL_ASSEMBLY); - MatAssemblyEnd(it.MatA, MAT_FINAL_ASSEMBLY); + BOUT_DO_PETSC(MatAssemblyBegin(it.MatA, MAT_FINAL_ASSEMBLY)); + BOUT_DO_PETSC(MatAssemblyEnd(it.MatA, MAT_FINAL_ASSEMBLY)); } // Increase reuse count @@ -714,9 +715,9 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { // Copy matrix into preconditioner if (coefs_set) { // Preconditioner already set - MatDestroy(&it.MatP); + BOUT_DO_PETSC(MatDestroy(&it.MatP)); } - MatConvert(it.MatA, MATSAME, MAT_INITIAL_MATRIX, &it.MatP); + BOUT_DO_PETSC(MatConvert(it.MatA, MATSAME, MAT_INITIAL_MATRIX, &it.MatP)); } // Set operators @@ -724,9 +725,9 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { // Note: This is a hack to force update of the preconditioner matrix #if PETSC_VERSION_GE(3, 5, 0) - KSPSetOperators(it.ksp, it.MatA, it.MatP); + KSPSetOperators(it.ksp, it.BOUT_DO_PETSC(MatA, it.MatP)); #else - KSPSetOperators(it.ksp, it.MatA, it.MatP, SAME_NONZERO_PATTERN); + KSPSetOperators(it.ksp, it.BOUT_DO_PETSC(MatA, it.MatP, SAME_NONZERO_PATTERN)); #endif } } else { @@ -736,7 +737,7 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { #if PETSC_VERSION_GE(3, 5, 0) KSPSetReusePreconditioner(it.ksp, PETSC_TRUE); #else - KSPSetOperators(it.ksp, it.MatA, it.MatP, SAME_PRECONDITIONER); + KSPSetOperators(it.ksp, it.BOUT_DO_PETSC(MatA, it.MatP, SAME_PRECONDITIONER)); #endif } } @@ -785,7 +786,7 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { ////////////////////////// // Load initial guess x0 into xs and rhs into bs int Istart, Iend; - MatGetOwnershipRange(it.MatA, &Istart, &Iend); + BOUT_DO_PETSC(MatGetOwnershipRange(it.MatA, &Istart, &Iend)); // Starting index int ind = Istart; From 6fb890c7d517350fb37f864fde4a7767348912f3 Mon Sep 17 00:00:00 2001 From: David Bold Date: Tue, 3 Nov 2026 10:51:18 +0100 Subject: [PATCH 3/6] Check more return values --- .../laplacexy/impls/petsc/laplacexy-petsc.cxx | 92 +++++++++---------- 1 file changed, 46 insertions(+), 46 deletions(-) diff --git a/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx b/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx index 540647d01d..7b2234f65c 100644 --- a/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx +++ b/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx @@ -90,15 +90,15 @@ LaplaceXYpetsc::LaplaceXYpetsc(Mesh* m, Options* opt, const CELL_LOC loc) const int localN = localSize(); // Create Vectors - VecCreate(comm, &xs); - VecSetSizes(xs, localN, PETSC_DETERMINE); - VecSetFromOptions(xs); - VecDuplicate(xs, &bs); + BOUT_DO_PETSC(VecCreate(comm, &xs)); + BOUT_DO_PETSC(VecSetSizes(xs, localN, PETSC_DETERMINE)); + BOUT_DO_PETSC(VecSetFromOptions(xs)); + BOUT_DO_PETSC(VecDuplicate(xs, &bs)); // Set size of Matrix on each processor to localN x localN - MatCreate(comm, &MatA); - MatSetSizes(MatA, localN, localN, PETSC_DETERMINE, PETSC_DETERMINE); - MatSetFromOptions(MatA); + BOUT_DO_PETSC(MatCreate(comm, &MatA)); + BOUT_DO_PETSC(MatSetSizes(MatA, localN, localN, PETSC_DETERMINE, PETSC_DETERMINE)); + BOUT_DO_PETSC(MatSetFromOptions(MatA)); ////////////////////////////////////////////////// // Specify local indices. This creates a mapping @@ -187,8 +187,8 @@ LaplaceXYpetsc::LaplaceXYpetsc(Mesh* m, Options* opt, const CELL_LOC loc) setPreallocationFiniteDifference(d_nnz, o_nnz); } // Pre-allocate - MatMPIAIJSetPreallocation(MatA, 0, d_nnz, 0, o_nnz); - MatSetUp(MatA); + BOUT_DO_PETSC(MatMPIAIJSetPreallocation(MatA, 0, d_nnz, 0, o_nnz)); + BOUT_DO_PETSC(MatSetUp(MatA)); BOUT_DO_PETSC(PetscFree(d_nnz)); BOUT_DO_PETSC(PetscFree(o_nnz)); @@ -196,7 +196,7 @@ LaplaceXYpetsc::LaplaceXYpetsc(Mesh* m, Options* opt, const CELL_LOC loc) // Determine which row/columns of the matrix are locally owned int Istart = 0; int Iend = 0; - MatGetOwnershipRange(MatA, &Istart, &Iend); + BOUT_DO_PETSC(MatGetOwnershipRange(MatA, &Istart, &Iend)); // Convert indexXY from local index to global index indexXY += Istart; @@ -209,14 +209,14 @@ LaplaceXYpetsc::LaplaceXYpetsc(Mesh* m, Options* opt, const CELL_LOC loc) // Set up KSP // Declare KSP Context - KSPCreate(comm, &ksp); + BOUT_DO_PETSC(KSPCreate(comm, &ksp)); // Configure Linear Solver const bool direct = (*opt)["direct"].doc("Use a direct LU solver").withDefault(false); if (direct) { - KSPGetPC(ksp, &pc); + BOUT_DO_PETSC(KSPGetPC(ksp, &pc)); PCSetType(pc, PCLU); #if PETSC_VERSION_GE(3, 9, 0) PCFactorSetMatSolverType(pc, "mumps"); @@ -248,12 +248,12 @@ LaplaceXYpetsc::LaplaceXYpetsc(Mesh* m, Options* opt, const CELL_LOC loc) const std::string pctype = (*opt)["pctype"].doc("Preconditioner type").withDefault("none"); - KSPSetType(ksp, ksptype.c_str()); - KSPSetTolerances(ksp, rtol, atol, dtol, maxits); + BOUT_DO_PETSC(KSPSetType(ksp, ksptype.c_str())); + BOUT_DO_PETSC(KSPSetTolerances(ksp, rtol, atol, dtol, maxits)); - KSPSetInitialGuessNonzero(ksp, static_cast(true)); + BOUT_DO_PETSC(KSPSetInitialGuessNonzero(ksp, static_cast(true))); - KSPGetPC(ksp, &pc); + BOUT_DO_PETSC(KSPGetPC(ksp, &pc)); PCSetType(pc, pctype.c_str()); if (pctype == "shell") { @@ -264,9 +264,9 @@ LaplaceXYpetsc::LaplaceXYpetsc(Mesh* m, Options* opt, const CELL_LOC loc) const bool rightprec = (*opt)["rightprec"].doc("Use right preconditioning?").withDefault(true); if (rightprec) { - KSPSetPCSide(ksp, PC_RIGHT); // Right preconditioning + BOUT_DO_PETSC(KSPSetPCSide(ksp, PC_RIGHT)); // Right preconditioning } else { - KSPSetPCSide(ksp, PC_LEFT); // Left preconditioning + BOUT_DO_PETSC(KSPSetPCSide(ksp, PC_LEFT)); // Left preconditioning } } } @@ -578,10 +578,10 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { for (int y = localmesh->ystart; y <= localmesh->yend; y++) { int row = globalIndex(localmesh->xstart - 1, y); PetscScalar val = 0.5; - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = globalIndex(localmesh->xstart, y); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); // Preconditioner bcoef(y - localmesh->ystart, 0) = 0.5; @@ -594,11 +594,11 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { for (int y = localmesh->ystart; y <= localmesh->yend; y++) { int row = globalIndex(localmesh->xstart - 1, y); PetscScalar val = 1.0; - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = globalIndex(localmesh->xstart, y); val = -1.0; - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); // Preconditioner bcoef(y - localmesh->ystart, 0) = 1.0; @@ -612,10 +612,10 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { for (int y = localmesh->ystart; y <= localmesh->yend; y++) { int row = globalIndex(localmesh->xend + 1, y); PetscScalar val = 0.5; - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = globalIndex(localmesh->xend, y); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); // Preconditioner acoef(y - localmesh->ystart, localmesh->xend + 1 - xstart) = 0.5; @@ -632,10 +632,10 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { } int row = globalIndex(it.ind, localmesh->ystart - 1); PetscScalar val = 0.5; - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = globalIndex(it.ind, localmesh->ystart); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } for (RangeIterator it = localmesh->iterateBndryUpperY(); !it.isDone(); it++) { @@ -645,10 +645,10 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { } int row = globalIndex(it.ind, localmesh->yend + 1); PetscScalar val = 0.5; - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); int col = globalIndex(it.ind, localmesh->yend); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } } else if (y_bndry == "neumann") { // Neumann on Y boundaries @@ -659,12 +659,12 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { } int row = globalIndex(it.ind, localmesh->ystart - 1); PetscScalar val = 1.0; - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); val = -1.0; int col = globalIndex(it.ind, localmesh->ystart); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } for (RangeIterator it = localmesh->iterateBndryUpperY(); !it.isDone(); it++) { @@ -674,11 +674,11 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { } int row = globalIndex(it.ind, localmesh->yend + 1); PetscScalar val = 1.0; - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); val = -1.0; int col = globalIndex(it.ind, localmesh->yend); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } } else if (y_bndry == "free_o3") { // 'free_o3' extrapolating boundary condition on Y boundaries @@ -689,19 +689,19 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { } int row = globalIndex(it.ind, localmesh->ystart - 1); PetscScalar val = 1.0; - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); val = -3.0; int col = globalIndex(it.ind, localmesh->ystart); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); val = 3.0; col = globalIndex(it.ind, localmesh->ystart + 1); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); val = -1.0; col = globalIndex(it.ind, localmesh->ystart + 2); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } for (RangeIterator it = localmesh->iterateBndryUpperY(); !it.isDone(); it++) { @@ -711,19 +711,19 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { } int row = globalIndex(it.ind, localmesh->yend + 1); PetscScalar val = 1.0; - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); val = -3.0; int col = globalIndex(it.ind, localmesh->yend); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); val = 3.0; col = globalIndex(it.ind, localmesh->yend - 1); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); val = -1.0; col = globalIndex(it.ind, localmesh->yend - 2); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } } else { throw BoutException("Unsupported option for y_bndry"); @@ -741,11 +741,11 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { // f(xs-1,ys-1) = f(xs-1,ys) PetscScalar val = 1.0; int row = globalIndex(localmesh->xstart - 1, localmesh->ystart - 1); - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); val = -1.0; int col = globalIndex(localmesh->xstart - 1, localmesh->ystart); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } else if (y_bndry == "free_o3" or y_bndry == "dirichlet") { // 'free_o3' extrapolating boundary condition on Y boundaries // f(xs-1,ys-1) = 3*f(xs-1,ys) - 3*f(xs-1,ys+1) + f(xs-1,ys+2) @@ -754,19 +754,19 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { // what value to pass for the corner PetscScalar val = 1.0; int row = globalIndex(localmesh->xstart - 1, localmesh->ystart - 1); - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); val = -3.0; int col = globalIndex(localmesh->xstart - 1, localmesh->ystart); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); val = 3.0; col = globalIndex(localmesh->xstart - 1, localmesh->ystart + 1); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); val = -1.0; col = globalIndex(localmesh->xstart - 1, localmesh->ystart + 2); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } else { throw BoutException("Unsupported option for y_bndry"); } From 6778df3cdefacc021bdcd92624069c5275000810 Mon Sep 17 00:00:00 2001 From: David Bold Date: Tue, 3 Nov 2026 10:54:01 +0100 Subject: [PATCH 4/6] Check more return values from PETSc --- .../laplacexy/impls/petsc/laplacexy-petsc.cxx | 244 +++++++++--------- .../laplacexz/impls/petsc/laplacexz-petsc.cxx | 81 +++--- 2 files changed, 165 insertions(+), 160 deletions(-) diff --git a/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx b/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx index 7b2234f65c..f578d7851c 100644 --- a/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx +++ b/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx @@ -777,11 +777,11 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { // f(xe+1,ys-1) = f(xe+1,ys) PetscScalar val = 1.0; int row = globalIndex(localmesh->xend + 1, localmesh->ystart - 1); - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); val = -1.0; int col = globalIndex(localmesh->xend + 1, localmesh->ystart); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } else if (y_bndry == "free_o3" or y_bndry == "dirichlet") { // 'free_o3' extrapolating boundary condition on Y boundaries // f(xe+1,ys-1) = 3*f(xe+1,ys) - 3*f(xe+1,ys+1) + f(xe+1,ys+2) @@ -790,19 +790,19 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { // what value to pass for the corner PetscScalar val = 1.0; int row = globalIndex(localmesh->xend + 1, localmesh->ystart - 1); - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); val = -3.0; int col = globalIndex(localmesh->xend + 1, localmesh->ystart); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); val = 3.0; col = globalIndex(localmesh->xend + 1, localmesh->ystart + 1); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); val = -1.0; col = globalIndex(localmesh->xend + 1, localmesh->ystart + 2); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } else { throw BoutException("Unsupported option for y_bndry"); } @@ -815,11 +815,11 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { // f(xs-1,ys-1) = f(xs-1,ys) PetscScalar val = 1.0; int row = globalIndex(localmesh->xstart - 1, localmesh->yend + 1); - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); val = -1.0; int col = globalIndex(localmesh->xstart - 1, localmesh->yend); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } else if (y_bndry == "free_o3" or y_bndry == "dirichlet") { // 'free_o3' extrapolating boundary condition on Y boundaries // f(xs-1,ys-1) = 3*f(xs-1,ys) - 3*f(xs-1,ys+1) + f(xs-1,ys+2) @@ -828,19 +828,19 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { // what value to pass for the corner PetscScalar val = 1.0; int row = globalIndex(localmesh->xstart - 1, localmesh->yend + 1); - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); val = -3.0; int col = globalIndex(localmesh->xstart - 1, localmesh->yend); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); val = 3.0; col = globalIndex(localmesh->xstart - 1, localmesh->yend - 1); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); val = -1.0; col = globalIndex(localmesh->xstart - 1, localmesh->yend - 2); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } else { throw BoutException("Unsupported option for y_bndry"); } @@ -851,11 +851,11 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { // f(xe+1,ys-1) = f(xe+1,ys) PetscScalar val = 1.0; int row = globalIndex(localmesh->xend + 1, localmesh->yend + 1); - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); val = -1.0; int col = globalIndex(localmesh->xend + 1, localmesh->yend); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } else if (y_bndry == "free_o3" or y_bndry == "dirichlet") { // 'free_o3' extrapolating boundary condition on Y boundaries // f(xe+1,ys-1) = 3*f(xe+1,ys) - 3*f(xe+1,ys+1) + f(xe+1,ys+2) @@ -864,19 +864,19 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { // what value to pass for the corner PetscScalar val = 1.0; int row = globalIndex(localmesh->xend + 1, localmesh->yend + 1); - MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &val, INSERT_VALUES)); val = -3.0; int col = globalIndex(localmesh->xend + 1, localmesh->yend); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); val = 3.0; col = globalIndex(localmesh->xend + 1, localmesh->yend - 1); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); val = -1.0; col = globalIndex(localmesh->xend + 1, localmesh->yend - 2); - MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &val, INSERT_VALUES)); } else { throw BoutException("Unsupported option for y_bndry"); } @@ -885,14 +885,15 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { } // Assemble Matrix - MatAssemblyBegin(MatA, MAT_FINAL_ASSEMBLY); - MatAssemblyEnd(MatA, MAT_FINAL_ASSEMBLY); + BOUT_DO_PETSC(MatAssemblyBegin(MatA, MAT_FINAL_ASSEMBLY)); + BOUT_DO_PETSC(MatAssemblyEnd(MatA, MAT_FINAL_ASSEMBLY)); // Set the operator #if PETSC_VERSION_GE(3, 5, 0) - KSPSetOperators(ksp, MatA, MatA); + BOUT_DO_PETSC(KSPSetOperators(ksp, BOUT_DO_PETSC(MatA, MatA))); #else - KSPSetOperators(ksp, MatA, MatA, DIFFERENT_NONZERO_PATTERN); + BOUT_DO_PETSC( + KSPSetOperators(ksp, BOUT_DO_PETSC(MatA, MatA, DIFFERENT_NONZERO_PATTERN))); #endif // Set coefficients for preconditioner @@ -981,24 +982,24 @@ void LaplaceXYpetsc::setMatrixElementsFiniteVolume(const Field2D& A, const Field int row = globalIndex(x, y); // Set the centre (diagonal) - MatSetValues(MatA, 1, &row, 1, &row, &c, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &c, INSERT_VALUES)); // X + 1 int col = globalIndex(x + 1, y); - MatSetValues(MatA, 1, &row, 1, &col, &xp, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &xp, INSERT_VALUES)); // X - 1 col = globalIndex(x - 1, y); - MatSetValues(MatA, 1, &row, 1, &col, &xm, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &xm, INSERT_VALUES)); if (include_y_derivs) { // Y + 1 col = globalIndex(x, y + 1); - MatSetValues(MatA, 1, &row, 1, &col, &yp, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &yp, INSERT_VALUES)); // Y - 1 col = globalIndex(x, y - 1); - MatSetValues(MatA, 1, &row, 1, &col, &ym, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &ym, INSERT_VALUES)); } } } @@ -1113,40 +1114,40 @@ void LaplaceXYpetsc::setMatrixElementsFiniteDifference(const Field2D& A, int row = globalIndex(x, y); // Set the centre (diagonal) - MatSetValues(MatA, 1, &row, 1, &row, &c, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &row, &c, INSERT_VALUES)); // X + 1 int col = globalIndex(x + 1, y); - MatSetValues(MatA, 1, &row, 1, &col, &xp, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &xp, INSERT_VALUES)); // X - 1 col = globalIndex(x - 1, y); - MatSetValues(MatA, 1, &row, 1, &col, &xm, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &xm, INSERT_VALUES)); if (include_y_derivs) { // Y + 1 col = globalIndex(x, y + 1); - MatSetValues(MatA, 1, &row, 1, &col, &yp, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &yp, INSERT_VALUES)); // Y - 1 col = globalIndex(x, y - 1); - MatSetValues(MatA, 1, &row, 1, &col, &ym, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &ym, INSERT_VALUES)); // X + 1, Y + 1 col = globalIndex(x + 1, y + 1); - MatSetValues(MatA, 1, &row, 1, &col, &xpyp, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &xpyp, INSERT_VALUES)); // X + 1, Y - 1 col = globalIndex(x + 1, y - 1); - MatSetValues(MatA, 1, &row, 1, &col, &xpym, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &xpym, INSERT_VALUES)); // X - 1, Y + 1 col = globalIndex(x - 1, y + 1); - MatSetValues(MatA, 1, &row, 1, &col, &xmyp, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &xmyp, INSERT_VALUES)); // X - 1, Y - 1 col = globalIndex(x - 1, y - 1); - MatSetValues(MatA, 1, &row, 1, &col, &xmym, INSERT_VALUES); + BOUT_DO_PETSC(MatSetValues(MatA, 1, &row, 1, &col, &xmym, INSERT_VALUES)); } } } @@ -1158,12 +1159,12 @@ LaplaceXYpetsc::~LaplaceXYpetsc() { if (!is_finalised) { // PetscFinalize may already have destroyed this object - KSPDestroy(&ksp); + BOUT_DO_PETSC(KSPDestroy(&ksp)); } - VecDestroy(&xs); - VecDestroy(&bs); - MatDestroy(&MatA); + BOUT_DO_PETSC(VecDestroy(&xs)); + BOUT_DO_PETSC(VecDestroy(&bs)); + BOUT_DO_PETSC(MatDestroy(&MatA)); } Field2D LaplaceXYpetsc::solve(const Field2D& rhs, const Field2D& x0) { @@ -1181,10 +1182,10 @@ Field2D LaplaceXYpetsc::solve(const Field2D& rhs, const Field2D& x0) { int ind = globalIndex(x, y); PetscScalar val = x0(x, y); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); val = rhs(x, y); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } @@ -1195,22 +1196,23 @@ Field2D LaplaceXYpetsc::solve(const Field2D& rhs, const Field2D& x0) { } // Assemble RHS Vector - VecAssemblyBegin(bs); - VecAssemblyEnd(bs); + BOUT_DO_PETSC(VecAssemblyBegin(bs)); + BOUT_DO_PETSC(VecAssemblyEnd(bs)); // Assemble Trial Solution Vector - VecAssemblyBegin(xs); - VecAssemblyEnd(xs); + BOUT_DO_PETSC(VecAssemblyBegin(xs)); + BOUT_DO_PETSC(VecAssemblyEnd(xs)); // Solve the system - KSPSolve(ksp, bs, xs); + BOUT_DO_PETSC(KSPSolve(ksp, bs, xs)); - KSPConvergedReason reason; - KSPGetConvergedReason(ksp, &reason); + BOUT_DO_PETSC(KSPConvergedReason reason); + BOUT_DO_PETSC(KSPGetConvergedReason(ksp, &reason)); if (reason <= 0) { - throw BoutException("LaplaceXYpetsc failed to converge. Reason {} ({:d})", - KSPConvergedReasons[reason], static_cast(reason)); + throw BoutException( + "LaplaceXYpetsc failed to converge. Reason {} ({:d})", + BOUT_DO_PETSC(KSPConvergedReasons[reason], static_cast(reason))); } if (save_performance) { @@ -1218,7 +1220,7 @@ Field2D LaplaceXYpetsc::solve(const Field2D& rhs, const Field2D& x0) { n_calls++; int iterations = 0; - KSPGetIterationNumber(ksp, &iterations); + BOUT_DO_PETSC(KSPGetIterationNumber(ksp, &iterations)); average_iterations = BoutReal(n_calls - 1) / BoutReal(n_calls) * average_iterations + BoutReal(iterations) / BoutReal(n_calls); @@ -1236,7 +1238,7 @@ Field2D LaplaceXYpetsc::solve(const Field2D& rhs, const Field2D& x0) { int ind = globalIndex(x, y); PetscScalar val; - VecGetValues(xs, 1, &ind, &val); + BOUT_DO_PETSC(VecGetValues(xs, 1, &ind, &val)); result(x, y) = val; } } @@ -1246,7 +1248,7 @@ Field2D LaplaceXYpetsc::solve(const Field2D& rhs, const Field2D& x0) { for (int y = localmesh->ystart; y <= localmesh->yend; y++) { int ind = globalIndex(localmesh->xstart - 1, y); PetscScalar val; - VecGetValues(xs, 1, &ind, &val); + BOUT_DO_PETSC(VecGetValues(xs, 1, &ind, &val)); for (int x = localmesh->xstart - 1; x >= 0; x--) { result(x, y) = val; } @@ -1258,7 +1260,7 @@ Field2D LaplaceXYpetsc::solve(const Field2D& rhs, const Field2D& x0) { for (int y = localmesh->ystart; y <= localmesh->yend; y++) { int ind = globalIndex(localmesh->xend + 1, y); PetscScalar val; - VecGetValues(xs, 1, &ind, &val); + BOUT_DO_PETSC(VecGetValues(xs, 1, &ind, &val)); for (int x = localmesh->xend + 1; x < localmesh->LocalNx; x++) { result(x, y) = val; } @@ -1278,7 +1280,7 @@ Field2D LaplaceXYpetsc::solve(const Field2D& rhs, const Field2D& x0) { } int ind = globalIndex(it.ind, localmesh->ystart - 1); PetscScalar val; - VecGetValues(xs, 1, &ind, &val); + BOUT_DO_PETSC(VecGetValues(xs, 1, &ind, &val)); for (int y = localmesh->ystart - 1; y >= 0; y--) { result(it.ind, y) = val; } @@ -1297,7 +1299,7 @@ Field2D LaplaceXYpetsc::solve(const Field2D& rhs, const Field2D& x0) { } int ind = globalIndex(it.ind, localmesh->yend + 1); PetscScalar val; - VecGetValues(xs, 1, &ind, &val); + BOUT_DO_PETSC(VecGetValues(xs, 1, &ind, &val)); for (int y = localmesh->yend + 1; y < localmesh->LocalNy; y++) { result(it.ind, y) = val; } @@ -1315,10 +1317,10 @@ void LaplaceXYpetsc::solveFiniteVolume(const Field2D& x0) { int ind = globalIndex(localmesh->xstart - 1, y); PetscScalar val = x0(localmesh->xstart - 1, y); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); val = 0.5 * (x0(localmesh->xstart - 1, y) + x0(localmesh->xstart, y)); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } else { // Inner X boundary (Neumann) @@ -1326,10 +1328,10 @@ void LaplaceXYpetsc::solveFiniteVolume(const Field2D& x0) { int ind = globalIndex(localmesh->xstart - 1, y); PetscScalar val = x0(localmesh->xstart - 1, y); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); val = 0.0; // x0(localmesh->xstart-1,y) - x0(localmesh->xstart,y); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } } @@ -1340,10 +1342,10 @@ void LaplaceXYpetsc::solveFiniteVolume(const Field2D& x0) { int ind = globalIndex(localmesh->xend + 1, y); PetscScalar val = x0(localmesh->xend + 1, y); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); val = 0.5 * (x0(localmesh->xend, y) + x0(localmesh->xend + 1, y)); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } @@ -1356,10 +1358,10 @@ void LaplaceXYpetsc::solveFiniteVolume(const Field2D& x0) { int ind = globalIndex(it.ind, localmesh->ystart - 1); PetscScalar val = x0(it.ind, localmesh->ystart - 1); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); val = 0.5 * (x0(it.ind, localmesh->ystart - 1) + x0(it.ind, localmesh->ystart)); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } for (RangeIterator it = localmesh->iterateBndryUpperY(); !it.isDone(); it++) { @@ -1370,10 +1372,10 @@ void LaplaceXYpetsc::solveFiniteVolume(const Field2D& x0) { int ind = globalIndex(it.ind, localmesh->yend + 1); PetscScalar val = x0(it.ind, localmesh->yend + 1); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); val = 0.5 * (x0(it.ind, localmesh->yend + 1) + x0(it.ind, localmesh->yend)); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } else if (y_bndry == "neumann" or y_bndry == "free_o3") { // Y boundaries Neumann @@ -1385,10 +1387,10 @@ void LaplaceXYpetsc::solveFiniteVolume(const Field2D& x0) { int ind = globalIndex(it.ind, localmesh->ystart - 1); PetscScalar val = x0(it.ind, localmesh->ystart - 1); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } for (RangeIterator it = localmesh->iterateBndryUpperY(); !it.isDone(); it++) { @@ -1399,10 +1401,10 @@ void LaplaceXYpetsc::solveFiniteVolume(const Field2D& x0) { int ind = globalIndex(it.ind, localmesh->yend + 1); PetscScalar val = x0(it.ind, localmesh->yend + 1); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } else { throw BoutException("Unsupported option for y_bndry"); @@ -1422,11 +1424,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { // For the boundary value of the initial guess, use the value that would be set by // applying the boundary condition to the initial guess PetscScalar val = 2. * x0(localmesh->xstart - 1, y) - x0(localmesh->xstart, y); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Pass the value from boundary cell of x0 as the boundary condition to the rhs val = x0(localmesh->xstart - 1, y); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } else { // Inner X boundary (Neumann) @@ -1436,10 +1438,10 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { // For the boundary value of the initial guess, use the value that would be set by // applying the boundary condition to the initial guess PetscScalar val = x0(localmesh->xstart, y); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } } @@ -1452,11 +1454,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { // For the boundary value of the initial guess, use the value that would be set by // applying the boundary condition to the initial guess PetscScalar val = 2. * x0(localmesh->xend + 1, y) - x0(localmesh->xend, y); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Pass the value from boundary cell of x0 as the boundary condition to the rhs val = x0(localmesh->xend + 1, y); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } @@ -1472,11 +1474,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { // applying the boundary condition to the initial guess PetscScalar val = 2. * x0(it.ind, localmesh->ystart - 1) - x0(it.ind, localmesh->ystart); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Pass the value from boundary cell of x0 as the boundary condition to the rhs val = x0(it.ind, localmesh->ystart - 1); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } for (RangeIterator it = localmesh->iterateBndryUpperY(); !it.isDone(); it++) { @@ -1490,11 +1492,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { // applying the boundary condition to the initial guess PetscScalar val = 2. * x0(it.ind, localmesh->yend + 1) - x0(it.ind, localmesh->yend); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Pass the value from boundary cell of x0 as the boundary condition to the rhs val = x0(it.ind, localmesh->yend + 1); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } // Use free_o3 for the corner boundary cells @@ -1507,11 +1509,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { PetscScalar val = 3. * x0(localmesh->xstart - 1, localmesh->ystart) - 3. * x0(localmesh->xstart - 1, localmesh->ystart + 1) + x0(localmesh->xstart - 1, localmesh->ystart + 2); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } if (localmesh->lastX()) { int ind = globalIndex(localmesh->xend + 1, localmesh->ystart - 1); @@ -1521,11 +1523,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { PetscScalar val = 3. * x0(localmesh->xend + 1, localmesh->ystart) - 3. * x0(localmesh->xend + 1, localmesh->ystart + 1) + x0(localmesh->xend + 1, localmesh->ystart + 2); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } if (localmesh->hasBndryUpperY()) { @@ -1537,11 +1539,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { PetscScalar val = 3. * x0(localmesh->xstart - 1, localmesh->yend) - 3. * x0(localmesh->xstart - 1, localmesh->yend - 1) + x0(localmesh->xstart - 1, localmesh->yend - 2); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } if (localmesh->lastX()) { int ind = globalIndex(localmesh->xend + 1, localmesh->yend + 1); @@ -1551,11 +1553,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { PetscScalar val = 3. * x0(localmesh->xend + 1, localmesh->yend) - 3. * x0(localmesh->xend + 1, localmesh->yend - 1) + x0(localmesh->xend + 1, localmesh->yend - 2); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } } else if (y_bndry == "neumann") { @@ -1570,11 +1572,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { // For the boundary value of the initial guess, use the value that would be set by // applying the boundary condition to the initial guess PetscScalar val = x0(it.ind, localmesh->ystart); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } if (localmesh->hasBndryLowerY()) { if (localmesh->firstX()) { @@ -1583,11 +1585,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { // For the boundary value of the initial guess, use the value that would be set by // applying the boundary condition to the initial guess PetscScalar val = x0(localmesh->xstart - 1, localmesh->ystart); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } if (localmesh->lastX()) { int ind = globalIndex(localmesh->xend + 1, localmesh->ystart - 1); @@ -1595,11 +1597,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { // For the boundary value of the initial guess, use the value that would be set by // applying the boundary condition to the initial guess PetscScalar val = x0(localmesh->xend + 1, localmesh->ystart); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } @@ -1613,11 +1615,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { // For the boundary value of the initial guess, use the value that would be set by // applying the boundary condition to the initial guess PetscScalar val = x0(it.ind, localmesh->yend); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } if (localmesh->hasBndryUpperY()) { if (localmesh->firstX()) { @@ -1626,11 +1628,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { // For the boundary value of the initial guess, use the value that would be set by // applying the boundary condition to the initial guess PetscScalar val = x0(localmesh->xstart - 1, localmesh->yend); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } if (localmesh->lastX()) { int ind = globalIndex(localmesh->xend + 1, localmesh->yend + 1); @@ -1638,11 +1640,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { // For the boundary value of the initial guess, use the value that would be set by // applying the boundary condition to the initial guess PetscScalar val = x0(localmesh->xend + 1, localmesh->yend); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } } else if (y_bndry == "free_o3") { @@ -1659,11 +1661,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { PetscScalar val = 3. * x0(it.ind, localmesh->ystart) - 3. * x0(it.ind, localmesh->ystart + 1) + x0(it.ind, localmesh->ystart + 2); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } if (localmesh->hasBndryLowerY()) { if (localmesh->firstX()) { @@ -1674,11 +1676,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { PetscScalar val = 3. * x0(localmesh->xstart - 1, localmesh->ystart) - 3. * x0(localmesh->xstart - 1, localmesh->ystart + 1) + x0(localmesh->xstart - 1, localmesh->ystart + 2); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } if (localmesh->lastX()) { int ind = globalIndex(localmesh->xend + 1, localmesh->ystart - 1); @@ -1688,11 +1690,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { PetscScalar val = 3. * x0(localmesh->xend + 1, localmesh->ystart) - 3. * x0(localmesh->xend + 1, localmesh->ystart + 1) + x0(localmesh->xend + 1, localmesh->ystart + 2); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } @@ -1708,11 +1710,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { PetscScalar val = 3. * x0(it.ind, localmesh->yend) - 3. * x0(it.ind, localmesh->yend - 1) + x0(it.ind, localmesh->yend - 2); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } if (localmesh->hasBndryUpperY()) { if (localmesh->firstX()) { @@ -1723,11 +1725,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { PetscScalar val = 3. * x0(localmesh->xstart - 1, localmesh->yend) - 3. * x0(localmesh->xstart - 1, localmesh->yend - 1) + x0(localmesh->xstart - 1, localmesh->yend - 2); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } if (localmesh->lastX()) { int ind = globalIndex(localmesh->xend + 1, localmesh->yend + 1); @@ -1737,11 +1739,11 @@ void LaplaceXYpetsc::solveFiniteDifference(const Field2D& x0) { PetscScalar val = 3. * x0(localmesh->xend + 1, localmesh->yend) - 3. * x0(localmesh->xend + 1, localmesh->yend - 1) + x0(localmesh->xend + 1, localmesh->yend - 2); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Set the value for the rhs of the boundary condition to zero val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); } } } else { @@ -1764,8 +1766,8 @@ int LaplaceXYpetsc::precon(Vec input, Vec result) { } const int ind = globalIndex(itdwn.ind, localmesh->ystart - 1); PetscScalar val; - VecGetValues(input, 1, &ind, &val); - VecSetValues(result, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecGetValues(input, 1, &ind, &val)); + BOUT_DO_PETSC(VecSetValues(result, 1, &ind, &val, INSERT_VALUES)); } for (auto itup = localmesh->iterateBndryUpperY(); !itup.isDone(); itup++) { // Should not go into corner cells, LaplaceXYpetsc stencil does not include them @@ -1774,8 +1776,8 @@ int LaplaceXYpetsc::precon(Vec input, Vec result) { } const int ind = globalIndex(itup.ind, localmesh->yend + 1); PetscScalar val; - VecGetValues(input, 1, &ind, &val); - VecSetValues(result, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecGetValues(input, 1, &ind, &val)); + BOUT_DO_PETSC(VecSetValues(result, 1, &ind, &val, INSERT_VALUES)); } // Load vector x into bvals array @@ -1783,7 +1785,7 @@ int LaplaceXYpetsc::precon(Vec input, Vec result) { for (int y = localmesh->ystart; y <= localmesh->yend; y++) { const int ind = globalIndex(x, y); PetscScalar val; - VecGetValues(input, 1, &ind, &val); + BOUT_DO_PETSC(VecGetValues(input, 1, &ind, &val)); bvals(y - localmesh->ystart, x - xstart) = val; } } @@ -1796,11 +1798,11 @@ int LaplaceXYpetsc::precon(Vec input, Vec result) { for (int y = localmesh->ystart; y <= localmesh->yend; y++) { const int ind = globalIndex(x, y); PetscScalar val = xvals(y - localmesh->ystart, x - xstart); - VecSetValues(result, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(result, 1, &ind, &val, INSERT_VALUES)); } } - VecAssemblyBegin(result); - VecAssemblyEnd(result); + BOUT_DO_PETSC(VecAssemblyBegin(result)); + BOUT_DO_PETSC(VecAssemblyEnd(result)); return 0; } diff --git a/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx b/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx index 7538c9d3f2..0b10802937 100644 --- a/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx +++ b/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx @@ -167,10 +167,10 @@ LaplaceXZpetsc::LaplaceXZpetsc(Mesh* m, Options* opt, const CELL_LOC loc) } // Create Vectors - VecCreate(comm, &xs); - VecSetSizes(xs, localN, PETSC_DETERMINE); - VecSetFromOptions(xs); - VecDuplicate(xs, &bs); + BOUT_DO_PETSC(VecCreate(comm, &xs)); + BOUT_DO_PETSC(VecSetSizes(xs, localN, PETSC_DETERMINE)); + BOUT_DO_PETSC(VecSetFromOptions(xs)); + BOUT_DO_PETSC(VecDuplicate(xs, &bs)); for (int y = localmesh->ystart; y <= localmesh->yend; y++) { YSlice data; @@ -725,9 +725,10 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { // Note: This is a hack to force update of the preconditioner matrix #if PETSC_VERSION_GE(3, 5, 0) - KSPSetOperators(it.ksp, it.BOUT_DO_PETSC(MatA, it.MatP)); + BOUT_DO_PETSC(KSPSetOperators(it.ksp, it.BOUT_DO_PETSC(MatA, it.MatP))); #else - KSPSetOperators(it.ksp, it.BOUT_DO_PETSC(MatA, it.MatP, SAME_NONZERO_PATTERN)); + BOUT_DO_PETSC( + KSPSetOperators(it.ksp, it.BOUT_DO_PETSC(MatA, it.MatP, SAME_NONZERO_PATTERN))); #endif } } else { @@ -735,9 +736,10 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { /// Reuse the preconditioner, even if the operator changes #if PETSC_VERSION_GE(3, 5, 0) - KSPSetReusePreconditioner(it.ksp, PETSC_TRUE); + BOUT_DO_PETSC(KSPSetReusePreconditioner(it.ksp, PETSC_TRUE)); #else - KSPSetOperators(it.ksp, it.BOUT_DO_PETSC(MatA, it.MatP, SAME_PRECONDITIONER)); + BOUT_DO_PETSC( + KSPSetOperators(it.ksp, it.BOUT_DO_PETSC(MatA, it.MatP, SAME_PRECONDITIONER))); #endif } } @@ -781,7 +783,7 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { int y = it.yindex; /// Specify non-zero starting guess for solution (from input x0) - KSPSetInitialGuessNonzero(it.ksp, PETSC_TRUE); + BOUT_DO_PETSC(KSPSetInitialGuessNonzero(it.ksp, PETSC_TRUE)); ////////////////////////// // Load initial guess x0 into xs and rhs into bs @@ -798,11 +800,11 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { for (int z = localmesh->zstart; z <= localmesh->zend; z++) { // Setting the initial guess x0 PetscScalar val = x0(localmesh->xstart - 1, y, z); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Setting the solution b val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); ind++; } } else if (inner_boundary_flags & INVERT_SET) { @@ -810,11 +812,11 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { for (int z = localmesh->zstart; z <= localmesh->zend; z++) { // Setting the initial guess x0 PetscScalar val = x0(localmesh->xstart - 1, y, z); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Setting the solution b val = x0(localmesh->xstart, y, z); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); ind++; } } else if (inner_boundary_flags & INVERT_RHS) { @@ -822,11 +824,11 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { for (int z = localmesh->zstart; z <= localmesh->zend; z++) { // Setting the initial guess x0 PetscScalar val = x0(localmesh->xstart - 1, y, z); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Setting the solution b val = b(localmesh->xstart, y, z); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); ind++; } } else { @@ -834,11 +836,11 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { for (int z = localmesh->zstart; z <= localmesh->zend; z++) { // Setting the initial guess x0 PetscScalar val = x0(localmesh->xstart - 1, y, z); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Setting the solution b val = x0(localmesh->xstart - 1, y, z) - x0(localmesh->xstart, y, z); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); ind++; } } @@ -848,10 +850,10 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { for (int x = localmesh->xstart; x <= localmesh->xend; x++) { for (int z = localmesh->zstart; z <= localmesh->zend; z++) { PetscScalar val = x0(x, y, z); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); val = b(x, y, z); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); ind++; } } @@ -863,11 +865,11 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { for (int z = localmesh->zstart; z <= localmesh->zend; z++) { // Setting the initial guess x0 PetscScalar val = x0(localmesh->xend + 1, y, z); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Setting the solution b val = 0.0; - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); ind++; } @@ -876,11 +878,11 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { for (int z = localmesh->zstart; z <= localmesh->zend; z++) { // Setting the initial guess x0 PetscScalar val = x0(localmesh->xend + 1, y, z); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Setting the solution b val = x0(localmesh->xend + 1, y, z); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); ind++; } @@ -889,11 +891,11 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { for (int z = localmesh->zstart; z <= localmesh->zend; z++) { // Setting the initial guess x0 PetscScalar val = x0(localmesh->xend + 1, y, z); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Setting the solution b val = b(localmesh->xend + 1, y, z); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); ind++; } @@ -902,11 +904,11 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { for (int z = localmesh->zstart; z <= localmesh->zend; z++) { // Setting the initial guess x0 PetscScalar val = x0(localmesh->xend + 1, y, z); - VecSetValues(xs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(xs, 1, &ind, &val, INSERT_VALUES)); // Setting the solution b val = 0.5 * (x0(localmesh->xend, y, z) + x0(localmesh->xend + 1, y, z)); - VecSetValues(bs, 1, &ind, &val, INSERT_VALUES); + BOUT_DO_PETSC(VecSetValues(bs, 1, &ind, &val, INSERT_VALUES)); ind++; } @@ -916,25 +918,26 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { ASSERT1(ind == Iend); // Reached end of range // Assemble RHS Vector - VecAssemblyBegin(bs); - VecAssemblyEnd(bs); + BOUT_DO_PETSC(VecAssemblyBegin(bs)); + BOUT_DO_PETSC(VecAssemblyEnd(bs)); // Assemble Trial Solution Vector - VecAssemblyBegin(xs); - VecAssemblyEnd(xs); + BOUT_DO_PETSC(VecAssemblyBegin(xs)); + BOUT_DO_PETSC(VecAssemblyEnd(xs)); ////////////////////////// // Solve the system - KSPSolve(it.ksp, bs, xs); + BOUT_DO_PETSC(KSPSolve(it.ksp, bs, xs)); // Check if the solve converged - KSPConvergedReason reason; - KSPGetConvergedReason(it.ksp, &reason); + BOUT_DO_PETSC(KSPConvergedReason reason); + BOUT_DO_PETSC(KSPGetConvergedReason(it.ksp, &reason)); if (reason <= 0) { - throw BoutException("LaplaceXZ failed to converge. Reason {} ({:d})", - KSPConvergedReasons[reason], static_cast(reason)); + throw BoutException( + "LaplaceXZ failed to converge. Reason {} ({:d})", + BOUT_DO_PETSC(KSPConvergedReasons[reason], static_cast(reason))); } ////////////////////////// @@ -945,7 +948,7 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { if (localmesh->firstX()) { for (int z = localmesh->zstart; z <= localmesh->zend; z++) { PetscScalar val; - VecGetValues(xs, 1, &ind, &val); + BOUT_DO_PETSC(VecGetValues(xs, 1, &ind, &val)); for (int x = localmesh->xstart - 1; x >= 0; --x) { result(x, y, z) = val; } @@ -956,7 +959,7 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { for (int x = localmesh->xstart; x <= localmesh->xend; x++) { for (int z = localmesh->zstart; z <= localmesh->zend; z++) { PetscScalar val; - VecGetValues(xs, 1, &ind, &val); + BOUT_DO_PETSC(VecGetValues(xs, 1, &ind, &val)); result(x, y, z) = val; ind++; } @@ -966,7 +969,7 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { if (localmesh->lastX()) { for (int z = localmesh->zstart; z <= localmesh->zend; z++) { PetscScalar val; - VecGetValues(xs, 1, &ind, &val); + BOUT_DO_PETSC(VecGetValues(xs, 1, &ind, &val)); for (int x = localmesh->xend + 1; x < localmesh->LocalNx; ++x) { result(x, y, z) = val; } From 4e3c8e03fb752340cdb4cd33dad72fd58eb44902 Mon Sep 17 00:00:00 2001 From: David Bold Date: Tue, 3 Nov 2026 13:45:08 +0100 Subject: [PATCH 5/6] Fixup excessive BOUT_DO_PETSC --- .../laplacexy/impls/petsc/laplacexy-petsc.cxx | 12 +++++------- .../laplacexz/impls/petsc/laplacexz-petsc.cxx | 15 ++++++--------- 2 files changed, 11 insertions(+), 16 deletions(-) diff --git a/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx b/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx index f578d7851c..b3554a05c6 100644 --- a/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx +++ b/src/invert/laplacexy/impls/petsc/laplacexy-petsc.cxx @@ -890,10 +890,9 @@ void LaplaceXYpetsc::setCoefs(const Field2D& A, const Field2D& B) { // Set the operator #if PETSC_VERSION_GE(3, 5, 0) - BOUT_DO_PETSC(KSPSetOperators(ksp, BOUT_DO_PETSC(MatA, MatA))); + BOUT_DO_PETSC(KSPSetOperators(ksp, MatA, MatA)); #else - BOUT_DO_PETSC( - KSPSetOperators(ksp, BOUT_DO_PETSC(MatA, MatA, DIFFERENT_NONZERO_PATTERN))); + BOUT_DO_PETSC(KSPSetOperators(ksp, MatA, MatA, DIFFERENT_NONZERO_PATTERN)); #endif // Set coefficients for preconditioner @@ -1206,13 +1205,12 @@ Field2D LaplaceXYpetsc::solve(const Field2D& rhs, const Field2D& x0) { // Solve the system BOUT_DO_PETSC(KSPSolve(ksp, bs, xs)); - BOUT_DO_PETSC(KSPConvergedReason reason); + KSPConvergedReason reason; BOUT_DO_PETSC(KSPGetConvergedReason(ksp, &reason)); if (reason <= 0) { - throw BoutException( - "LaplaceXYpetsc failed to converge. Reason {} ({:d})", - BOUT_DO_PETSC(KSPConvergedReasons[reason], static_cast(reason))); + throw BoutException("LaplaceXYpetsc failed to converge. Reason {} ({:d})", + KSPConvergedReasons[reason], static_cast(reason)); } if (save_performance) { diff --git a/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx b/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx index 0b10802937..a6f8115025 100644 --- a/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx +++ b/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx @@ -725,10 +725,9 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { // Note: This is a hack to force update of the preconditioner matrix #if PETSC_VERSION_GE(3, 5, 0) - BOUT_DO_PETSC(KSPSetOperators(it.ksp, it.BOUT_DO_PETSC(MatA, it.MatP))); + BOUT_DO_PETSC(KSPSetOperators(it.ksp, it.MatA, it.MatP)); #else - BOUT_DO_PETSC( - KSPSetOperators(it.ksp, it.BOUT_DO_PETSC(MatA, it.MatP, SAME_NONZERO_PATTERN))); + BOUT_DO_PETSC(KSPSetOperators(it.ksp, it.MatA, it.MatP, SAME_NONZERO_PATTERN)); #endif } } else { @@ -738,8 +737,7 @@ void LaplaceXZpetsc::setCoefs(const Field3D& Ain, const Field3D& Bin) { #if PETSC_VERSION_GE(3, 5, 0) BOUT_DO_PETSC(KSPSetReusePreconditioner(it.ksp, PETSC_TRUE)); #else - BOUT_DO_PETSC( - KSPSetOperators(it.ksp, it.BOUT_DO_PETSC(MatA, it.MatP, SAME_PRECONDITIONER))); + BOUT_DO_PETSC(KSPSetOperators(it.ksp, it.MatA, it.MatP, SAME_PRECONDITIONER)); #endif } } @@ -931,13 +929,12 @@ Field3D LaplaceXZpetsc::solve(const Field3D& bin, const Field3D& x0in) { BOUT_DO_PETSC(KSPSolve(it.ksp, bs, xs)); // Check if the solve converged - BOUT_DO_PETSC(KSPConvergedReason reason); + KSPConvergedReason reason; BOUT_DO_PETSC(KSPGetConvergedReason(it.ksp, &reason)); if (reason <= 0) { - throw BoutException( - "LaplaceXZ failed to converge. Reason {} ({:d})", - BOUT_DO_PETSC(KSPConvergedReasons[reason], static_cast(reason))); + throw BoutException("LaplaceXZ failed to converge. Reason {} ({:d})", + KSPConvergedReasons[reason], static_cast(reason)); } ////////////////////////// From 5d906b547ca1c31c7adcd5b148bb7e8bc3df4832 Mon Sep 17 00:00:00 2001 From: David Bold Date: Mon, 9 Nov 2026 14:37:50 +0100 Subject: [PATCH 6/6] Resolve clang-tidy review issues --- src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx | 1 + src/solver/impls/imex-bdf2/imex-bdf2.cxx | 1 + 2 files changed, 2 insertions(+) diff --git a/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx b/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx index a6f8115025..93e2004720 100644 --- a/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx +++ b/src/invert/laplacexz/impls/petsc/laplacexz-petsc.cxx @@ -9,6 +9,7 @@ */ #include "bout/build_defines.hxx" +#include "bout/petsclib.hxx" #include "laplacexz-petsc.hxx" diff --git a/src/solver/impls/imex-bdf2/imex-bdf2.cxx b/src/solver/impls/imex-bdf2/imex-bdf2.cxx index e2f76572bb..4bf235feea 100644 --- a/src/solver/impls/imex-bdf2/imex-bdf2.cxx +++ b/src/solver/impls/imex-bdf2/imex-bdf2.cxx @@ -12,6 +12,7 @@ #include #include #include +#include #include #include