Block-Structured AMR Software Framework
 
Loading...
Searching...
No Matches
AMReX_PETSc.H
Go to the documentation of this file.
1#ifndef AMREX_PETSC_H_
2#define AMREX_PETSC_H_
3#include <AMReX_Config.H>
4
5#include <AMReX_Geometry.H>
6#include <AMReX_MultiFab.H>
7#include <AMReX_LayoutData.H>
8#include <AMReX_BndryData.H>
9
10// Do NOT include any other PETSc headers here. They should be included in
11// .cpp files before AMReX_PETSc.H is included.
12#include <petscsystypes.h>
13
14namespace amrex {
15
16class amrex_KSP;
17class amrex_Mat;
18class amrex_Vec;
19
21{
22public:
23
24 PETScABecLap (const BoxArray& grids, const DistributionMapping& dmap,
25 const Geometry& geom_, MPI_Comm comm_);
26
28
29 PETScABecLap (PETScABecLap const&) = delete;
30 PETScABecLap (PETScABecLap &&) noexcept = delete;
31 PETScABecLap& operator= (PETScABecLap const&) = delete;
32 PETScABecLap& operator= (PETScABecLap &&) noexcept = delete;
33
34 void setScalars (Real sa, Real sb);
35 void setACoeffs (const MultiFab& alpha);
36 void setBCoeffs (const Array<const MultiFab*,BL_SPACEDIM>& beta);
37 void setVerbose (int _verbose);
38 void solve (MultiFab& soln, const MultiFab& rhs, Real rel_tol, Real abs_tol,
39 int max_iter, const BndryData& bndry, int max_bndry_order);
40
41#ifdef AMREX_USE_EB
42 void setEBDirichlet (MultiFab const* beb) { m_eb_b_coeffs = beb; }
43#endif
44
45private:
46
47 static constexpr PetscInt regular_stencil_size = 2*AMREX_SPACEDIM + 1;
48 static constexpr PetscInt eb_stencil_size = AMREX_D_TERM(3, *3, *3);
49
52
53 int verbose = 0;
54
58
60
62 BndryData const* m_bndry = nullptr;
63 int m_maxorder = -1;
64
65 std::unique_ptr<amrex_KSP> solver;
66 std::unique_ptr<amrex_Mat> A;
67 std::unique_ptr<amrex_Vec> b;
68 std::unique_ptr<amrex_Vec> x;
69
73
74 MultiFab const* m_eb_b_coeffs = nullptr;
75
76public: // for cuda
77 void prepareSolver ();
78 void loadVectors (MultiFab& soln, const MultiFab& rhs);
79 void getSolution (MultiFab& soln);
80};
81
82[[nodiscard]] std::unique_ptr<PETScABecLap>
83makePetsc (const BoxArray& grids, const DistributionMapping& dmap,
84 const Geometry& geom, MPI_Comm comm_);
85
86}
87
88#endif
#define AMREX_D_TERM(a, b, c)
Definition AMReX_SPACE.H:129
#define BL_SPACEDIM
Definition AMReX_SPACE.H:15
int MPI_Comm
Definition AMReX_ccse-mpi.H:47
static constexpr int MPI_COMM_NULL
Definition AMReX_ccse-mpi.H:55
A BndryData stores and manipulates boundary data information on each side of each box in a BoxArray.
Definition AMReX_BndryData.H:41
A collection of Boxes stored in an Array.
Definition AMReX_BoxArray.H:550
Calculates the distribution of FABs to MPI processes.
Definition AMReX_DistributionMapping.H:41
An Array of FortranArrayBox(FAB)-like Objects.
Definition AMReX_FabArray.H:344
Definition AMReX_FabFactory.H:50
Rectangular problem domain geometry.
Definition AMReX_Geometry.H:73
a one-thingy-per-box distributed object
Definition AMReX_LayoutData.H:13
A collection (stored as an array) of FArrayBox objects.
Definition AMReX_MultiFab.H:38
Definition AMReX_PETSc.H:21
int verbose
Definition AMReX_PETSc.H:53
MultiFab acoefs
Definition AMReX_PETSc.H:55
FabArray< BaseFab< PetscInt > > cell_id
Definition AMReX_PETSc.H:71
MPI_Comm comm
Definition AMReX_PETSc.H:50
void loadVectors(MultiFab &soln, const MultiFab &rhs)
Definition AMReX_PETSc.cpp:608
void getSolution(MultiFab &soln)
Definition AMReX_PETSc.cpp:723
MultiFab diaginv
Definition AMReX_PETSc.H:59
void solve(MultiFab &soln, const MultiFab &rhs, Real rel_tol, Real abs_tol, int max_iter, const BndryData &bndry, int max_bndry_order)
Definition AMReX_PETSc.cpp:139
FabFactory< FArrayBox > const * m_factory
Definition AMReX_PETSc.H:61
static constexpr PetscInt eb_stencil_size
Definition AMReX_PETSc.H:48
std::unique_ptr< amrex_Mat > A
Definition AMReX_PETSc.H:66
void setVerbose(int _verbose)
Definition AMReX_PETSc.cpp:133
BndryData const * m_bndry
Definition AMReX_PETSc.H:62
PETScABecLap(PETScABecLap &&) noexcept=delete
LayoutData< PetscInt > ncells_grid
Definition AMReX_PETSc.H:70
std::unique_ptr< amrex_Vec > x
Definition AMReX_PETSc.H:68
Real scalar_a
Definition AMReX_PETSc.H:57
~PETScABecLap()
Definition AMReX_PETSc.cpp:96
std::unique_ptr< amrex_KSP > solver
Definition AMReX_PETSc.H:65
void setBCoeffs(const Array< const MultiFab *, BL_SPACEDIM > &beta)
Definition AMReX_PETSc.cpp:124
std::unique_ptr< amrex_Vec > b
Definition AMReX_PETSc.H:67
Real scalar_b
Definition AMReX_PETSc.H:57
int m_maxorder
Definition AMReX_PETSc.H:63
MultiFab const * m_eb_b_coeffs
Definition AMReX_PETSc.H:74
Array< MultiFab, AMREX_SPACEDIM > bcoefs
Definition AMReX_PETSc.H:56
void prepareSolver()
Definition AMReX_PETSc.cpp:178
FabArray< BaseFab< PetscInt > > cell_id_vec
Definition AMReX_PETSc.H:72
PETScABecLap(PETScABecLap const &)=delete
void setScalars(Real sa, Real sb)
Definition AMReX_PETSc.cpp:111
void setACoeffs(const MultiFab &alpha)
Definition AMReX_PETSc.cpp:118
static constexpr PetscInt regular_stencil_size
Definition AMReX_PETSc.H:47
Geometry geom
Definition AMReX_PETSc.H:51
Definition AMReX_Amr.cpp:49
std::unique_ptr< PETScABecLap > makePetsc(const BoxArray &grids, const DistributionMapping &dmap, const Geometry &geom, MPI_Comm comm_)
Definition AMReX_PETSc.cpp:54
std::array< T, N > Array
Definition AMReX_Array.H:24