Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

General parser function for external fields #5349

Merged
merged 10 commits into from
Oct 2, 2024
Merged

General parser function for external fields #5349

Show file tree
Hide file tree
Changes from 8 commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
6af1b0d
general parser function for external fields
Sep 30, 2024
ff36f1b
correct the parser dimension and variables
Oct 1, 2024
e216ed4
remove num_arguments and replace string with char
Oct 1, 2024
47a802a
fix bug
Oct 1, 2024
574d19f
remove const qualifier in header
Oct 1, 2024
9981ec7
replace overload function with optional inputs
Oct 1, 2024
4321988
fix indexing error
Oct 2, 2024
585a311
style changes
Oct 2, 2024
25b28e4
fix typo
Oct 2, 2024
2b9436a
set const qualifier
Oct 2, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
8 changes: 1 addition & 7 deletions Source/FieldSolver/FiniteDifferenceSolver/HybridPICModel/HybridPICModel.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -53,13 +53,7 @@ public:
* external current multifab. Note the external current can be a function
* of time and therefore this should be re-evaluated at every step.
*/
void GetCurrentExternal (
ablastr::fields::MultiLevelVectorField const& edge_lengths
);
void GetCurrentExternal (
ablastr::fields::VectorField const& edge_lengths,
int lev
);
void GetCurrentExternal ();

/**
* \brief
Expand Down
169 changes: 17 additions & 152 deletions Source/FieldSolver/FiniteDifferenceSolver/HybridPICModel/HybridPICModel.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -221,167 +221,32 @@ void HybridPICModel::InitData ()
// if the current is time dependent which is what needs to be done to
// write time independent fields on the first step.
for (int lev = 0; lev <= warpx.finestLevel(); ++lev) {
auto edge_lengths = std::array<std::unique_ptr<amrex::MultiFab>, 3>();
#ifdef AMREX_USE_EB
if (EB::enabled()) {
using ablastr::fields::Direction;
auto const & edge_lengths_x = *warpx.m_fields.get(FieldType::edge_lengths, Direction{0}, lev);
auto const & edge_lengths_y = *warpx.m_fields.get(FieldType::edge_lengths, Direction{1}, lev);
auto const & edge_lengths_z = *warpx.m_fields.get(FieldType::edge_lengths, Direction{2}, lev);

edge_lengths = std::array< std::unique_ptr<amrex::MultiFab>, 3 >{
std::make_unique<amrex::MultiFab>(
edge_lengths_x, amrex::make_alias, 0, edge_lengths_x.nComp()),
std::make_unique<amrex::MultiFab>(
edge_lengths_y, amrex::make_alias, 0, edge_lengths_y.nComp()),
std::make_unique<amrex::MultiFab>(
edge_lengths_z, amrex::make_alias, 0, edge_lengths_z.nComp())
};
}
#endif
GetCurrentExternal(ablastr::fields::a2m(edge_lengths), lev);
warpx.ComputeExternalFieldOnGridUsingParser(
FieldType::hybrid_current_fp_external,
m_J_external[0],
m_J_external[1],
m_J_external[2],
lev, PatchType::fine, 'e',
warpx.m_fields.get_alldirs(FieldType::edge_lengths, lev),
warpx.m_fields.get_alldirs(FieldType::face_areas, lev));
}
}

void HybridPICModel::GetCurrentExternal (
ablastr::fields::MultiLevelVectorField const& edge_lengths)
void HybridPICModel::GetCurrentExternal ()
{
if (!m_external_field_has_time_dependence) { return; }

auto& warpx = WarpX::GetInstance();
for (int lev = 0; lev <= warpx.finestLevel(); ++lev)
{
GetCurrentExternal(edge_lengths[lev], lev);
}
}


void HybridPICModel::GetCurrentExternal (
ablastr::fields::VectorField const& edge_lengths,
int lev)
{
// This logic matches closely to WarpX::InitializeExternalFieldsOnGridUsingParser
// except that the parsers include time dependence.
auto & warpx = WarpX::GetInstance();

auto t = warpx.gett_new(lev);

auto dx_lev = warpx.Geom(lev).CellSizeArray();
const RealBox& real_box = warpx.Geom(lev).ProbDomain();

using ablastr::fields::Direction;
amrex::MultiFab * mfx = warpx.m_fields.get(FieldType::hybrid_current_fp_external, Direction{0}, lev);
amrex::MultiFab * mfy = warpx.m_fields.get(FieldType::hybrid_current_fp_external, Direction{1}, lev);
amrex::MultiFab * mfz = warpx.m_fields.get(FieldType::hybrid_current_fp_external, Direction{2}, lev);

const amrex::IntVect x_nodal_flag = mfx->ixType().toIntVect();
const amrex::IntVect y_nodal_flag = mfy->ixType().toIntVect();
const amrex::IntVect z_nodal_flag = mfz->ixType().toIntVect();

// avoid implicit lambda capture
auto Jx_external = m_J_external[0];
auto Jy_external = m_J_external[1];
auto Jz_external = m_J_external[2];

for ( MFIter mfi(*mfx, TilingIfNotGPU()); mfi.isValid(); ++mfi)
{
const amrex::Box& tbx = mfi.tilebox( x_nodal_flag, mfx->nGrowVect() );
const amrex::Box& tby = mfi.tilebox( y_nodal_flag, mfy->nGrowVect() );
const amrex::Box& tbz = mfi.tilebox( z_nodal_flag, mfz->nGrowVect() );

auto const& mfxfab = mfx->array(mfi);
auto const& mfyfab = mfy->array(mfi);
auto const& mfzfab = mfz->array(mfi);

amrex::Array4<amrex::Real> lx, ly, lz;
if (EB::enabled()) {
lx = edge_lengths[0]->array(mfi);
ly = edge_lengths[1]->array(mfi);
lz = edge_lengths[2]->array(mfi);
}

amrex::ParallelFor (tbx, tby, tbz,
[=] AMREX_GPU_DEVICE (int i, int j, int k) {
// skip if node is covered by an embedded boundary
if (lx && lx(i, j, k) <= 0) { return; }

// Shift required in the x-, y-, or z- position
// depending on the index type of the multifab
#if defined(WARPX_DIM_1D_Z)
const amrex::Real x = 0._rt;
const amrex::Real y = 0._rt;
const amrex::Real fac_z = (1._rt - x_nodal_flag[0]) * dx_lev[0] * 0.5_rt;
const amrex::Real z = j*dx_lev[0] + real_box.lo(0) + fac_z;
#elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
const amrex::Real fac_x = (1._rt - x_nodal_flag[0]) * dx_lev[0] * 0.5_rt;
const amrex::Real x = i*dx_lev[0] + real_box.lo(0) + fac_x;
const amrex::Real y = 0._rt;
const amrex::Real fac_z = (1._rt - x_nodal_flag[1]) * dx_lev[1] * 0.5_rt;
const amrex::Real z = j*dx_lev[1] + real_box.lo(1) + fac_z;
#else
const amrex::Real fac_x = (1._rt - x_nodal_flag[0]) * dx_lev[0] * 0.5_rt;
const amrex::Real x = i*dx_lev[0] + real_box.lo(0) + fac_x;
const amrex::Real fac_y = (1._rt - x_nodal_flag[1]) * dx_lev[1] * 0.5_rt;
const amrex::Real y = j*dx_lev[1] + real_box.lo(1) + fac_y;
const amrex::Real fac_z = (1._rt - x_nodal_flag[2]) * dx_lev[2] * 0.5_rt;
const amrex::Real z = k*dx_lev[2] + real_box.lo(2) + fac_z;
#endif
// Initialize the x-component of the field.
mfxfab(i,j,k) = Jx_external(x,y,z,t);
},
[=] AMREX_GPU_DEVICE (int i, int j, int k) {
// skip if node is covered by an embedded boundary
if (ly && ly(i, j, k) <= 0) { return; }

#if defined(WARPX_DIM_1D_Z)
const amrex::Real x = 0._rt;
const amrex::Real y = 0._rt;
const amrex::Real fac_z = (1._rt - y_nodal_flag[0]) * dx_lev[0] * 0.5_rt;
const amrex::Real z = j*dx_lev[0] + real_box.lo(0) + fac_z;
#elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
const amrex::Real fac_x = (1._rt - y_nodal_flag[0]) * dx_lev[0] * 0.5_rt;
const amrex::Real x = i*dx_lev[0] + real_box.lo(0) + fac_x;
const amrex::Real y = 0._rt;
const amrex::Real fac_z = (1._rt - y_nodal_flag[1]) * dx_lev[1] * 0.5_rt;
const amrex::Real z = j*dx_lev[1] + real_box.lo(1) + fac_z;
#elif defined(WARPX_DIM_3D)
const amrex::Real fac_x = (1._rt - y_nodal_flag[0]) * dx_lev[0] * 0.5_rt;
const amrex::Real x = i*dx_lev[0] + real_box.lo(0) + fac_x;
const amrex::Real fac_y = (1._rt - y_nodal_flag[1]) * dx_lev[1] * 0.5_rt;
const amrex::Real y = j*dx_lev[1] + real_box.lo(1) + fac_y;
const amrex::Real fac_z = (1._rt - y_nodal_flag[2]) * dx_lev[2] * 0.5_rt;
const amrex::Real z = k*dx_lev[2] + real_box.lo(2) + fac_z;
#endif
// Initialize the y-component of the field.
mfyfab(i,j,k) = Jy_external(x,y,z,t);
},
[=] AMREX_GPU_DEVICE (int i, int j, int k) {
// skip if node is covered by an embedded boundary
if (lz && lz(i, j, k) <= 0) { return; }

#if defined(WARPX_DIM_1D_Z)
const amrex::Real x = 0._rt;
const amrex::Real y = 0._rt;
const amrex::Real fac_z = (1._rt - z_nodal_flag[0]) * dx_lev[0] * 0.5_rt;
const amrex::Real z = j*dx_lev[0] + real_box.lo(0) + fac_z;
#elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
const amrex::Real fac_x = (1._rt - z_nodal_flag[0]) * dx_lev[0] * 0.5_rt;
const amrex::Real x = i*dx_lev[0] + real_box.lo(0) + fac_x;
const amrex::Real y = 0._rt;
const amrex::Real fac_z = (1._rt - z_nodal_flag[1]) * dx_lev[1] * 0.5_rt;
const amrex::Real z = j*dx_lev[1] + real_box.lo(1) + fac_z;
#elif defined(WARPX_DIM_3D)
const amrex::Real fac_x = (1._rt - z_nodal_flag[0]) * dx_lev[0] * 0.5_rt;
const amrex::Real x = i*dx_lev[0] + real_box.lo(0) + fac_x;
const amrex::Real fac_y = (1._rt - z_nodal_flag[1]) * dx_lev[1] * 0.5_rt;
const amrex::Real y = j*dx_lev[1] + real_box.lo(1) + fac_y;
const amrex::Real fac_z = (1._rt - z_nodal_flag[2]) * dx_lev[2] * 0.5_rt;
const amrex::Real z = k*dx_lev[2] + real_box.lo(2) + fac_z;
#endif
// Initialize the z-component of the field.
mfzfab(i,j,k) = Jz_external(x,y,z,t);
}
);
warpx.ComputeExternalFieldOnGridUsingParser(
FieldType::hybrid_current_fp_external,
m_J_external[0],
m_J_external[1],
m_J_external[2],
lev, PatchType::fine, 'e',
warpx.m_fields.get_alldirs(FieldType::edge_lengths, lev),
warpx.m_fields.get_alldirs(FieldType::face_areas, lev));
}
}

Expand Down
3 changes: 1 addition & 2 deletions Source/FieldSolver/WarpXPushFieldsHybridPIC.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -68,8 +68,7 @@ void WarpX::HybridPICEvolveFields ()
const int sub_steps = m_hybrid_pic_model->m_substeps;

// Get the external current
m_hybrid_pic_model->GetCurrentExternal(
m_fields.get_mr_levels_alldirs(FieldType::edge_lengths, finest_level));
m_hybrid_pic_model->GetCurrentExternal();

// Reference hybrid-PIC multifabs
ablastr::fields::MultiLevelScalarField rho_fp_temp = m_fields.get_mr_levels(FieldType::hybrid_rho_fp_temp, finest_level);
Expand Down
17 changes: 8 additions & 9 deletions Source/Fluids/WarpXFluidContainer.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1010,24 +1010,23 @@ void WarpXFluidContainer::GatherAndPush (
// External field parsers
external_e_fields = (m_E_ext_s == "parse_e_ext_function");
external_b_fields = (m_B_ext_s == "parse_b_ext_function");

amrex::ParserExecutor<4> Exfield_parser;
amrex::ParserExecutor<4> Eyfield_parser;
amrex::ParserExecutor<4> Ezfield_parser;
amrex::ParserExecutor<4> Bxfield_parser;
amrex::ParserExecutor<4> Byfield_parser;
amrex::ParserExecutor<4> Bzfield_parser;

if (external_e_fields){
constexpr int num_arguments = 4; //x,y,z,t
Exfield_parser = m_Ex_parser->compile<num_arguments>();
Eyfield_parser = m_Ey_parser->compile<num_arguments>();
Ezfield_parser = m_Ez_parser->compile<num_arguments>();
Exfield_parser = m_Ex_parser->compile<4>();
Eyfield_parser = m_Ey_parser->compile<4>();
Ezfield_parser = m_Ez_parser->compile<4>();
}

if (external_b_fields){
constexpr int num_arguments = 4; //x,y,z,t
Bxfield_parser = m_Bx_parser->compile<num_arguments>();
Byfield_parser = m_By_parser->compile<num_arguments>();
Bzfield_parser = m_Bz_parser->compile<num_arguments>();
Bxfield_parser = m_Bx_parser->compile<4>();
Byfield_parser = m_By_parser->compile<4>();
Bzfield_parser = m_Bz_parser->compile<4>();
}


Expand Down
12 changes: 6 additions & 6 deletions Source/Initialization/ExternalField.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -127,11 +127,11 @@ ExternalFieldParams::ExternalFieldParams(const amrex::ParmParse& pp_warpx)
str_Bz_ext_grid_function);

Bxfield_parser = std::make_unique<amrex::Parser>(
utils::parser::makeParser(str_Bx_ext_grid_function,{"x","y","z"}));
utils::parser::makeParser(str_Bx_ext_grid_function,{"x","y","z","t"}));
Byfield_parser = std::make_unique<amrex::Parser>(
utils::parser::makeParser(str_By_ext_grid_function,{"x","y","z"}));
utils::parser::makeParser(str_By_ext_grid_function,{"x","y","z","t"}));
Bzfield_parser = std::make_unique<amrex::Parser>(
utils::parser::makeParser(str_Bz_ext_grid_function,{"x","y","z"}));
utils::parser::makeParser(str_Bz_ext_grid_function,{"x","y","z","t"}));
}
//___________________________________________________________________________

Expand Down Expand Up @@ -163,11 +163,11 @@ ExternalFieldParams::ExternalFieldParams(const amrex::ParmParse& pp_warpx)
str_Ez_ext_grid_function);

Exfield_parser = std::make_unique<amrex::Parser>(
utils::parser::makeParser(str_Ex_ext_grid_function,{"x","y","z"}));
utils::parser::makeParser(str_Ex_ext_grid_function,{"x","y","z","t"}));
Eyfield_parser = std::make_unique<amrex::Parser>(
utils::parser::makeParser(str_Ey_ext_grid_function,{"x","y","z"}));
utils::parser::makeParser(str_Ey_ext_grid_function,{"x","y","z","t"}));
Ezfield_parser = std::make_unique<amrex::Parser>(
utils::parser::makeParser(str_Ez_ext_grid_function,{"x","y","z"}));
utils::parser::makeParser(str_Ez_ext_grid_function,{"x","y","z","t"}));
}
//___________________________________________________________________________

Expand Down
Loading
Loading