|
Alamo
|
Go to the documentation of this file.
13 #ifndef INTEGRATOR_ALLENCAHN_H // Include guards
14 #define INTEGRATOR_ALLENCAHN_H //
23 #include "AMReX_ParallelDescriptor.H"
24 #include "AMReX_ParmParse.H"
36 #include "Numeric/Stencil.H"
74 std::string type =
"sphere";
76 pp_query_validate(
"alpha.ic.type", type, {
"sphere",
"constant",
"expression",
"bmp",
"random"});
77 if (type ==
"sphere") value.
ic =
new IC::Sphere(value.geom, pp,
"alpha.ic.sphere");
78 else if (type ==
"constant") value.
ic =
new IC::Constant(value.geom, pp,
"alpha.ic.constant");
79 else if (type ==
"expression") value.
ic =
new IC::Expression(value.geom, pp,
"alpha.ic.expression");
80 else if (type ==
"bmp") value.
ic =
new IC::BMP(value.geom, pp,
"alpha.ic.bmp");
81 else if (type ==
"random") value.
ic =
new IC::Random(value.geom, pp,
"alpha.ic.random");
109 const Set::Scalar *DX = this->geom[lev].CellSize();
111 for (amrex::MFIter mfi(*
alpha_mf[lev],
true); mfi.isValid(); ++mfi)
113 const amrex::Box& bx = mfi.tilebox();
114 amrex::Array4<Set::Scalar>
const& alpha_new = (*
alpha_mf[lev]).array(mfi);
115 amrex::Array4<const Set::Scalar>
const& alpha = (*
alpha_old_mf[lev]).array(mfi);
117 amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(
int i,
int j,
int k)
121 Set::Scalar alpha_2 = alpha(i, j, k) * alpha(i, j, k);
123 driving_force += (
ch.chempot /
ch.eps) * (2.0 * alpha(i, j, k) - 6.0 * alpha_2 + 4.0 * alpha_3);
126 driving_force -=
ch.eps *
ch.grad * alpha_lap;
129 alpha(i, j, k) -
ch.L *
dt * driving_force;
140 Set::Scalar dr = sqrt(AMREX_D_TERM(DX[0] * DX[0], +DX[1] * DX[1], +DX[2] * DX[2]));
143 for (amrex::MFIter mfi(*
alpha_mf[lev], amrex::TilingIfNotGPU()); mfi.isValid(); ++mfi)
146 const amrex::Box& bx = mfi.tilebox();
147 amrex::Array4<char>
const& tags = a_tags.array(mfi);
148 amrex::Array4<Set::Scalar>
const& alpha = (*
alpha_mf[lev]).array(mfi);
151 amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(
int i,
int j,
int k)
159 tags(i, j, k) = amrex::TagBox::SET;
void Initialize(const int &a_lev, Set::Field< Set::Scalar > &a_field, Set::Scalar a_time=0.0)
Pure abstract IC object from which all other IC objects inherit.
int number_of_ghost_cells
BC::BC< Set::Scalar > * bc
static void Parse(AllenCahn &value, IO::ParmParse &pp)
void TagCellsForRefinement(int lev, amrex::TagBoxArray &a_tags, Set::Scalar, int)
Eigen::Matrix< amrex::Real, AMREX_SPACEDIM, 1 > Vector
AMREX_FORCE_INLINE Set::Vector Gradient(const amrex::Array4< const Set::Scalar > &f, const int &i, const int &j, const int &k, const int &m, const Set::Scalar dx[AMREX_SPACEDIM], std::array< StencilType, AMREX_SPACEDIM > stencil=DefaultType)
void RegisterNewFab(Set::Field< Set::Scalar > &new_fab, BC::BC< Set::Scalar > *new_bc, int ncomp, int nghost, std::string name, bool writeout)
#define pp_queryclass(...)
AllenCahn(IO::ParmParse &pp)
AMREX_FORCE_INLINE Set::Scalar Laplacian(const amrex::Array4< const Set::Scalar > &f, const int &i, const int &j, const int &k, const int &m, const Set::Scalar dx[AMREX_SPACEDIM])
Set::Scalar refinement_threshold
#define pp_query_default(...)
void Abort(const char *msg)
Set::Field< Set::Scalar > alpha_mf
Collection of numerical integrator objects.
struct Integrator::AllenCahn::@0 ch
void Advance(int lev, Set::Scalar, Set::Scalar dt)
#define pp_query_validate(...)
Set::Field< Set::Scalar > alpha_old_mf
amrex::Vector< amrex::Real > dt
Timesteps for each level of refinement.