Skip to content

Commit

Permalink
Collect common solver functionality in base classes (#228)
Browse files Browse the repository at this point in the history
* Skeleton for EqnSysBase. Only templates on Nektar eqsys, not particle sys yet.

* Add TimeEvoEqnSysBase class.

* Move equation sys base classes to nektar_interface dir.

* m_field_to_index, m_required_flds, validate_fields, m_int_fld_names.

* Remove virtual inheritance.

* Move load_params to EqnSysBase.

* Add v_InitObject to check required fields, load parameters,  set time-integrated fields.

* Missing utils header.

* Add a base particle system, and a template parameter of TimeEvoEqnSysBase that must inherit from it.

* Add missing v_InitObject() pass-through in TimeEvoEqnSysBase.

* Convert SimpleSOL to use TimeEvoEqnSysBase, PartSysBase.

* Use SolverRunner for SimpleSOL entrypoint.

* Remove debug output.

* Move cell indices map into PartSysBase.

* Move num_parts_tot to PartSysBase and set its value from session file params at construction.

* Correct a param string.

* Correct order of calls in PartSysBase::free().

* Init particle_sys in EqnSysBase.

* Free particle sys memory after fluid solver loop has finished, rather than in destructor.

* Cope with particles being turned off.

* Move solver base classes to nektar_interface subdir.

* Consistent use of this-> in base classes.

* Remove some variables and functions in SimpleSOL/NeutralParticleSystem that duplicate base class members.

* Remove some variables and functions in SOLSystem that duplicate base class members.

* m_ => this->  in SimpleSOL

* Report some parameters commont to all particle systems.

* Add [] operator to NektarFieldIndexMap.

* Add write() functionality to PartSysBase.

* Remove superfluous destructor definition.

* Parameter name tweaks and further docstrings.

* Add zero_array_of_arrays util function to EqnSysBase.

* missing 'virtual'

* Correct v_DoInitialise signature for [email protected].

* Missing 'override'.

* Add an options struct for PartSysBase.

* Call extend_halos_fixed_offset in PartSysBase after constructing the particle-mesh interface.
  • Loading branch information
oparry-ukaea authored May 8, 2024
1 parent db278a4 commit c03be8d
Show file tree
Hide file tree
Showing 13 changed files with 557 additions and 237 deletions.
3 changes: 3 additions & 0 deletions CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -147,6 +147,7 @@ set(HEADER_FILES
${INC_DIR}/nektar_interface/bounding_box_intersection.hpp
${INC_DIR}/nektar_interface/cell_id_translation.hpp
${INC_DIR}/nektar_interface/coordinate_mapping.hpp
${INC_DIR}/nektar_interface/solver_base/eqnsys_base.hpp
${INC_DIR}/nektar_interface/expansion_looping/basis_evaluate_base.hpp
${INC_DIR}/nektar_interface/expansion_looping/expansion_looping.hpp
${INC_DIR}/nektar_interface/expansion_looping/geom_to_expansion_builder.hpp
Expand Down Expand Up @@ -176,6 +177,7 @@ set(HEADER_FILES
${INC_DIR}/nektar_interface/geometry_transport/remote_geom_3d.hpp
${INC_DIR}/nektar_interface/geometry_transport/shape_mapping.hpp
${INC_DIR}/nektar_interface/parameter_store.hpp
${INC_DIR}/nektar_interface/solver_base/partsys_base.hpp
${INC_DIR}/nektar_interface/particle_boundary_conditions.hpp
${INC_DIR}/nektar_interface/particle_cell_mapping/coarse_lookup_map.hpp
${INC_DIR}/nektar_interface/particle_cell_mapping/coarse_mappers_base.hpp
Expand Down Expand Up @@ -209,6 +211,7 @@ set(HEADER_FILES
${INC_DIR}/nektar_interface/particle_interface.hpp
${INC_DIR}/nektar_interface/particle_mesh_interface.hpp
${INC_DIR}/nektar_interface/special_functions.hpp
${INC_DIR}/nektar_interface/solver_base/time_evolved_eqnsys_base.hpp
${INC_DIR}/nektar_interface/utilities.hpp
${INC_DIR}/nektar_interface/utility_mesh_cartesian.hpp
${INC_DIR}/nektar_interface/utility_mesh_plotting.hpp
Expand Down
142 changes: 142 additions & 0 deletions include/nektar_interface/solver_base/eqnsys_base.hpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,142 @@
#ifndef __EQNSYS_BASE_H_
#define __EQNSYS_BASE_H_

#include <SolverUtils/EquationSystem.h>

#include <type_traits>

#include "nektar_interface/solver_base/partsys_base.hpp"
#include "nektar_interface/utilities.hpp"

namespace LU = Nektar::LibUtilities;
namespace SD = Nektar::SpatialDomains;
namespace SU = Nektar::SolverUtils;

namespace NESO::Solvers {

// Fwd declare NEKEQNSYS, PARTSYS as classes
class NEKEQNSYS;
class PARTSYS;

/**
* @brief Base class for Nektar++ equation systems, coupled to a NESO-Particles
* particle system derived from NESO::Solvers::PartSysBase.
*/
template <typename NEKEQNSYS, typename PARTSYS>
class EqnSysBase : public NEKEQNSYS {
// Template param must derive from Nektar's EquationSystem base class
static_assert(std::is_base_of<SU::EquationSystem, NEKEQNSYS>(),
"Template arg to EqnSysBase must derive from "
"Nektar::SolverUtils::EquationSystem");
// Particle system must derive from PartSysBase
static_assert(
std::is_base_of<PartSysBase, PARTSYS>(),
"PARTSYS template arg to EqnSysBase must derive from PartSysBase");

protected:
EqnSysBase(const LU::SessionReaderSharedPtr &session,
const SD::MeshGraphSharedPtr &graph)
: NEKEQNSYS(session, graph), field_to_index(session->GetVariables()),
required_fld_names() {

// If number of particles / number per cell was set in config; construct the
// particle system
int num_parts_per_cell, num_parts_tot;
session->LoadParameter(PartSysBase::NUM_PARTS_TOT_STR, num_parts_tot, -1);
session->LoadParameter(PartSysBase::NUM_PARTS_PER_CELL_STR,
num_parts_per_cell, -1);
this->particles_enabled = num_parts_tot > 0 || num_parts_per_cell > 0;
if (this->particles_enabled) {
this->particle_sys = std::make_shared<PARTSYS>(session, graph);
}
}

/// Field name => index mapper
NESO::NektarFieldIndexMap field_to_index;

/// Particle system
std::shared_ptr<PARTSYS> particle_sys;

/// Flag identifying whether particles were enabled in the config file
bool particles_enabled;

/// List of field names required by the solver
std::vector<std::string> required_fld_names;

/// Placeholder for subclasses to override; called in v_InitObject()
virtual void load_params(){};

/**
* @brief Check that all required fields are defined. All fields must have the
* same number of quad points for now.
*/
void validate_fields() {
int npts_exp = NEKEQNSYS::GetNpoints();
for (auto &fld_name : this->required_fld_names) {
int idx = this->field_to_index.get_idx(fld_name);
// Check field exists

std::string err_msg = "Required field [" + fld_name + "] is not defined.";
NESOASSERT(idx >= 0, err_msg.c_str());

// Check fields all have the same number of quad points
int npts = this->m_fields[idx]->GetNpoints();
err_msg = "Expecting " + std::to_string(npts_exp) +
" quad points, but field '" + fld_name + "' has " +
std::to_string(npts) +
". Check NUMMODES is the same for all required fields.";
NESOASSERT(npts == npts_exp, err_msg.c_str());
}
}

/** @brief Write particle params to stdout on task 0. Ensures they appear just
* after fluid params are written by nektar.
*
* */
virtual void v_DoInitialise() override {
if (this->m_session->GetComm()->TreatAsRankZero() &&
this->particles_enabled) {
particle_sys->add_params_report();
}
NEKEQNSYS::v_DoInitialise();
}

/**
* @brief Free particle system memory after solver loop has finished.
* Prevent further overrides to guarantee that subclasses do the same.
*/
virtual void v_DoSolve() override final {
NEKEQNSYS::v_DoSolve();
if (this->particle_sys) {
this->particle_sys->free();
}
}

/**
* @brief Initialise the equation system, then check required fields are set
* and load parameters.
*/
virtual void v_InitObject(bool create_fields) override {
NEKEQNSYS::v_InitObject(create_fields);

// Ensure that the session file defines all required variables
validate_fields();

// Load parameters
load_params();
}

/**
* @brief Utility function to zero a Nektar++ array of arrays
*
* @param arr Array-of-arrays to zero
*/
void zero_array_of_arrays(Array<OneD, Array<OneD, NekDouble>> &arr) {
for (auto ii = 0; ii < arr.size(); ii++) {
Vmath::Zero(arr[ii].size(), arr[ii], 1);
}
}
};

} // namespace NESO::Solvers
#endif
Loading

0 comments on commit c03be8d

Please sign in to comment.