Serac  0.1
Serac is an implicit thermal strucural mechanics simulation code.
Classes | Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | List of all members
serac::BasePhysics Class Referenceabstract

This is the abstract base class for a generic forward solver. More...

#include <base_physics.hpp>

Collaboration diagram for serac::BasePhysics:
Collaboration graph
[legend]

Classes

struct  ParameterInfo
 The information needed for the physics parameters stored as Finite Element State fields. More...
 

Public Member Functions

 BasePhysics (std::string physics_name, std::string mesh_tag, int cycle=0, double time=0.0, bool checkpoint_to_disk=false)
 Empty constructor. More...
 
 BasePhysics (BasePhysics &&other)=default
 Construct a new Base Physics object (copy constructor) More...
 
virtual double time () const
 Get the current forward-solution time. More...
 
virtual int cycle () const
 Get the current forward-solution cycle iteration number. More...
 
virtual double maxTime () const
 Get the maximum time reached by the forward solver. More...
 
virtual double minTime () const
 Get the initial time used by the forward solver. More...
 
virtual int maxCycle () const
 The maximum cycle (timestep iteration number) reached by the forward solver. More...
 
virtual int minCycle () const
 Get the initial cycle (timestep iteration number) used by the forward solver. More...
 
bool isQuasistatic () const
 Check if the physics is setup as quasistatic. More...
 
virtual std::vector< double > timesteps () const
 Get a vector of the timestep sizes (i.e. \(\Delta t\)s) taken by the forward solver. More...
 
virtual void resetStates (int cycle=0, double time=0.0)=0
 Base method to reset physics states to zero. This does not reset design parameters or shape. More...
 
virtual void completeSetup ()=0
 Complete the setup and allocate the necessary data structures. More...
 
virtual const FiniteElementStatestate (const std::string &state_name) const =0
 Accessor for getting named finite element state primal solution from the physics modules. More...
 
virtual void setState (const std::string &, const FiniteElementState &)=0
 Set the primal solution field values of the underlying physics solver.
 
virtual std::vector< std::string > stateNames () const =0
 Get a vector of the finite element state primal solution names. More...
 
virtual const FiniteElementStateadjoint (const std::string &adjoint_name) const =0
 Accessor for getting named finite element state adjoint solution from the physics modules. More...
 
virtual std::vector< std::string > adjointNames () const
 Get a vector of the finite element state adjoint solution names. More...
 
const FiniteElementStateshapeDisplacement () const
 Accessor for getting the shape displacement field from the physics modules. More...
 
const FiniteElementStateparameter (const std::string &parameter_name) const
 Accessor for getting named finite element state parameter fields from the physics modules. More...
 
const FiniteElementStateparameter (std::size_t parameter_index) const
 Accessor for getting indexed finite element state parameter fields from the physics modules. More...
 
std::vector< std::string > parameterNames ()
 Get a vector of the finite element state parameter names. More...
 
void setParameter (const size_t parameter_index, const FiniteElementState &parameter_state)
 Deep copy a parameter field into the internally-owned parameter used for simulations. More...
 
void setShapeDisplacement (const FiniteElementState &shape_displacement)
 Set the current shape displacement for the underlying mesh. More...
 
virtual const FiniteElementDualcomputeTimestepSensitivity (size_t)
 Compute the implicit sensitivity of the quantity of interest used in defining the adjoint load with respect to the parameter field (d QOI/d state * d state/d parameter). More...
 
virtual const FiniteElementDualcomputeTimestepShapeSensitivity ()
 Compute the implicit sensitivity of the quantity of interest used in defining the adjoint load with respect to the shape displacement field (d QOI/d state * d state/d shape displacement). More...
 
virtual const std::unordered_map< std::string, const serac::FiniteElementDual & > computeInitialConditionSensitivity ()
 Compute the implicit sensitivity of the quantity of interest with respect to the initial condition fields. More...
 
virtual void advanceTimestep (double dt)=0
 Advance the state variables according to the chosen time integrator. More...
 
virtual void setAdjointLoad (std::unordered_map< std::string, const serac::FiniteElementDual & >)
 Set the loads for the adjoint reverse timestep solve.
 
virtual void reverseAdjointTimestep ()
 Solve the adjoint reverse timestep problem. More...
 
virtual void outputStateToDisk (std::optional< std::string > paraview_output_dir={}) const
 Output the current state of the PDE fields in Sidre format and optionally in Paraview format if paraview_output_dir is given. More...
 
FiniteElementState loadCheckpointedState (const std::string &state_name, int cycle) const
 Accessor for getting a single named finite element state primal solution from the physics modules at a given checkpointed cycle index. More...
 
virtual double getCheckpointedTimestep (int cycle) const
 Get a timestep increment which has been previously checkpointed at the give cycle. More...
 
virtual void initializeSummary (axom::sidre::DataStore &datastore, const double t_final, const double dt) const
 Initializes the Sidre structure for simulation summary data. More...
 
virtual void saveSummary (axom::sidre::DataStore &datastore, const double t) const
 Saves the summary data to the Sidre Datastore. More...
 
virtual ~BasePhysics ()=default
 Destroy the Base Solver object.
 
const mfem::ParMesh & mesh () const
 Returns a reference to the mesh object.
 
mfem::ParMesh & mesh ()
 This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.
 

Protected Member Functions

void CreateParaviewDataCollection () const
 Create a paraview data collection for the physics package if requested.
 
void UpdateParaviewDataCollection (const std::string &paraview_output_dir) const
 Update the paraview states, duals, parameters, and metadata (cycle, time) in preparation for output. More...
 
void initializeBasePhysicsStates (int cycle, double time)
 Protected, non-virtual method to reset physics states to zero. This does not reset design parameters or shape. More...
 
virtual std::unordered_map< std::string, FiniteElementStategetCheckpointedStates (int cycle) const
 Accessor for getting all of the primal solutions from the physics modules at a given checkpointed cycle index. More...
 

Protected Attributes

std::string name_ = {}
 Name of the physics module.
 
std::string mesh_tag_ = {}
 ID of the corresponding MFEMSidreDataCollection (denoting a mesh)
 
mfem::ParMesh & mesh_
 The primary mesh.
 
MPI_Comm comm_
 The MPI communicator.
 
std::vector< const serac::FiniteElementState * > states_
 List of finite element primal states associated with this physics module.
 
std::vector< const serac::FiniteElementState * > adjoints_
 List of finite element adjoint states associated with this physics module.
 
std::vector< const serac::FiniteElementDual * > duals_
 List of finite element duals associated with this physics module.
 
std::vector< ParameterInfoparameters_
 A vector of the parameters associated with this physics module.
 
FiniteElementStateshape_displacement_
 The parameter info associated with the shape displacement field. More...
 
std::unique_ptr< FiniteElementDualshape_displacement_sensitivity_
 Sensitivity with respect to the shape displacement field. More...
 
std::unordered_map< std::string, std::vector< serac::FiniteElementState > > checkpoint_states_
 A map containing optionally in-memory checkpointed primal states for transient adjoint solvers.
 
std::unordered_map< std::string, serac::FiniteElementStatecached_checkpoint_states_
 A container relating a checkpointed cycle and the associated finite element state fields. More...
 
std::optional< int > cached_checkpoint_cycle_
 An optional int for disk-based checkpointing containing the cycle number of the last retrieved checkpoint.
 
bool is_quasistatic_ = true
 Whether the simulation is time-independent.
 
double time_
 Current time for the forward pass.
 
double max_time_
 The maximum time reached for the forward solver.
 
double min_time_
 The time the forward solver was initialized to.
 
std::vector< double > timesteps_
 A vector of the timestep sizes (i.e. \(\Delta t\)) taken by the forward solver.
 
int cycle_
 Current cycle (forward pass time iteration count)
 
int max_cycle_
 The maximum cycle (forward pass iteration count) reached by the forward solver.
 
int min_cycle_
 The cycle the forward solver was initialized to.
 
double ode_time_point_
 The value of time at which the ODE solver wants to evaluate the residual.
 
int mpi_rank_
 MPI rank.
 
int mpi_size_
 MPI size.
 
std::unique_ptr< mfem::ParaViewDataCollection > paraview_dc_
 DataCollection pointer for optional paraview output.
 
std::unordered_map< std::string, std::unique_ptr< mfem::ParGridFunction > > paraview_dual_grid_functions_
 A optional map of the dual names and duals in grid function form for paraview output.
 
std::unique_ptr< mfem::ParGridFunction > shape_sensitivity_grid_function_
 A optional view of the shape sensitivity in grid function form for paraview output.
 
BoundaryConditionManager bcs_
 Boundary condition manager instance.
 
bool checkpoint_to_disk_
 A flag denoting whether to save the state to disk or memory as needed for dynamic adjoint solves.
 

Static Protected Attributes

static constexpr int FLOAT_PRECISION_ = 8
 Number of significant figures to output for floating-point.
 

Detailed Description

This is the abstract base class for a generic forward solver.

Definition at line 51 of file base_physics.hpp.

Constructor & Destructor Documentation

◆ BasePhysics() [1/2]

serac::BasePhysics::BasePhysics ( std::string  physics_name,
std::string  mesh_tag,
int  cycle = 0,
double  time = 0.0,
bool  checkpoint_to_disk = false 
)

Empty constructor.

Parameters
[in]physics_nameName of the physics module instance
[in]mesh_tagThe tag for the mesh in the StateManager to construct the physics module on
[in]cycleThe simulation cycle (i.e. timestep iteration) to intialize the physics module to
[in]timeThe simulation time to initialize the physics module to
[in]checkpoint_to_diskA flag to save the transient states on disk instead of memory for the transient adjoint solves

Definition at line 21 of file base_physics.cpp.

◆ BasePhysics() [2/2]

serac::BasePhysics::BasePhysics ( BasePhysics &&  other)
default

Construct a new Base Physics object (copy constructor)

Parameters
otherThe other base physics to copy from

Member Function Documentation

◆ adjoint()

virtual const FiniteElementState& serac::BasePhysics::adjoint ( const std::string &  adjoint_name) const
pure virtual

Accessor for getting named finite element state adjoint solution from the physics modules.

Parameters
adjoint_nameThe name of the Finite Element State adjoint solution to retrieve
Returns
The named adjoint Finite Element State

Implemented in serac::Thermomechanics< order, dim, parameter_space >, serac::SolidMechanics< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >, and serac::HeatTransfer< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >.

◆ adjointNames()

virtual std::vector<std::string> serac::BasePhysics::adjointNames ( ) const
inlinevirtual

Get a vector of the finite element state adjoint solution names.

Returns
The adjoint solution names

Reimplemented in serac::SolidMechanics< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >.

Definition at line 177 of file base_physics.hpp.

◆ advanceTimestep()

virtual void serac::BasePhysics::advanceTimestep ( double  dt)
pure virtual

◆ completeSetup()

virtual void serac::BasePhysics::completeSetup ( )
pure virtual

◆ computeInitialConditionSensitivity()

virtual const std::unordered_map<std::string, const serac::FiniteElementDual&> serac::BasePhysics::computeInitialConditionSensitivity ( )
inlinevirtual

Compute the implicit sensitivity of the quantity of interest with respect to the initial condition fields.

Returns
Fields states corresponding to the sensitivities with respect to the initial condition fields
Precondition
reverseAdjointTimestep with an appropriate adjoint load must be called prior to this method as many times as the forward advance is called.

Reimplemented in serac::SolidMechanics< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >, and serac::HeatTransfer< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >.

Definition at line 303 of file base_physics.hpp.

◆ computeTimestepSensitivity()

virtual const FiniteElementDual& serac::BasePhysics::computeTimestepSensitivity ( size_t  )
inlinevirtual

Compute the implicit sensitivity of the quantity of interest used in defining the adjoint load with respect to the parameter field (d QOI/d state * d state/d parameter).

Returns
The sensitivity of the QOI (given implicitly by the adjoint load) with respect to the parameter
Precondition
completeSetup(), advanceTimestep(), and solveAdjoint() must be called prior to this method.

Reimplemented in serac::SolidMechanics< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >, and serac::HeatTransfer< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >.

Definition at line 275 of file base_physics.hpp.

◆ computeTimestepShapeSensitivity()

virtual const FiniteElementDual& serac::BasePhysics::computeTimestepShapeSensitivity ( )
inlinevirtual

Compute the implicit sensitivity of the quantity of interest used in defining the adjoint load with respect to the shape displacement field (d QOI/d state * d state/d shape displacement).

Returns
The sensitivity with respect to the shape displacement
Precondition
completeSetup(), advanceTimestep(), and solveAdjoint() must be called prior to this method.

Reimplemented in serac::SolidMechanics< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >, and serac::HeatTransfer< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >.

Definition at line 289 of file base_physics.hpp.

◆ cycle()

int serac::BasePhysics::cycle ( ) const
virtual

Get the current forward-solution cycle iteration number.

Returns
The current forward-solution cycle iteration number

Definition at line 50 of file base_physics.cpp.

◆ getCheckpointedStates()

std::unordered_map< std::string, FiniteElementState > serac::BasePhysics::getCheckpointedStates ( int  cycle) const
protectedvirtual

Accessor for getting all of the primal solutions from the physics modules at a given checkpointed cycle index.

Parameters
cycleThe cycle to retrieve state from
Returns
A map containing the primal field names and their associated FiniteElementStates at the requested cycle

Reimplemented in serac::SolidMechanics< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >, and serac::HeatTransfer< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >.

Definition at line 361 of file base_physics.cpp.

◆ getCheckpointedTimestep()

double serac::BasePhysics::getCheckpointedTimestep ( int  cycle) const
virtual

Get a timestep increment which has been previously checkpointed at the give cycle.

Parameters
cycleThe previous 'timestep' number where the timestep increment is requested
Returns
The timestep increment

Definition at line 369 of file base_physics.cpp.

◆ initializeBasePhysicsStates()

void serac::BasePhysics::initializeBasePhysicsStates ( int  cycle,
double  time 
)
protected

Protected, non-virtual method to reset physics states to zero. This does not reset design parameters or shape.

Parameters
[in]cycleThe simulation cycle (i.e. timestep iteration) to intialize the physics module to
[in]timeThe simulation time to initialize the physics module to

Definition at line 62 of file base_physics.cpp.

◆ initializeSummary()

void serac::BasePhysics::initializeSummary ( axom::sidre::DataStore &  datastore,
const double  t_final,
const double  dt 
) const
virtual

Initializes the Sidre structure for simulation summary data.

Parameters
[in]datastoreSidre DataStore where data are saved
[in]t_finalFinal time of the simulation
[in]dtThe time step

Definition at line 219 of file base_physics.cpp.

◆ isQuasistatic()

bool serac::BasePhysics::isQuasistatic ( ) const
inline

Check if the physics is setup as quasistatic.

Returns
true if quasistatic, false if transient

Definition at line 119 of file base_physics.hpp.

◆ loadCheckpointedState()

FiniteElementState serac::BasePhysics::loadCheckpointedState ( const std::string &  state_name,
int  cycle 
) const

Accessor for getting a single named finite element state primal solution from the physics modules at a given checkpointed cycle index.

Parameters
cycleThe cycle to retrieve state from
state_nameThe name of the state to retrieve (e.g. "temperature", "displacement")
Returns
The named primal Finite Element State

Definition at line 336 of file base_physics.cpp.

◆ maxCycle()

int serac::BasePhysics::maxCycle ( ) const
virtual

The maximum cycle (timestep iteration number) reached by the forward solver.

Returns
The maximum cycle reached by the forward solver

Definition at line 54 of file base_physics.cpp.

◆ maxTime()

double serac::BasePhysics::maxTime ( ) const
virtual

Get the maximum time reached by the forward solver.

Returns
The maximum time reached by the forward solver

Definition at line 52 of file base_physics.cpp.

◆ minCycle()

int serac::BasePhysics::minCycle ( ) const
virtual

Get the initial cycle (timestep iteration number) used by the forward solver.

Returns
The initial cycle used by the forward solver

Definition at line 58 of file base_physics.cpp.

◆ minTime()

double serac::BasePhysics::minTime ( ) const
virtual

Get the initial time used by the forward solver.

Returns
The initial time used by the forward solver

Definition at line 56 of file base_physics.cpp.

◆ outputStateToDisk()

void serac::BasePhysics::outputStateToDisk ( std::optional< std::string >  paraview_output_dir = {}) const
virtual

Output the current state of the PDE fields in Sidre format and optionally in Paraview format if paraview_output_dir is given.

Parameters
[in]paraview_output_dirOptional output directory for paraview visualization files

Definition at line 183 of file base_physics.cpp.

◆ parameter() [1/2]

const FiniteElementState& serac::BasePhysics::parameter ( const std::string &  parameter_name) const
inline

Accessor for getting named finite element state parameter fields from the physics modules.

Parameters
parameter_nameThe name of the Finite Element State parameter to retrieve
Returns
The named parameter Finite Element State
Note
The input parameter name should not contain the base physics name. It should be identical to what is in the physics module constructor argument list.

Definition at line 195 of file base_physics.hpp.

◆ parameter() [2/2]

const FiniteElementState& serac::BasePhysics::parameter ( std::size_t  parameter_index) const
inline

Accessor for getting indexed finite element state parameter fields from the physics modules.

Parameters
parameter_indexThe index of the Finite Element State parameter to retrieve
Returns
The indexed parameter Finite Element State

Definition at line 217 of file base_physics.hpp.

◆ parameterNames()

std::vector<std::string> serac::BasePhysics::parameterNames ( )
inline

Get a vector of the finite element state parameter names.

Returns
The parameter names

Definition at line 232 of file base_physics.hpp.

◆ resetStates()

virtual void serac::BasePhysics::resetStates ( int  cycle = 0,
double  time = 0.0 
)
pure virtual

Base method to reset physics states to zero. This does not reset design parameters or shape.

Parameters
[in]cycleThe simulation cycle (i.e. timestep iteration) to intialize the physics module to
[in]timeThe simulation time to initialize the physics module to

Implemented in serac::Thermomechanics< order, dim, parameter_space >, serac::SolidMechanics< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >, and serac::HeatTransfer< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >.

◆ reverseAdjointTimestep()

virtual void serac::BasePhysics::reverseAdjointTimestep ( )
inlinevirtual

Solve the adjoint reverse timestep problem.

Precondition
It is expected that the forward analysis is complete and the current states are valid

Reimplemented in serac::SolidMechanics< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >, and serac::HeatTransfer< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >.

Definition at line 328 of file base_physics.hpp.

◆ saveSummary()

void serac::BasePhysics::saveSummary ( axom::sidre::DataStore &  datastore,
const double  t 
) const
virtual

Saves the summary data to the Sidre Datastore.

Parameters
[in]datastoreSidre DataStore where curves are saved
[in]tThe current time of the simulation

Definition at line 272 of file base_physics.cpp.

◆ setParameter()

void serac::BasePhysics::setParameter ( const size_t  parameter_index,
const FiniteElementState parameter_state 
)

Deep copy a parameter field into the internally-owned parameter used for simulations.

Parameters
parameter_indexthe index of the parameter
parameter_statethe values to use for the specified parameter
Precondition
The discretization space and mesh for this finite element state must be consistent with the arguments provided in the physics module constructor.

The physics module constructs its own parameter FiniteElementState in the physics module constructor. This call sets the internally-owned parameter object by value (i.e. deep copies) from the given argument.

Definition at line 81 of file base_physics.cpp.

◆ setShapeDisplacement()

void serac::BasePhysics::setShapeDisplacement ( const FiniteElementState shape_displacement)

Set the current shape displacement for the underlying mesh.

Parameters
shape_displacementThe shape displacement to copy for use in the physics module

This updates the shape displacement field associated with the underlying mesh. Note that the input FiniteElementState is deep copied into the shape displacement object owned by the StateManager.

Definition at line 101 of file base_physics.cpp.

◆ shapeDisplacement()

const FiniteElementState& serac::BasePhysics::shapeDisplacement ( ) const
inline

Accessor for getting the shape displacement field from the physics modules.

Returns
The shape displacement finite element state

Definition at line 184 of file base_physics.hpp.

◆ state()

virtual const FiniteElementState& serac::BasePhysics::state ( const std::string &  state_name) const
pure virtual

Accessor for getting named finite element state primal solution from the physics modules.

Parameters
state_nameThe name of the Finite Element State primal solution to retrieve
Returns
The named primal Finite Element State

Implemented in serac::Thermomechanics< order, dim, parameter_space >, serac::SolidMechanics< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >, and serac::HeatTransfer< order, dim, Parameters< parameter_space... >, std::integer_sequence< int, parameter_indices... > >.

◆ stateNames()

virtual std::vector<std::string> serac::BasePhysics::stateNames ( ) const
pure virtual

◆ time()

double serac::BasePhysics::time ( ) const
virtual

Get the current forward-solution time.

Returns
The current forward-solution time

Definition at line 48 of file base_physics.cpp.

◆ timesteps()

std::vector< double > serac::BasePhysics::timesteps ( ) const
virtual

Get a vector of the timestep sizes (i.e. \(\Delta t\)s) taken by the forward solver.

Returns
The vector of timestep sizes taken by the foward solver

Definition at line 60 of file base_physics.cpp.

◆ UpdateParaviewDataCollection()

void serac::BasePhysics::UpdateParaviewDataCollection ( const std::string &  paraview_output_dir) const
protected

Update the paraview states, duals, parameters, and metadata (cycle, time) in preparation for output.

Parameters
paraview_output_dirThe directory to write the paraview output

Definition at line 160 of file base_physics.cpp.

Member Data Documentation

◆ cached_checkpoint_states_

std::unordered_map<std::string, serac::FiniteElementState> serac::BasePhysics::cached_checkpoint_states_
mutableprotected

A container relating a checkpointed cycle and the associated finite element state fields.

Note
This is only used when the disk-based checkpointing is used. It avoids thrashing the disk IO by performing a file open/read/close for every separate retrieval of state checkpoint data for the various primal fields.

Definition at line 506 of file base_physics.hpp.

◆ shape_displacement_

FiniteElementState& serac::BasePhysics::shape_displacement_
protected

The parameter info associated with the shape displacement field.

Note
This is owned by the State Manager since it is associated with the mesh

Definition at line 488 of file base_physics.hpp.

◆ shape_displacement_sensitivity_

std::unique_ptr<FiniteElementDual> serac::BasePhysics::shape_displacement_sensitivity_
protected

Sensitivity with respect to the shape displacement field.

Note
This quantity is also called the vector-Jacobian product during back propagation in data science.
This is owned by the physics instance as the sensitivity is with respect to a certain PDE residual (i.e. physics module)

Definition at line 494 of file base_physics.hpp.


The documentation for this class was generated from the following files: