Conversation
Read operators from the mesh in CSR format. PETSc handles communication pattern for MPI decomposition.
Handles the mapping between mesh indices and PETSc global indices.
Simple wrapper functions around the Options functionality, so that Array types can be read from the mesh file.
Reads row, column and weights from Mesh, constructs Petsc Mat objects. Working on composition and operation on Field3D objects.
Map from Field3D, including yup/ydown boundaries, to PETSc Vec, perform matrix-vector product, and map back to a Field3D.
Docstrings to explain PetscMapping, PetscOperator and PetscOperators classes. Extend `BOUT_DO_PETSC` macro to print error message, file and line number.
Alternative to initializer_list constructor that avoids
confusion with e.g. `Array{2}` : Is it an array of size 2,
or an array of size 1 containing value 2?
`std::vector` uses `vector(2)` means size 2; `vector{2}`
means value 2.
With this commit, both `Array{2}` and `Array(2)` mean size 2;
`Array<T>::fromValues({2})` means size 1, value 2. It's more verbose
but less prone to typos.
Some care needed with yup/down and yp/ym indices.
`diagonal` creates diagonal matrices from Field3D. `transpose` creates a new matrix that is a transpose of a PetscOperator. This is useful for converting gradients into divergences.
Allows input operators to contain entries in boundary rows. This is important so that the gradient and divergence operators depend on boundary cell values.
![github-actions[bot]](https://avatars.githubusercontent.com/in/15368?s=60&v=4){kind=small}
| void next() final { ++bndry_position; } | ||
| bool isDone() final { return (bndry_position == end(bndry_points)); } | ||
|
|
||
| std::size_t size() { return this->bndry_points.size(); } |
There was a problem hiding this comment.
warning: no header providing "std::size_t" is directly included [misc-include-cleaner]
include/bout/parallel_boundary_region.hxx:5:
- #include <vector>
+ #include <cstddef>
+ #include <vector>
include/bout/petsc_operators.hxx
Outdated
| /// | ||
| /// A pair of PETSc matrices is constructed internally to reorder vectors | ||
| /// between mesh-global and PETSc-global numbering. | ||
| class PetscMapping { |
There was a problem hiding this comment.
warning: class 'PetscMapping' defines a destructor but does not define a copy constructor, a copy assignment operator, a move constructor or a move assignment operator [cppcoreguidelines-special-member-functions]
class PetscMapping {
^| /// The local index is the offset into the packed local PETSc vector. | ||
| template <typename Function> | ||
| void map_local_yup(Function func) const { | ||
| PetscInt row = evolving_region.size() + xin_region.size() + xout_region.size(); |
There was a problem hiding this comment.
warning: narrowing conversion from 'unsigned int' to signed type 'PetscInt' (aka 'int') is implementation-defined [bugprone-narrowing-conversions]
PetscInt row = evolving_region.size() + xin_region.size() + xout_region.size();
^| /// The local index is the offset into the packed local PETSc vector. | ||
| template <typename Function> | ||
| void map_local_ydown(Function func) const { | ||
| PetscInt row = evolving_region.size() + xin_region.size() + xout_region.size() |
There was a problem hiding this comment.
warning: narrowing conversion from 'unsigned int' to signed type 'PetscInt' (aka 'int') is implementation-defined [bugprone-narrowing-conversions]
PetscInt row = evolving_region.size() + xin_region.size() + xout_region.size()
^| /// - applied to a `Field3D`; | ||
| /// - composed with other operators; and | ||
| /// - added or subtracted. | ||
| class PetscOperator { |
There was a problem hiding this comment.
warning: class 'PetscOperator' defines a destructor but does not define a copy constructor, a copy assignment operator, a move constructor or a move assignment operator [cppcoreguidelines-special-member-functions]
class PetscOperator {
^
src/mesh/petsc_operators.cxx
Outdated
| return PetscOperator(this->mapping, mat); | ||
| } | ||
|
|
||
| Field3D PetscOperators::meshGetField3D(Mesh* mesh, const std::string& name) const { |
There was a problem hiding this comment.
warning: method 'meshGetField3D' can be made static [readability-convert-member-functions-to-static]
| Field3D PetscOperators::meshGetField3D(Mesh* mesh, const std::string& name) const { | |
| Field3D PetscOperators::meshGetField3D(Mesh* mesh, const std::string& name) { |
include/bout/petsc_operators.hxx:388:
- Field3D meshGetField3D(Mesh* mesh, const std::string& name) const;
+ static Field3D meshGetField3D(Mesh* mesh, const std::string& name) ;
src/mesh/petsc_operators.cxx
Outdated
| meshGetField3D(mesh, "forward_cell_number"), | ||
| meshGetField3D(mesh, "backward_cell_number"))) { | ||
|
|
||
| int mesh_total_cells; |
There was a problem hiding this comment.
warning: variable 'mesh_total_cells' is not initialized [cppcoreguidelines-init-variables]
| int mesh_total_cells; | |
| int mesh_total_cells = 0; |
src/mesh/petsc_operators.cxx
Outdated
| int mesh_total_cells; | ||
| if (mesh->get(mesh_total_cells, "total_cells") == 0) { | ||
| // Check total number of cells | ||
| if (mesh_total_cells != mapping->size()) { |
There was a problem hiding this comment.
warning: comparison of integers of different signs: 'int' and 'unsigned int' [clang-diagnostic-sign-compare]
if (mesh_total_cells != mapping->size()) {
^| IMEXBDF2* s; | ||
| ierr = PCShellGetContext(pc, reinterpret_cast<void**>(&s)); | ||
| CHKERRQ(ierr); | ||
| PetscCall(PCShellGetContext(pc, reinterpret_cast<void**>(&s))); |
There was a problem hiding this comment.
warning: do not use reinterpret_cast [cppcoreguidelines-pro-type-reinterpret-cast]
PetscCall(PCShellGetContext(pc, reinterpret_cast<void**>(&s)));
^| IMEXBDF2* s; | ||
| ierr = PCShellGetContext(pc, reinterpret_cast<void**>(&s)); | ||
| CHKERRQ(ierr); | ||
| PetscCall(PCShellGetContext(pc, reinterpret_cast<void**>(&s))); |
There was a problem hiding this comment.
warning: multilevel pointer conversion from 'void **' to 'void *', please use explicit cast [bugprone-multi-level-implicit-pointer-conversion]
PetscCall(PCShellGetContext(pc, reinterpret_cast<void**>(&s)));
^Not available unless configured with PETSc.
src/mesh/petsc_operators.cxx
Outdated
| return Region<Ind3D>(xout_indices); | ||
| } | ||
|
|
||
| PetscMapping::PetscMapping(const Field3D& cell_number, const Field3D& forward_cell_number, |
There was a problem hiding this comment.
warning: constructor does not initialize these fields: row_start, row_end, mat_mesh_to_petsc, mat_petsc_to_mesh [cppcoreguidelines-pro-type-member-init]
include/bout/petsc_operators.hxx:236:
- PetscInt row_start, row_end;
+ PetscInt row_start{}, row_end{};include/bout/petsc_operators.hxx:239:
- Mat mat_mesh_to_petsc;
+ Mat mat_mesh_to_petsc{};include/bout/petsc_operators.hxx:242:
- Mat mat_petsc_to_mesh;
+ Mat mat_petsc_to_mesh{};
src/mesh/petsc_operators.cxx
Outdated
| MatDestroy(&this->mat_petsc_to_mesh); | ||
| } | ||
|
|
||
| PetscOperator::PetscOperator(PetscMappingPtr mapping, Array<int> rows, Array<int> cols, |
There was a problem hiding this comment.
warning: constructor does not initialize these fields: mat_operator [cppcoreguidelines-pro-type-member-init]
include/bout/petsc_operators.hxx:335:
- Mat mat_operator;
+ Mat mat_operator{};
src/mesh/petsc_operators.cxx
Outdated
| BOUT_DO_PETSC(MatDestroy(&mat)); | ||
| } | ||
|
|
||
| void PetscOperator::copyToVec(PetscMappingPtr mapping, const Field3D& f, Vec vec) { |
There was a problem hiding this comment.
warning: the parameter 'mapping' is copied for each invocation but only used as a const reference; consider making it a const reference [performance-unnecessary-value-param]
| void PetscOperator::copyToVec(PetscMappingPtr mapping, const Field3D& f, Vec vec) { | |
| void PetscOperator::copyToVec(const PetscMappingPtr& mapping, const Field3D& f, Vec vec) { |
include/bout/petsc_operators.hxx:345:
- static void copyToVec(PetscMappingPtr mapping, const Field3D& f, Vec vec);
+ static void copyToVec(const PetscMappingPtr& mapping, const Field3D& f, Vec vec);
src/mesh/petsc_operators.cxx
Outdated
| return PetscOperator(this->mapping, mat); | ||
| } | ||
|
|
||
| Field3D PetscOperators::meshGetField3D(Mesh* mesh, const std::string& name) const { |
There was a problem hiding this comment.
warning: method 'meshGetField3D' can be made static [readability-convert-member-functions-to-static]
| Field3D PetscOperators::meshGetField3D(Mesh* mesh, const std::string& name) const { | |
| Field3D PetscOperators::meshGetField3D(Mesh* mesh, const std::string& name) { |
include/bout/petsc_operators.hxx:392:
- Field3D meshGetField3D(Mesh* mesh, const std::string& name) const;
+ static Field3D meshGetField3D(Mesh* mesh, const std::string& name) ;Check number of cells against global mapping size, not local size.
Handle non-standard CSR inputs: - Negative indices for missing rows - Missing final entry in rows array
include/bout/petsc_operators.hxx
Outdated
| /// PETSc matrix algebra such as composition, addition, or subtraction. | ||
| PetscOperator(PetscMappingPtr mapping, Mat mat) | ||
| : mapping(std::move(mapping)), mat_operator(mat), | ||
| rhs_vec(createVec(this->mapping->localSize())), |
There was a problem hiding this comment.
warning: narrowing conversion from 'unsigned int' to signed type 'PetscInt' (aka 'int') is implementation-defined [bugprone-narrowing-conversions]
rhs_vec(createVec(this->mapping->localSize())),
^
include/bout/petsc_operators.hxx
Outdated
| PetscOperator(PetscMappingPtr mapping, Mat mat) | ||
| : mapping(std::move(mapping)), mat_operator(mat), | ||
| rhs_vec(createVec(this->mapping->localSize())), | ||
| result_vec(createVec(this->mapping->localSize())) {} |
There was a problem hiding this comment.
warning: narrowing conversion from 'unsigned int' to signed type 'PetscInt' (aka 'int') is implementation-defined [bugprone-narrowing-conversions]
result_vec(createVec(this->mapping->localSize())) {}
^
src/mesh/petsc_operators.cxx
Outdated
| return Region<Ind3D>(xout_indices); | ||
| } | ||
|
|
||
| PetscMapping::PetscMapping(const Field3D& cell_number, const Field3D& forward_cell_number, |
There was a problem hiding this comment.
warning: constructor does not initialize these fields: row_start, row_end, mat_mesh_to_petsc, mat_petsc_to_mesh [cppcoreguidelines-pro-type-member-init]
include/bout/petsc_operators.hxx:245:
- PetscInt row_start, row_end;
+ PetscInt row_start{}, row_end{};include/bout/petsc_operators.hxx:248:
- Mat mat_mesh_to_petsc;
+ Mat mat_mesh_to_petsc{};include/bout/petsc_operators.hxx:251:
- Mat mat_petsc_to_mesh;
+ Mat mat_petsc_to_mesh{};
src/mesh/petsc_operators.cxx
Outdated
| MatDestroy(&this->mat_petsc_to_mesh); | ||
| } | ||
|
|
||
| PetscOperator::PetscOperator(PetscMappingPtr mapping, Array<int> rows, Array<int> cols, |
There was a problem hiding this comment.
warning: constructor does not initialize these fields: mat_operator [cppcoreguidelines-pro-type-member-init]
include/bout/petsc_operators.hxx:344:
- Mat mat_operator;
+ Mat mat_operator{};
src/mesh/petsc_operators.cxx
Outdated
|
|
||
| PetscOperator::PetscOperator(PetscMappingPtr mapping, Array<int> rows, Array<int> cols, | ||
| Array<BoutReal> weights) | ||
| : mapping(std::move(mapping)), rhs_vec(createVec(this->mapping->localSize())), |
There was a problem hiding this comment.
warning: narrowing conversion from 'unsigned int' to signed type 'PetscInt' (aka 'int') is implementation-defined [bugprone-narrowing-conversions]
: mapping(std::move(mapping)), rhs_vec(createVec(this->mapping->localSize())),
^
src/mesh/petsc_operators.cxx
Outdated
| BOUT_DO_PETSC(MatDestroy(&mat)); | ||
| } | ||
|
|
||
| void PetscOperator::copyToVec(PetscMappingPtr mapping, const Field3D& f, Vec vec) { |
There was a problem hiding this comment.
warning: the parameter 'mapping' is copied for each invocation but only used as a const reference; consider making it a const reference [performance-unnecessary-value-param]
| void PetscOperator::copyToVec(PetscMappingPtr mapping, const Field3D& f, Vec vec) { | |
| void PetscOperator::copyToVec(const PetscMappingPtr& mapping, const Field3D& f, Vec vec) { |
include/bout/petsc_operators.hxx:354:
- static void copyToVec(PetscMappingPtr mapping, const Field3D& f, Vec vec);
+ static void copyToVec(const PetscMappingPtr& mapping, const Field3D& f, Vec vec);
src/mesh/petsc_operators.cxx
Outdated
| } | ||
|
|
||
| PetscOperator PetscOperator::diagonal(PetscMappingPtr mapping, const Field3D& f) { | ||
| Vec diag{createVec(mapping->localSize())}; |
There was a problem hiding this comment.
warning: narrowing conversion from 'unsigned int' to signed type 'PetscInt' (aka 'int') is implementation-defined [bugprone-narrowing-conversions]
Vec diag{createVec(mapping->localSize())};
^
src/mesh/petsc_operators.cxx
Outdated
| // -> Create MPIAIJ | ||
| BOUT_DO_PETSC(MatCreate(BoutComm::get(), &mat)); | ||
| BOUT_DO_PETSC(MatSetType(mat, MATMPIAIJ)); | ||
| const int nlocal = mapping->localSize(); |
There was a problem hiding this comment.
warning: narrowing conversion from 'unsigned int' to signed type 'int' is implementation-defined [bugprone-narrowing-conversions]
const int nlocal = mapping->localSize();
^
src/mesh/petsc_operators.cxx
Outdated
| return PetscOperator(this->mapping, mat); | ||
| } | ||
|
|
||
| Field3D PetscOperators::meshGetField3D(Mesh* mesh, const std::string& name) const { |
There was a problem hiding this comment.
warning: method 'meshGetField3D' can be made static [readability-convert-member-functions-to-static]
| Field3D PetscOperators::meshGetField3D(Mesh* mesh, const std::string& name) const { | |
| Field3D PetscOperators::meshGetField3D(Mesh* mesh, const std::string& name) { |
include/bout/petsc_operators.hxx:401:
- Field3D meshGetField3D(Mesh* mesh, const std::string& name) const;
+ static Field3D meshGetField3D(Mesh* mesh, const std::string& name) ;
src/mesh/petsc_operators.cxx
Outdated
| int mesh_total_cells; | ||
| if (mesh->get(mesh_total_cells, "total_cells") == 0) { | ||
| // Check total number of cells | ||
| if (mesh_total_cells != mapping->globalSize()) { |
There was a problem hiding this comment.
warning: comparison of integers of different signs: 'int' and 'unsigned int' [clang-diagnostic-sign-compare]
if (mesh_total_cells != mapping->globalSize()) {
^Changes PetscMapping::size() to PetscMapping::localSize()
Defined in `bout/petsc_interface.hxx`, `bout::petsc::UniqueVec and `UniqueMat` are `std::unique_ptr` wrappers around PETSc Vec and Mat pointers, with custom deleters. PetscOperator copy constructor and assignment deleted, move assignment set to default. This enables some optimizations to operate in-place on temporaries.
|
This looks amazing, thanks @bendudson! I've not finished going through this yet, but how much is this a replacement for |
| struct VecDeleter { | ||
| void operator()(Vec* vec_ptr) const { | ||
| VecDestroy(vec_ptr); | ||
| delete vec_ptr; |
There was a problem hiding this comment.
warning: deleting a pointer through a type that is not marked 'gsl::owner<>'; consider using a smart pointer instead [cppcoreguidelines-owning-memory]
delete vec_ptr;
^Additional context
include/bout/petsc_interface.hxx:67: variable declared here
void operator()(Vec* vec_ptr) const {
^| struct MatDeleter { | ||
| void operator()(Mat* mat_ptr) const { | ||
| MatDestroy(mat_ptr); | ||
| delete mat_ptr; |
There was a problem hiding this comment.
warning: deleting a pointer through a type that is not marked 'gsl::owner<>'; consider using a smart pointer instead [cppcoreguidelines-owning-memory]
delete mat_ptr;
^Additional context
include/bout/petsc_interface.hxx:80: variable declared here
void operator()(Mat* mat_ptr) const {
^
src/mesh/petsc_operators.cxx
Outdated
| return Region<Ind3D>(xout_indices); | ||
| } | ||
|
|
||
| PetscMapping::PetscMapping(const Field3D& cell_number, const Field3D& forward_cell_number, |
There was a problem hiding this comment.
warning: constructor does not initialize these fields: row_start, row_end, mat_mesh_to_petsc, mat_petsc_to_mesh [cppcoreguidelines-pro-type-member-init]
include/bout/petsc_operators.hxx:246:
- PetscInt row_start, row_end;
+ PetscInt row_start{}, row_end{};include/bout/petsc_operators.hxx:249:
- Mat mat_mesh_to_petsc;
+ Mat mat_mesh_to_petsc{};include/bout/petsc_operators.hxx:252:
- Mat mat_petsc_to_mesh;
+ Mat mat_petsc_to_mesh{};
src/mesh/petsc_operators.cxx
Outdated
| PetscOperator::PetscOperator(PetscMappingPtr mapping, Array<int> rows, Array<int> cols, | ||
| Array<BoutReal> weights) | ||
| : mapping(std::move(mapping)), mat_operator(new Mat{nullptr}), | ||
| rhs_vec(createVec(this->mapping->localSize())), |
There was a problem hiding this comment.
warning: narrowing conversion from 'unsigned int' to signed type 'PetscInt' (aka 'int') is implementation-defined [bugprone-narrowing-conversions]
rhs_vec(createVec(this->mapping->localSize())),
^
src/mesh/petsc_operators.cxx
Outdated
| BOUT_DO_PETSC(MatDestroy(&mat)); | ||
| } | ||
|
|
||
| void PetscOperator::copyToVec(PetscMappingPtr mapping, const Field3D& f, Vec vec) { |
There was a problem hiding this comment.
warning: the parameter 'mapping' is copied for each invocation but only used as a const reference; consider making it a const reference [performance-unnecessary-value-param]
| void PetscOperator::copyToVec(PetscMappingPtr mapping, const Field3D& f, Vec vec) { | |
| void PetscOperator::copyToVec(const PetscMappingPtr& mapping, const Field3D& f, Vec vec) { |
include/bout/petsc_operators.hxx:382:
- static void copyToVec(PetscMappingPtr mapping, const Field3D& f, Vec vec);
+ static void copyToVec(const PetscMappingPtr& mapping, const Field3D& f, Vec vec);
src/mesh/petsc_operators.cxx
Outdated
| } | ||
|
|
||
| PetscOperator PetscOperator::diagonal(PetscMappingPtr mapping, const Field3D& f) { | ||
| const UniqueVec diag{createVec(mapping->localSize())}; |
There was a problem hiding this comment.
warning: narrowing conversion from 'unsigned int' to signed type 'PetscInt' (aka 'int') is implementation-defined [bugprone-narrowing-conversions]
const UniqueVec diag{createVec(mapping->localSize())};
^
src/mesh/petsc_operators.cxx
Outdated
| // -> Create MPIAIJ | ||
| BOUT_DO_PETSC(MatCreate(BoutComm::get(), mat.get())); | ||
| BOUT_DO_PETSC(MatSetType(*mat, MATMPIAIJ)); | ||
| const PetscInt nlocal = mapping->localSize(); |
There was a problem hiding this comment.
warning: narrowing conversion from 'unsigned int' to signed type 'PetscInt' (aka 'int') is implementation-defined [bugprone-narrowing-conversions]
const PetscInt nlocal = mapping->localSize();
^Shows how to use the PetscOperators, and compares against the yup/ydown operators. Includes a script to generate the mesh, which currently requires a branch of Zoidberg to run. The code to calculate Grad_par, Div_par and Div_par_Grad_par from forward and backward operators is now in `PetscOperators::getParallel()`. That calculates and returns a collection of `Parallel` operators. This example could be turned into a test: - `Grad_par` should be close between the two methods - Volume integrals of `Div_par` and `Div_par_Grad_par` results should be near zero (flux conserved).
|
Thanks for going through this @ZedThree ! Sorry there's quite a lot, though most is self-contained. The implementation is independent of I've added an example of how these operators are used (https://github.com/boutproject/BOUT-dev/pull/3330/changes#diff-d24773f587abce8257e893cfcfec54002aafaf875fb9cf93f54076b5c60ec7abR16). There is some amount of setup needed to calculate differential operators from the forward / backward maps, that I've put into a helper function that returns a struct of operators (https://github.com/boutproject/BOUT-dev/pull/3330/changes#diff-62e947100b5e999ca6d6d57a6f1281449adc7e2db9e8b7fc5e9273fea887ad88R302). Not sure this is the right API; suggestions welcome! |
I had thought that's essentially what |
The |
|
Ahhh, I see! Thanks! |
This is a function rather than a PetscOperator because it requires at least two matrix-vector products. This implementation requires six Mat-Vec products that could probably be optimized.
Trace multiple field lines per cell, to construct smoother maps. A complication is cells that partially intersect a boundary.
Mesh not available in repo, only generation script.
| dump["div_par_grad_par_op"] = parallel.Div_par_Grad_par(f_neumann); | ||
| dump["div_par_grad_par_yud"] = Div_par_K_Grad_par(1.0, f_neumann); | ||
|
|
||
| bout::writeDefaultOutputFile(dump); |
There was a problem hiding this comment.
warning: no header providing "bout::writeDefaultOutputFile" is directly included [misc-include-cleaner]
examples/fci-operators/fci_operators_example.cxx:7:
- #include "bout/petsc_operators.hxx"
+ #include "bout/options_io.hxx"
+ #include "bout/petsc_operators.hxx"| // Create a parallel gradient operator by combining the parallel | ||
| // length dl = dy * sqrt(g_22) with forward & backward operators | ||
| auto* coords = this->mesh->getCoordinates(); | ||
| Field3D dl = coords->dy * sqrt(coords->g_22); |
There was a problem hiding this comment.
warning: no header providing "sqrt" is directly included [misc-include-cleaner]
src/mesh/petsc_operators.cxx:1:
+ #include <cmath>
src/mesh/petsc_operators.cxx
Outdated
| meshGetField3D(this->mesh, "backward_boundary_fraction"); | ||
|
|
||
| BOUT_FOR(i, dV.getRegion("RGN_NOBNDRY")) { | ||
| dV.yup()[i.yp()] *= forward_boundary_fraction[i]; |
There was a problem hiding this comment.
warning: no viable overloaded operator[] for type 'const Field3D' [clang-diagnostic-error]
dV.yup()[i.yp()] *= forward_boundary_fraction[i];
^Additional context
include/bout/field3d.hxx:360: candidate function not viable: no known conversion from 'SpecificInd<1>' to 'const SpecificInd<IND_TYPE::IND_3D aka 0>' for 1st argument
const BoutReal& operator[](const Ind3D& d) const { return data[d.ind]; }
^include/bout/field3d.hxx:359: candidate function not viable: 'this' argument has type 'const Field3D', but method is not marked const
BoutReal& operator[](const Ind3D& d) { return data[d.ind]; }
^
src/mesh/petsc_operators.cxx
Outdated
|
|
||
| BOUT_FOR(i, dV.getRegion("RGN_NOBNDRY")) { | ||
| dV.yup()[i.yp()] *= forward_boundary_fraction[i]; | ||
| dV.ydown()[i.ym()] *= backward_boundary_fraction[i]; |
There was a problem hiding this comment.
warning: no viable overloaded operator[] for type 'const Field3D' [clang-diagnostic-error]
dV.ydown()[i.ym()] *= backward_boundary_fraction[i];
^Additional context
include/bout/field3d.hxx:360: candidate function not viable: no known conversion from 'SpecificInd<1>' to 'const SpecificInd<IND_TYPE::IND_3D aka 0>' for 1st argument
const BoutReal& operator[](const Ind3D& d) const { return data[d.ind]; }
^include/bout/field3d.hxx:359: candidate function not viable: 'this' argument has type 'const Field3D', but method is not marked const
BoutReal& operator[](const Ind3D& d) { return data[d.ind]; }
^
src/mesh/petsc_operators.cxx
Outdated
| auto Div_par_Grad_par = ((Div_minus * Grad_plus) + (Div_plus * Grad_minus)) * 0.5; | ||
|
|
||
| return Parallel{std::move(Grad_par), std::move(Div_par), | ||
| std::move(Div_par_Grad_par), std::move(dV), |
There was a problem hiding this comment.
warning: passing result of std::move() as a const reference argument; no move will actually happen [performance-move-const-arg]
| std::move(Div_par_Grad_par), std::move(dV), | |
| std::move(Div_par_Grad_par), dV, |
Fix for when metrics are 2D.
More consistent treatment of where values are defined: On cells (C), forward (L+) and backward (L-) legs. Uses typed vectors and operators to ensure that quantities are used consistently. Coauthor: ChatGPT.
Use VecCopy to make a copy of the contents. Rename function from `duplicate` to `copy`.
Still debugging
|
|
||
| dump["grad_par_op"] = parallel.Grad_par(f, f); | ||
| Field3D forward{0.0}; | ||
| BOUT_FOR(i, forward.getRegion("RGN_NOBNDRY")) { forward[i] = f.yup()[i.yp()]; } |
There was a problem hiding this comment.
warning: no header providing "BOUT_FOR" is directly included [misc-include-cleaner]
examples/fci-operators/fci_operators_example.cxx:8:
+ #include "bout/region.hxx"| #include "bout/field3d.hxx" | ||
| #include "bout/mesh.hxx" | ||
| #include "bout/output.hxx" | ||
| #include "bout/output_bout_types.hxx" |
There was a problem hiding this comment.
warning: included header output.hxx is not used directly [misc-include-cleaner]
| #include "bout/output_bout_types.hxx" | |
| #include "bout/output_bout_types.hxx" |
| Mat getPetscToStored() const { return mat_petsc_to_stored; } | ||
|
|
||
| protected: | ||
| PetscLib lib; |
There was a problem hiding this comment.
warning: member variable 'lib' has protected visibility [cppcoreguidelines-non-private-member-variables-in-classes]
PetscLib lib;
^
include/bout/petsc_operators.hxx
Outdated
| void buildPermutation(PetscInt nlocal, PetscInt nglobal, | ||
| const std::vector<int>& stored_indices); | ||
|
|
||
| PetscInt global_size{0}; |
There was a problem hiding this comment.
warning: member variable 'global_size' has protected visibility [cppcoreguidelines-non-private-member-variables-in-classes]
PetscInt global_size{0};
^
include/bout/petsc_operators.hxx
Outdated
| const std::vector<int>& stored_indices); | ||
|
|
||
| PetscInt global_size{0}; | ||
| PetscInt row_start{0}; |
There was a problem hiding this comment.
warning: member variable 'row_start' has protected visibility [cppcoreguidelines-non-private-member-variables-in-classes]
PetscInt row_start{0};
^|
|
||
| /// Wrapper around a PETSc Vec with a compile-time space tag. | ||
| template <typename SpaceTag> | ||
| class PetscSpaceVector { |
There was a problem hiding this comment.
warning: class 'PetscSpaceVector' defines a copy constructor, a copy assignment operator, a move constructor and a move assignment operator but does not define a destructor [cppcoreguidelines-special-member-functions]
class PetscSpaceVector {
^| private: | ||
| void assembleFromCSR(const Array<int>& rows, const Array<int>& cols, | ||
| const Array<BoutReal>& weights) { | ||
| UniqueMat temp{new Mat{nullptr}}; |
There was a problem hiding this comment.
warning: variable 'temp' of type 'UniqueMat' (aka 'unique_ptr<_p_Mat *, bout::petsc::MatDeleter>') can be declared 'const' [misc-const-correctness]
| UniqueMat temp{new Mat{nullptr}}; | |
| UniqueMat const temp{new Mat{nullptr}}; |
| #include "bout/boutexception.hxx" | ||
| #include "bout/field3d.hxx" | ||
| #include "bout/output.hxx" | ||
| #include "bout/output_bout_types.hxx" |
There was a problem hiding this comment.
warning: included header output.hxx is not used directly [misc-include-cleaner]
| #include "bout/output_bout_types.hxx" | |
| #include "bout/output_bout_types.hxx" |
| // Maps PETSc row (or column) indices to the global indices stored in | ||
| // the mesh. This is needed because the PETSc indices depend on the | ||
| // number of processors. | ||
| BOUT_DO_PETSC(MatCreate(BoutComm::get(), &mat_stored_to_petsc)); |
There was a problem hiding this comment.
warning: no header providing "BOUT_DO_PETSC" is directly included [misc-include-cleaner]
src/mesh/petsc_operators.cxx:1:
+ #include "bout/petsclib.hxx"| const Field3D& boundary_leg_number, | ||
| std::shared_ptr<const PetscLegMapping> leg_mapping) const { | ||
| using Op = PetscOperator<LegTag, CellSpaceTag>; | ||
| bout::petsc::UniqueMat mat{new Mat{nullptr}}; |
There was a problem hiding this comment.
warning: no header providing "bout::petsc::UniqueMat" is directly included [misc-include-cleaner]
src/mesh/petsc_operators.cxx:1:
+ #include "bout/petsc_interface.hxx"`applyToField` was ambiguous in meaning apply and convert to Field, and `operator()(Field3D, Field3D)` ended up applying the operator twice. Simplified the application to vectors and Field3D objects.
PetscSpaceVector operators on temporaries can modify in-place.
Coauthor: Claude
| /// @brief BOUT++ PETSc library handle; ensures PETSc is initialised. | ||
| PetscLib lib; | ||
|
|
||
| /// @brief Construct the permutation matrices and populate internal lookup structures. |
There was a problem hiding this comment.
warning: member variable 'lib' has protected visibility [cppcoreguidelines-non-private-member-variables-in-classes]
lib;
^| void buildPermutation(PetscInt nlocal, PetscInt nglobal, | ||
| const std::vector<int>& stored_indices); | ||
|
|
||
| PetscInt global_size{0}; ///< Total global row count. |
There was a problem hiding this comment.
warning: member variable 'global_size' has protected visibility [cppcoreguidelines-non-private-member-variables-in-classes]
es);
^| const std::vector<int>& stored_indices); | ||
|
|
||
| PetscInt global_size{0}; ///< Total global row count. | ||
| PetscInt row_start{0}; ///< First global PETSc row owned locally. |
There was a problem hiding this comment.
warning: member variable 'row_start' has protected visibility [cppcoreguidelines-non-private-member-variables-in-classes]
ount.
^|
|
||
| PetscInt global_size{0}; ///< Total global row count. | ||
| PetscInt row_start{0}; ///< First global PETSc row owned locally. | ||
| PetscInt row_end{0}; ///< One-past-last global PETSc row owned locally. |
There was a problem hiding this comment.
warning: member variable 'row_end' has protected visibility [cppcoreguidelines-non-private-member-variables-in-classes]
ally.
^| PetscInt global_size{0}; ///< Total global row count. | ||
| PetscInt row_start{0}; ///< First global PETSc row owned locally. | ||
| PetscInt row_end{0}; ///< One-past-last global PETSc row owned locally. | ||
| std::vector<int> local_stored_indices; ///< Stored indices in local PETSc row order. |
There was a problem hiding this comment.
warning: member variable 'local_stored_indices' has protected visibility [cppcoreguidelines-non-private-member-variables-in-classes]
ally.
^| std::vector<int> local_stored_indices; ///< Stored indices in local PETSc row order. | ||
| std::unordered_map<int, PetscInt> local_stored_to_petsc; ///< Stored → PETSc lookup. | ||
| Mat mat_stored_to_petsc{nullptr}; ///< Permutation: stored numbering → PETSc ordering. | ||
| Mat mat_petsc_to_stored{nullptr}; ///< Permutation: PETSc ordering → stored numbering. |
There was a problem hiding this comment.
warning: member variable 'mat_petsc_to_stored' has protected visibility [cppcoreguidelines-non-private-member-variables-in-classes]
ordering.
^| /// @param func Callback invoked for each locally stored forward boundary cell. | ||
| template <typename Function> | ||
| void mapLocalYup(Function func) const { | ||
| PetscInt row = evolving_region.size() + xin_region.size() + xout_region.size(); |
There was a problem hiding this comment.
warning: narrowing conversion from 'unsigned int' to signed type 'PetscInt' (aka 'int') is implementation-defined [bugprone-narrowing-conversions]
func) const {
^| /// @param func Callback invoked for each locally stored backward boundary cell. | ||
| template <typename Function> | ||
| void mapLocalYdown(Function func) const { | ||
| PetscInt row = evolving_region.size() + xin_region.size() + xout_region.size() |
There was a problem hiding this comment.
warning: narrowing conversion from 'unsigned int' to signed type 'PetscInt' (aka 'int') is implementation-defined [bugprone-narrowing-conversions]
func) const {
^| /// | ||
| /// @tparam SpaceTag One of CellSpaceTag, ForwardLegSpaceTag, BackwardLegSpaceTag. | ||
| template <typename SpaceTag> | ||
| class PetscSpaceVector { |
There was a problem hiding this comment.
warning: class 'PetscSpaceVector' defines a copy constructor, a copy assignment operator, a move constructor and a move assignment operator but does not define a destructor [cppcoreguidelines-special-member-functions]
ame SpaceTag>
^| /// @param weights CSR nonzero-value array. | ||
| void assembleFromCSR(const Array<int>& rows, const Array<int>& cols, | ||
| const Array<BoutReal>& weights) { | ||
| UniqueMat temp{new Mat{nullptr}}; |
There was a problem hiding this comment.
warning: variable 'temp' of type 'UniqueMat' (aka 'unique_ptr<_p_Mat *, bout::petsc::MatDeleter>') can be declared 'const' [misc-const-correctness]
| UniqueMat temp{new Mat{nullptr}}; | |
| al>& weights) {const |
|
@bendudson "arbitrary parallelisation" is quite a phrase! What exactly will it enable us to do? |
3 participants
Reads matrix operators from the mesh file in CSR format, creating PETSc matrices that apply to Field3Ds. Mainly intended to implement FCI parallel operators as an alternative to yup/down communication.
This enables:
snessolver needs to know which cells are connected, so that PETSc can perform coloring and calculate the system Jacobian efficiently. Knowing the forward and backward matrices will enable this to be done.There are a couple of downsides: