35 #ifndef IC_EXPRESSION_H_
36 #define IC_EXPRESSION_H_
40 #include "AMReX_Parser.H"
49 std::vector<amrex::ParserExecutor<4>>
f;
60 for (amrex::MFIter mfi(*a_field[lev], amrex::TilingIfNotGPU()); mfi.isValid(); ++mfi)
64 amrex::IndexType type = a_field[lev]->ixType();
65 if (type == amrex::IndexType::TheCellType()) bx = mfi.growntilebox();
66 else if (type == amrex::IndexType::TheNodeType()) bx = mfi.grownnodaltilebox();
69 amrex::Array4<Set::Scalar>
const& field = a_field[lev]->array(mfi);
70 for (
unsigned int n = 0; n <
f.size(); n++)
72 amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(
int i,
int j,
int k)
75 if (
coord == Expression::CoordSys::Cartesian)
77 #if AMREX_SPACEDIM == 1
78 field(i, j, k, n) =
f[n](x(0), 0.0, 0.0, a_time);
79 #elif AMREX_SPACEDIM == 2
80 field(i, j, k, n) =
f[n](x(0), x(1), 0.0, a_time);
81 #elif AMREX_SPACEDIM == 3
82 field(i, j, k, n) =
f[n](x(0), x(1), x(2), a_time);
86 else if (
coord == Expression::CoordSys::Polar)
88 field(i, j, k, n) =
f[n](sqrt(x(0)* x(0) + x(1) * x(1)), std::atan2(x(1), x(0)), x(2), a_time);
94 a_field[lev]->FillBoundary();
101 for (amrex::MFIter mfi(*a_field[lev], amrex::TilingIfNotGPU()); mfi.isValid(); ++mfi)
104 amrex::IndexType type = a_field[lev]->ixType();
105 if (type == amrex::IndexType::TheCellType()) bx = mfi.growntilebox();
106 else if (type == amrex::IndexType::TheNodeType()) bx = mfi.grownnodaltilebox();
110 for (
unsigned int n = 0; n < AMREX_SPACEDIM; n++)
112 amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(
int i,
int j,
int k)
115 if (
coord == Expression::CoordSys::Cartesian)
117 #if AMREX_SPACEDIM == 1
118 field(i, j, k)(n) =
f[n](x(0), 0.0, 0.0, a_time);
119 #elif AMREX_SPACEDIM == 2
120 field(i, j, k)(n) =
f[n](x(0), x(1), 0.0, a_time);
121 #elif AMREX_SPACEDIM == 3
122 field(i, j, k)(n) =
f[n](x(0), x(1), x(2), a_time);
126 else if (
coord == Expression::CoordSys::Polar)
128 field(i, j, k)(n) =
f[n](sqrt(x(0)* x(0) + x(1) * x(1)), std::atan2(x(1), x(0)), x(2), a_time);
134 a_field[lev]->FillBoundary();
139 for (
int i = 0;
true; i++)
141 std::string coordstr =
"cartesian";
145 if (coordstr ==
"cartesian") value.
coord = Expression::CoordSys::Cartesian;
146 else if (coordstr ==
"polar") value.
coord = Expression::CoordSys::Polar;
149 std::string func =
"0.0";
150 std::string name =
"region" + std::to_string(i);
152 if (!pp.
contains(name.data()))
break;
155 value.
parser.push_back(amrex::Parser(func));
161 std::set<std::string> entries = pp.getEntries(prefix +
".constant");
162 std::set<std::string>::iterator entry;
163 for (entry = entries.begin(); entry != entries.end(); entry++)
166 std::string fullname = *entry;
170 value.
parser.back().setConstant(name,val);
173 if (value.
coord == Expression::CoordSys::Cartesian)
175 value.
parser.back().registerVariables({
"x",
"y",
"z",
"t" });
176 value.
f.push_back(value.
parser.back().compile<4>());
178 else if (value.
coord == Expression::CoordSys::Polar)
180 value.
parser.back().registerVariables({
"r",
"theta",
"z",
"t" });
181 value.
f.push_back(value.
parser.back().compile<4>());