1 #ifndef AMREX_INTEGRATOR_BASE_H
2 #define AMREX_INTEGRATOR_BASE_H
3 #include <AMReX_Config.H>
8 #if defined(AMREX_PARTICLES)
13 #include <type_traits>
19 #if defined(AMREX_PARTICLES)
21 struct IntegratorOps<T, std::enable_if_t<std::is_base_of_v<amrex::ParticleContainerBase, T> > >
24 static void CreateLike (
amrex::Vector<std::unique_ptr<T> >& V,
const T& Other)
27 V.emplace_back(std::make_unique<T>(Other.Geom(0), Other.ParticleDistributionMap(0), Other.ParticleBoxArray(0)));
28 T& pc = *V[V.size()-1];
35 static void Copy (T& Y,
const T& Other)
38 const bool local =
true;
39 Y.copyParticles(Other, local);
42 static void Saxpy (T& Y,
const amrex::Real a, T& X)
52 auto checkValid = [&]() ->
bool {
53 bool pty_v = pty.isValid();
54 bool ptx_v = ptx.isValid();
56 return pty_v && ptx_v;
59 auto ptIncrement = [&](){ ++pty; ++ptx; };
62 #pragma omp parallel if (Gpu::notInLaunchRegion())
64 for (; checkValid(); ptIncrement())
66 const int npy = pty.numParticles();
67 const int npx = ptx.numParticles();
70 ParticleType* psy = &(pty.GetArrayOfStructs()[0]);
71 ParticleType* psx = &(ptx.GetArrayOfStructs()[0]);
73 auto particle_apply_rhs = T::particle_apply_rhs;
76 ParticleType& py = psy[i];
77 const ParticleType& px = psx[i];
78 particle_apply_rhs(py, a, px);
87 struct IntegratorOps<T, std::enable_if_t<std::is_same_v<amrex::Vector<amrex::MultiFab>, T> > >
93 V.emplace_back(std::make_unique<T>());
94 for (
auto const& other_mf : Other) {
96 V.back()->push_back(
amrex::MultiFab(other_mf.boxArray(), other_mf.DistributionMap(), other_mf.nComp(), nGrow));
103 const int size = Y.size();
104 bool specify_components = !scomp.empty() && ncomp.
size() == scomp.size();
105 for (
int i = 0; i <
size; ++i) {
107 const int iscomp = specify_components ? scomp[i] : 0;
108 const int incomp = specify_components ? ncomp[i] : Other[i].nComp();
118 const int size = Y.size();
119 bool specify_components = !scomp.empty() && ncomp.
size() == scomp.size();
120 for (
int i = 0; i <
size; ++i) {
122 const int iscomp = specify_components ? scomp[i] : 0;
123 const int incomp = specify_components ? ncomp[i] : X[i].nComp();
133 struct IntegratorOps<T, std::enable_if_t<std::is_same_v<amrex::MultiFab, T> > >
140 V.emplace_back(std::make_unique<T>(Other.boxArray(), Other.DistributionMap(), Other.nComp(), nGrow));
143 static void Copy (T& Y,
const T& Other,
const int scomp=0,
const int ncomp=-1,
bool Grow =
true)
147 const int mf_ncomp = ncomp > 0 ? ncomp : Other.nComp();
151 static void Saxpy (T& Y,
const amrex::Real a,
const T& X,
const int scomp=0,
const int ncomp=-1,
bool Grow =
false)
155 const int mf_ncomp = ncomp > 0 ? ncomp : X.nComp();
168 std::function<void(T& rhs, T& state,
const amrex::Real time)>
Rhs;
174 std::function<void(T& rhs, T& state,
const amrex::Real time)>
RhsIm;
180 std::function<void(T& rhs, T& state,
const amrex::Real time)>
RhsEx;
186 std::function<void(T& rhs, T& state,
const amrex::Real time)>
RhsFast;
280 void set_rhs (std::function<
void(T&, T&,
const amrex::Real)> F)
286 std::function<
void(T&, T&,
const amrex::Real)> Fe)
365 virtual amrex::Real
advance (T& S_old, T& S_new, amrex::Real time,
371 virtual void evolve (T& S_out,
const amrex::Real time_out) = 0;
374 amrex::Real timestep_fraction, T& data) = 0;
376 virtual void map_data (std::function<
void(T&)> Map) = 0;
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_ALWAYS_ASSERT(EX)
Definition: AMReX_BLassert.H:50
#define AMREX_GPU_DEVICE
Definition: AMReX_GpuQualifiers.H:18
Definition: AMReX_IntegratorBase.H:163
void set_post_step_action(std::function< void(T &, amrex::Real)> A)
Definition: AMReX_IntegratorBase.H:302
void set_rhs(std::function< void(T &, T &, const amrex::Real)> F)
Definition: AMReX_IntegratorBase.H:280
void set_adaptive_step()
Definition: AMReX_IntegratorBase.H:328
bool use_adaptive_fast_time_step
Flag to enable/disable adaptive time stepping at the fast time scale in multirate methods (bool)
Definition: AMReX_IntegratorBase.H:232
amrex::Real fast_rel_tol
Relative tolerance for adaptive time stepping at the fast time scale (Real)
Definition: AMReX_IntegratorBase.H:264
amrex::Real rel_tol
Relative tolerance for adaptive time stepping (Real)
Definition: AMReX_IntegratorBase.H:253
void set_fast_tolerances(amrex::Real rtol, amrex::Real atol)
Definition: AMReX_IntegratorBase.H:355
void set_post_fast_stage_action(std::function< void(T &, amrex::Real)> A)
Definition: AMReX_IntegratorBase.H:307
std::function< void(T &, amrex::Real)> post_fast_stage_action
The post_stage_action function is called by the integrator on the computed stage just after it is com...
Definition: AMReX_IntegratorBase.H:204
virtual void time_interpolate(const T &S_new, const T &S_old, amrex::Real timestep_fraction, T &data)=0
void set_adaptive_fast_step()
Definition: AMReX_IntegratorBase.H:339
virtual amrex::Real advance(T &S_old, T &S_new, amrex::Real time, amrex::Real dt)=0
Take a single time step from (time, S_old) to (time + dt, S_new) with the given step size.
amrex::Real get_time_step()
Definition: AMReX_IntegratorBase.H:317
amrex::Real fast_abs_tol
Absolute tolerance for adaptive time stepping at the fast time scale (Real)
Definition: AMReX_IntegratorBase.H:270
amrex::Real fast_time_step
Current integrator fast time scale time step size with multirate methods (Real)
Definition: AMReX_IntegratorBase.H:238
std::function< void(T &rhs, T &state, const amrex::Real time)> RhsEx
RhsEx is the explicit right-hand-side function an ImEx integrator will use.
Definition: AMReX_IntegratorBase.H:180
std::function< void(T &, amrex::Real)> post_step_action
The post_step_action function is called by the integrator on the computed state just after it is comp...
Definition: AMReX_IntegratorBase.H:198
void set_post_stage_action(std::function< void(T &, amrex::Real)> A)
Definition: AMReX_IntegratorBase.H:297
std::function< void(T &, amrex::Real)> post_fast_step_action
The post_step_action function is called by the integrator on the computed state just after it is comp...
Definition: AMReX_IntegratorBase.H:210
int max_steps
Max number of internal steps before an error is returned (Long)
Definition: AMReX_IntegratorBase.H:248
virtual void evolve(T &S_out, const amrex::Real time_out)=0
Evolve the current (internal) integrator state to time_out.
std::function< void(T &rhs, T &state, const amrex::Real time)> RhsIm
RhsIm is the implicit right-hand-side function an ImEx integrator will use.
Definition: AMReX_IntegratorBase.H:174
virtual void map_data(std::function< void(T &)> Map)=0
void set_post_fast_step_action(std::function< void(T &, amrex::Real)> A)
Definition: AMReX_IntegratorBase.H:312
void set_imex_rhs(std::function< void(T &, T &, const amrex::Real)> Fi, std::function< void(T &, T &, const amrex::Real)> Fe)
Definition: AMReX_IntegratorBase.H:285
void set_max_steps(int steps)
Definition: AMReX_IntegratorBase.H:344
void set_time_step(amrex::Real dt)
Definition: AMReX_IntegratorBase.H:322
void set_tolerances(amrex::Real rtol, amrex::Real atol)
Definition: AMReX_IntegratorBase.H:349
std::function< void(T &rhs, T &state, const amrex::Real time)> Rhs
Rhs is the right-hand-side function the integrator will use.
Definition: AMReX_IntegratorBase.H:168
virtual ~IntegratorBase()=default
amrex::Long num_steps
Number of integrator time steps (Long)
Definition: AMReX_IntegratorBase.H:243
bool use_adaptive_time_step
Flag to enable/disable adaptive time stepping in single rate methods or at the slow time scale in mul...
Definition: AMReX_IntegratorBase.H:216
void set_fast_rhs(std::function< void(T &, T &, const amrex::Real)> F)
Definition: AMReX_IntegratorBase.H:292
std::function< void(T &, amrex::Real)> post_stage_action
The post_stage_action function is called by the integrator on the computed stage just after it is com...
Definition: AMReX_IntegratorBase.H:192
void set_fast_time_step(amrex::Real dt)
Definition: AMReX_IntegratorBase.H:333
IntegratorBase(const T &)
Definition: AMReX_IntegratorBase.H:276
amrex::Real time_step
Current integrator time step size (Real)
Definition: AMReX_IntegratorBase.H:221
std::function< void(T &rhs, T &state, const amrex::Real time)> RhsFast
RhsFast is the fast timescale right-hand-side function a multirate integrator will use.
Definition: AMReX_IntegratorBase.H:186
amrex::Real previous_time_step
Step size of the last completed step (Real)
Definition: AMReX_IntegratorBase.H:226
amrex::Real abs_tol
Absolute tolerance for adaptive time stepping (Real)
Definition: AMReX_IntegratorBase.H:258
A collection (stored as an array) of FArrayBox objects.
Definition: AMReX_MultiFab.H:38
static void Saxpy(MultiFab &dst, Real a, const MultiFab &src, int srccomp, int dstcomp, int numcomp, int nghost)
dst += a*src
Definition: AMReX_MultiFab.cpp:345
static void Copy(MultiFab &dst, const MultiFab &src, int srccomp, int dstcomp, int numcomp, int nghost)
Copy from src to dst including nghost ghost cells. The two MultiFabs MUST have the same underlying Bo...
Definition: AMReX_MultiFab.cpp:193
Definition: AMReX_ParIter.H:113
This class is a thin wrapper around std::vector. Unlike vector, Vector::operator[] provides bound che...
Definition: AMReX_Vector.H:27
Long size() const noexcept
Definition: AMReX_Vector.H:50
AMREX_GPU_HOST_DEVICE Long size(T const &b) noexcept
integer version
Definition: AMReX_GpuRange.H:26
Definition: AMReX_Amr.cpp:49
std::enable_if_t< std::is_integral_v< T > > ParallelFor(TypeList< CTOs... > ctos, std::array< int, sizeof...(CTOs)> const &runtime_options, T N, F &&f)
Definition: AMReX_CTOParallelForImpl.H:200
void Saxpy(MF &dst, typename MF::value_type a, MF const &src, int scomp, int dcomp, int ncomp, IntVect const &nghost)
dst += a * src
Definition: AMReX_FabArrayUtility.H:1646
void Copy(FabArray< DFAB > &dst, FabArray< SFAB > const &src, int srccomp, int dstcomp, int numcomp, int nghost)
Definition: AMReX_FabArray.H:179
IntVectND< AMREX_SPACEDIM > IntVect
Definition: AMReX_BaseFwd.H:30
static void CreateLike(amrex::Vector< std::unique_ptr< T > > &V, const T &Other, bool Grow=false)
Definition: AMReX_IntegratorBase.H:136
static void Saxpy(T &Y, const amrex::Real a, const T &X, const int scomp=0, const int ncomp=-1, bool Grow=false)
Definition: AMReX_IntegratorBase.H:151
static void Copy(T &Y, const T &Other, const int scomp=0, const int ncomp=-1, bool Grow=true)
Definition: AMReX_IntegratorBase.H:143
static void Saxpy(T &Y, const amrex::Real a, const T &X, const Vector< int > &scomp={}, const Vector< int > &ncomp={}, bool Grow=false)
Definition: AMReX_IntegratorBase.H:115
static void CreateLike(amrex::Vector< std::unique_ptr< T > > &V, const T &Other, bool Grow=false)
Definition: AMReX_IntegratorBase.H:90
static void Copy(T &Y, const T &Other, const Vector< int > &scomp={}, const Vector< int > &ncomp={}, bool Grow=true)
Definition: AMReX_IntegratorBase.H:100
Definition: AMReX_IntegratorBase.H:17
The struct used to store particles.
Definition: AMReX_Particle.H:295