Cardinal
Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | List of all members
NekRSProblem Class Reference

Solve nekRS wrapped as a MOOSE app. More...

#include <NekRSProblem.h>

Inheritance diagram for NekRSProblem:
[legend]

Public Member Functions

 NekRSProblem (const InputParameters &params)
 
 ~NekRSProblem ()
 
int nPoints () const
 
bool isDataTransferHappening (ExternalProblem::Direction direction)
 
Transient * transientExecutioner () const
 
bool isUsrWrkSlotReservedForCoupling (const unsigned int &slot) const
 
void copyIndividualScratchSlot (const unsigned int &slot) const
 
void writeFieldFile (const Real &time, const int &step) const
 
virtual bool isOutputStep () const
 Whether nekRS should write an output file for the current time step. More...
 
virtual void initialSetup () override
 
virtual void externalSolve () override
 
virtual bool converged (unsigned int) override
 
virtual void addExternalVariables () override
 
virtual void syncSolutions (ExternalProblem::Direction direction) override
 
virtual bool nondimensional () const
 
virtual const bool hasMovingNekMesh () const
 
virtual bool synchronizeIn ()
 
virtual bool synchronizeOut ()
 
unsigned int nUsrWrkSlots () const
 
template<typename T >
void volumeSolution (const T &field, double *s)
 
template<typename T >
void boundarySolution (const T &field, double *s)
 
template<typename T >
void writeVolumeSolution (const int elem_id, const T &field, double *s, const std::vector< double > *add=nullptr)
 
template<typename T >
void writeBoundarySolution (const int elem_id, const T &field, double *s)
 
std::string casename () const
 
void mapFaceDataToNekFace (const unsigned int &e, const unsigned int &var_num, const Real &divisor, const Real &additive, double **outgoing_data)
 
void mapVolumeDataToNekVolume (const unsigned int &e, const unsigned int &var_num, const Real &divisor, const Real &additive, double **outgoing_data)
 
void mapFaceDataToNekVolume (const unsigned int &e, const unsigned int &var_num, const Real &divisor, const Real &additive, double **outgoing_data)
 Map nodal points on a MOOSE face element to the GLL points on a Nek volume element. More...
 
template<typename T >
void checkDuplicateEntries (const std::vector< T > &var, const std::string &name) const
 
void checkDuplicateVariableName (const std::string &name) const
 
bool stringHasEnding (std::string const &full, std::string const &ending) const
 

Static Public Member Functions

static InputParameters validParams ()
 

Protected Member Functions

void copyScratchToDevice ()
 Copy the data sent from MOOSE->Nek from host to device. More...
 
void interpolateBoundarySolutionToNek (double *incoming_moose_value, double *outgoing_nek_value)
 
void interpolateVolumeSolutionToNek (const int elem_id, double *incoming_moose_value, double *outgoing_nek_value)
 
void initializeInterpolationMatrices ()
 Initialize interpolation matrices for transfers in/out of nekRS. More...
 
std::string fieldFilePrefix (const int &number) const
 

Protected Attributes

std::unique_ptr< NumericVector< Number > > _serialized_solution
 
const std::string & _casename
 NekRS casename. More...
 
const bool & _write_fld_files
 
const bool & _disable_fld_file_output
 Whether to turn off all field file writing. More...
 
bool _nondimensional
 Whether a dimensionalization action has been added. More...
 
unsigned int _n_usrwrk_slots
 
const unsigned int & _constant_interval
 For constant synchronization intervals, the desired frequency (in units of Nek time steps) More...
 
const bool & _skip_final_field_file
 Whether to skip writing a field file on NekRS's last time steo. More...
 
int _n_surface_elems
 Number of surface elements in the data transfer mesh, across all processes. More...
 
int _n_vertices_per_surface
 Number of vertices per surface element of the transfer mesh. More...
 
int _n_volume_elems
 Number of volume elements in the data transfer mesh, across all processes. More...
 
int _n_vertices_per_volume
 Number of vertices per volume element of the transfer mesh. More...
 
double _start_time
 Start time of the simulation based on NekRS's .par file. More...
 
bool _is_output_step
 Whether the most recent time step was an output file writing step. More...
 
NekRSMesh_nek_mesh
 
NekTimeStepper_timestepper = nullptr
 The time stepper used for selection of time step size. More...
 
Transient * _transient_executioner = nullptr
 Underlying executioner. More...
 
bool _needs_interpolation
 
int _n_points
 Number of points for interpolated fields on the MOOSE mesh. More...
 
const PostprocessorValue * _transfer_in = nullptr
 Postprocessor containing the signal of when a synchronization has occurred. More...
 
double * _interpolation_outgoing = nullptr
 Vandermonde interpolation matrix (for outgoing transfers) More...
 
double * _interpolation_incoming = nullptr
 Vandermonde interpolation matrix (for incoming transfers) More...
 
int _moose_Nq
 For the MOOSE mesh, the number of quadrature points in each coordinate direction. More...
 
const std::vector< unsigned int > * _usrwrk_output = nullptr
 Slots in the nrs->o_usrwrk array to write to a field file. More...
 
const std::vector< std::string > * _usrwrk_output_prefix = nullptr
 Filename prefix to use for naming the field files containing the nrs->o_usrwrk array slots. More...
 
double _elapsedStepSum
 Sum of the elapsed time in NekRS solves. More...
 
double _elapsedTime
 Sum of the total elapsed time in NekRS solves. More...
 
double _tSolveStepMin
 Minimum step solve time. More...
 
double _tSolveStepMax
 Maximum step solve time. More...
 
synchronization::SynchronizationEnum _sync_interval
 When to synchronize the NekRS solution with the mesh mirror. More...
 
std::vector< FieldTransferBase * > _field_transfers
 All of the FieldTransfer objects which pass data in/out of NekRS. More...
 
std::vector< ScalarTransferBase * > _scalar_transfers
 All of the ScalarTransfer objecst which pass data in/out of NekRS. More...
 
std::set< unsigned int > _usrwrk_slots
 Usrwrk slots managed by Cardinal. More...
 

Static Protected Attributes

static bool _first
 flag to indicate whether this is the first pass to serialize the solution More...
 

Detailed Description

Solve nekRS wrapped as a MOOSE app.

This object controls all of the execution of and data transfers to/from nekRS. Any number of data transfers can be added using the [FieldTransfers] block.

Constructor & Destructor Documentation

◆ NekRSProblem()

NekRSProblem::NekRSProblem ( const InputParameters &  params)

◆ ~NekRSProblem()

NekRSProblem::~NekRSProblem ( )

Member Function Documentation

◆ addExternalVariables()

virtual void NekRSProblem::addExternalVariables ( )
overridevirtual

◆ boundarySolution()

template<typename T >
void NekRSProblem::boundarySolution ( const T &  field,
double *  s 
)
inline

Interpolate the NekRS boundary solution onto the boundary MOOSE mesh mirror (re2 -> mirror)

Parameters
[in]ffield to interpolate
[out]sinterpolated boundary value

◆ casename()

std::string NekRSProblem::casename ( ) const
inline

The casename (prefix) of the NekRS files

Returns
casename

◆ checkDuplicateEntries()

template<typename T >
void CardinalProblem::checkDuplicateEntries ( const std::vector< T > &  var,
const std::string &  name 
) const
inlineinherited

Check for duplicate entries in a 1-d vector

Parameters
[in]varinput vector
[in]namestring to use for printing error message

◆ checkDuplicateVariableName()

void CardinalProblem::checkDuplicateVariableName ( const std::string &  name) const
inherited

Check whether the user has already created a variable using one of the protected names that the wrapping is using.

Parameters
[in]namevariable name

◆ converged()

virtual bool NekRSProblem::converged ( unsigned int  )
inlineoverridevirtual

◆ copyIndividualScratchSlot()

void NekRSProblem::copyIndividualScratchSlot ( const unsigned int &  slot) const

Copy an individual slice in the usrwrk array from host to device

Parameters
[in]slotslot in usrwrk array

◆ copyScratchToDevice()

void NekRSProblem::copyScratchToDevice ( )
protected

Copy the data sent from MOOSE->Nek from host to device.

◆ externalSolve()

virtual void NekRSProblem::externalSolve ( )
overridevirtual

◆ fieldFilePrefix()

std::string NekRSProblem::fieldFilePrefix ( const int &  number) const
protected

Get a three-character prefix for use in writing output files for repeated Nek sibling apps.

Parameters
[in]numbermulti-app number

◆ hasMovingNekMesh()

virtual const bool NekRSProblem::hasMovingNekMesh ( ) const
inlinevirtual

Whether the mesh is moving

Returns
whether the mesh is moving

◆ initializeInterpolationMatrices()

void NekRSProblem::initializeInterpolationMatrices ( )
protected

Initialize interpolation matrices for transfers in/out of nekRS.

◆ initialSetup()

virtual void NekRSProblem::initialSetup ( )
overridevirtual

◆ interpolateBoundarySolutionToNek()

void NekRSProblem::interpolateBoundarySolutionToNek ( double *  incoming_moose_value,
double *  outgoing_nek_value 
)
protected

Interpolate the MOOSE mesh mirror solution onto the NekRS boundary mesh (mirror -> re2)

Parameters
[in]incoming_moose_valueMOOSE face values
[out]outgoing_nek_valueinterpolated MOOSE face values onto the NekRS boundary mesh

◆ interpolateVolumeSolutionToNek()

void NekRSProblem::interpolateVolumeSolutionToNek ( const int  elem_id,
double *  incoming_moose_value,
double *  outgoing_nek_value 
)
protected

Interpolate the MOOSE mesh mirror solution onto the NekRS volume mesh (mirror -> re2)

Parameters
[in]elem_idelement ID
[in]incoming_moose_valueMOOSE face values
[out]outgoing_nek_valueinterpolated MOOSE face values onto the NekRS volume mesh

◆ isDataTransferHappening()

bool NekRSProblem::isDataTransferHappening ( ExternalProblem::Direction  direction)

Whether a data transfer to/from Nek is occurring

Parameters
[in]directiondirection of data transfer
Returns
whether a data transfer to Nek is about to occur

◆ isOutputStep()

virtual bool NekRSProblem::isOutputStep ( ) const
virtual

Whether nekRS should write an output file for the current time step.

A nekRS output file (suffix .f000xx) is written if the time step is an integer multiple of the output writing interval or if the time step is the last time step.

Returns
whether to write a nekRS output file

◆ isUsrWrkSlotReservedForCoupling()

bool NekRSProblem::isUsrWrkSlotReservedForCoupling ( const unsigned int &  slot) const

Whether a given slot space is reserved for coupling

Parameters
[in]slotslot in usrwrk array
Returns
whether a usrwrk slot space is reserved by Cardinal

◆ mapFaceDataToNekFace()

void NekRSProblem::mapFaceDataToNekFace ( const unsigned int &  e,
const unsigned int &  var_num,
const Real &  divisor,
const Real &  additive,
double **  outgoing_data 
)

Map nodal points on a MOOSE face element to the GLL points on a Nek face element.

Parameters
[in]eMOOSE element ID
[in]var_numvariable index to fetch MOOSE data from
[in]divisornumber to divide MOOSE data by before sending to Nek (to non-dimensionalize it)
[in]additivenumber to subtract from MOOSE data, before dividing by divisor and sending to Nek (to non-dimensionalize)
[out]outgoing_datadata represented on Nek's GLL points, ready to be applied in Nek

◆ mapFaceDataToNekVolume()

void NekRSProblem::mapFaceDataToNekVolume ( const unsigned int &  e,
const unsigned int &  var_num,
const Real &  divisor,
const Real &  additive,
double **  outgoing_data 
)

Map nodal points on a MOOSE face element to the GLL points on a Nek volume element.

This function is to be used when MOOSE variables are defined over the entire volume (maybe the MOOSE transfer only sent meaningful values to the coupling boundaries), so we need to do a volume interpolation of the incoming MOOSE data into nrs->usrwrk, rather than a face interpolation. This could be optimized in the future to truly only just write the boundary values into the nekRS scratch space rather than the volume values, but it looks right now that our biggest expense occurs in the MOOSE transfer system, not these transfers internally to nekRS.

Parameters
[in]eMOOSE element ID
[in]var_numvariable index to fetch MOOSE data from
[in]divisornumber to divide MOOSE data by before sending to Nek (to non-dimensionalize it)
[in]additivenumber to subtract from MOOSE data, before dividing by divisor and sending to Nek (to non-dimensionalize)
[out]outgoing_datadata represented on Nek's GLL points, ready to be applied in Nek

◆ mapVolumeDataToNekVolume()

void NekRSProblem::mapVolumeDataToNekVolume ( const unsigned int &  e,
const unsigned int &  var_num,
const Real &  divisor,
const Real &  additive,
double **  outgoing_data 
)

Map nodal points on a MOOSE volume element to the GLL points on a Nek volume element.

Parameters
[in]eMOOSE element ID
[in]var_numvariable index to fetch MOOSE data from
[in]divisornumber to divide MOOSE data by before sending to Nek (to non-dimensionalize it)
[in]additivenumber to subtract from MOOSE data, before dividing by divisor and sending to Nek (to non-dimensionalize)
[out]outgoing_datadata represented on Nek's GLL points, ready to be applied in Nek

◆ nondimensional()

virtual bool NekRSProblem::nondimensional ( ) const
inlinevirtual

Whether the solve is in nondimensional form

Returns
whether solve is in nondimensional form

◆ nPoints()

int NekRSProblem::nPoints ( ) const
inline

The number of points (DOFs) on the mesh mirror for this rank; this is used to define the size of various allocated terms

Returns
number of points on this rank

◆ nUsrWrkSlots()

unsigned int NekRSProblem::nUsrWrkSlots ( ) const
inline

Get the number of usrwrk slots allocated

Returns
number of allocated usrwrk slots

◆ stringHasEnding()

bool CardinalProblem::stringHasEnding ( std::string const &  full,
std::string const &  ending 
) const
inherited

Whether a string ends in a particular sub-string

Parameters
[in]fullfull string
[in]endingsub-string ending
Returns
whether full string has ending

◆ synchronizeIn()

virtual bool NekRSProblem::synchronizeIn ( )
virtual

Whether data should be synchronized in to nekRS

Returns
whether inward data synchronization should occur

◆ synchronizeOut()

virtual bool NekRSProblem::synchronizeOut ( )
virtual

Whether data should be synchronized out of nekRS

Returns
whether outward data synchronization should occur

◆ syncSolutions()

virtual void NekRSProblem::syncSolutions ( ExternalProblem::Direction  direction)
overridevirtual

◆ transientExecutioner()

Transient* NekRSProblem::transientExecutioner ( ) const
inline

◆ validParams()

static InputParameters NekRSProblem::validParams ( )
static

◆ volumeSolution()

template<typename T >
void NekRSProblem::volumeSolution ( const T &  field,
double *  s 
)
inline

Interpolate the NekRS volume solution onto the volume MOOSE mesh mirror (re2 -> mirror)

Parameters
[in]ffield to interpolate
[out]sinterpolated volume value

◆ writeBoundarySolution()

template<typename T >
void NekRSProblem::writeBoundarySolution ( const int  elem_id,
const T &  field,
double *  s 
)
inline

Write into the NekRS solution space for coupling boundaries; for setting a mesh position in terms of a displacement, we need to add the displacement to the initial mesh coordinates.

Parameters
[in]elem_idelement ID
[in]fieldfield to write
[in]ssolution values to write for the field for the given element

◆ writeFieldFile()

void NekRSProblem::writeFieldFile ( const Real &  time,
const int &  step 
) const

Write NekRS solution field file

Parameters
[in]timesolution time in NekRS (if NekRS is non-dimensional, this will be non-dimensional)
[in]steptime step index

◆ writeVolumeSolution()

template<typename T >
void NekRSProblem::writeVolumeSolution ( const int  elem_id,
const T &  field,
double *  s,
const std::vector< double > *  add = nullptr 
)
inline

Write into the NekRS solution space for coupling volumes; for setting a mesh position in terms of a displacement, we need to add the displacement to the initial mesh coordinates. For this, the 'add' parameter lets you pass in a vector of values (in NekRS's mesh order, i.e. the re2 order) to add.

Parameters
[in]elem_idelement ID
[in]fieldfield to write
[in]ssolution values to write for the field for the given element
[in]addoptional vector of values to add to each value set on the NekRS end

Member Data Documentation

◆ _casename

const std::string& NekRSProblem::_casename
protected

NekRS casename.

◆ _constant_interval

const unsigned int& NekRSProblem::_constant_interval
protected

For constant synchronization intervals, the desired frequency (in units of Nek time steps)

◆ _disable_fld_file_output

const bool& NekRSProblem::_disable_fld_file_output
protected

Whether to turn off all field file writing.

◆ _elapsedStepSum

double NekRSProblem::_elapsedStepSum
protected

Sum of the elapsed time in NekRS solves.

◆ _elapsedTime

double NekRSProblem::_elapsedTime
protected

Sum of the total elapsed time in NekRS solves.

◆ _field_transfers

std::vector<FieldTransferBase *> NekRSProblem::_field_transfers
protected

All of the FieldTransfer objects which pass data in/out of NekRS.

◆ _first

bool NekRSProblem::_first
staticprotected

flag to indicate whether this is the first pass to serialize the solution

◆ _interpolation_incoming

double* NekRSProblem::_interpolation_incoming = nullptr
protected

Vandermonde interpolation matrix (for incoming transfers)

◆ _interpolation_outgoing

double* NekRSProblem::_interpolation_outgoing = nullptr
protected

Vandermonde interpolation matrix (for outgoing transfers)

◆ _is_output_step

bool NekRSProblem::_is_output_step
protected

Whether the most recent time step was an output file writing step.

◆ _moose_Nq

int NekRSProblem::_moose_Nq
protected

For the MOOSE mesh, the number of quadrature points in each coordinate direction.

◆ _n_points

int NekRSProblem::_n_points
protected

Number of points for interpolated fields on the MOOSE mesh.

◆ _n_surface_elems

int NekRSProblem::_n_surface_elems
protected

Number of surface elements in the data transfer mesh, across all processes.

◆ _n_usrwrk_slots

unsigned int NekRSProblem::_n_usrwrk_slots
protected

Number of slices/slots to allocate in nrs->usrwrk to hold fields for coupling (i.e. data going into NekRS, written by Cardinal), or used for custom user actions, but not for coupling. By default, we just allocate 7 slots (no inherent reason, just a fairly big amount). For memory-limited cases, you can reduce this number to just the bare minimum necessary for your use case.

◆ _n_vertices_per_surface

int NekRSProblem::_n_vertices_per_surface
protected

Number of vertices per surface element of the transfer mesh.

◆ _n_vertices_per_volume

int NekRSProblem::_n_vertices_per_volume
protected

Number of vertices per volume element of the transfer mesh.

◆ _n_volume_elems

int NekRSProblem::_n_volume_elems
protected

Number of volume elements in the data transfer mesh, across all processes.

◆ _needs_interpolation

bool NekRSProblem::_needs_interpolation
protected

Whether an interpolation needs to be performed on the nekRS temperature solution, or if we can just grab the solution at specified points

◆ _nek_mesh

NekRSMesh* NekRSProblem::_nek_mesh
protected

Underlying mesh object on which NekRS exchanges fields with MOOSE or extracts NekRS's solution for I/O features

◆ _nondimensional

bool NekRSProblem::_nondimensional
protected

Whether a dimensionalization action has been added.

◆ _scalar_transfers

std::vector<ScalarTransferBase *> NekRSProblem::_scalar_transfers
protected

All of the ScalarTransfer objecst which pass data in/out of NekRS.

◆ _serialized_solution

std::unique_ptr<NumericVector<Number> > NekRSProblem::_serialized_solution
protected

◆ _skip_final_field_file

const bool& NekRSProblem::_skip_final_field_file
protected

Whether to skip writing a field file on NekRS's last time steo.

◆ _start_time

double NekRSProblem::_start_time
protected

Start time of the simulation based on NekRS's .par file.

◆ _sync_interval

synchronization::SynchronizationEnum NekRSProblem::_sync_interval
protected

When to synchronize the NekRS solution with the mesh mirror.

This parameter determines when to synchronize the NekRS solution with the mesh mirror - this entails:

  • Mapping from the NekRS spectral element mesh to the finite element mesh mirror, to extract information from NekRS and make it available to MOOSE
  • Mapping from the finite element mesh mirror into the NekRS spectral element mesh, to send information from MOOSE into NekRS

Several options are available:

  • 'constant' will simply keep the NekRS solution and the mesh mirror entirely consistent with one another on a given constant frequency of time steps. By default, the 'constant_interval' is 1, so that NekRS and MOOSE communicate with each other on every single time step
  • 'parent_app' will only send data between NekRS and a parent application when (1) the main application has just sent "new" information to NekRS, and when (2) the main application is just about to run a new time step (with updated BCs/source terms from NekRS).

    nekRS is often subcycled relative to the application controlling it - that is, nekRS may be run with a time step 10x smaller than a conduction MOOSE app. If 'interpolate_transfers = false' in the master application, then the data going into nekRS is fixed for each of the subcycled time steps it takes, so these extra data transfers are completely unnecssary. This flag indicates that the information sent from MOOSE to NekRS should only be updated if the data from MOOSE is "new", and likewise whether the NekRS solution should only be interpolated to the mesh mirror once MOOSE is actually "ready" to solve a time step using it.

    NOTE: if 'interpolate_transfers = true' in the master application, then the data coming into nekRS is unique on each subcycled time step, so setting this to true will in effect override interpolate_transfers to be false. For the best performance, you should set interpolate_transfers to false so that you don't even bother computing the interpolated data, since it's not used if this parameter is set to true.

◆ _timestepper

NekTimeStepper* NekRSProblem::_timestepper = nullptr
protected

The time stepper used for selection of time step size.

◆ _transfer_in

const PostprocessorValue* NekRSProblem::_transfer_in = nullptr
protected

Postprocessor containing the signal of when a synchronization has occurred.

◆ _transient_executioner

Transient* NekRSProblem::_transient_executioner = nullptr
protected

Underlying executioner.

◆ _tSolveStepMax

double NekRSProblem::_tSolveStepMax
protected

Maximum step solve time.

◆ _tSolveStepMin

double NekRSProblem::_tSolveStepMin
protected

Minimum step solve time.

◆ _usrwrk_output

const std::vector<unsigned int>* NekRSProblem::_usrwrk_output = nullptr
protected

Slots in the nrs->o_usrwrk array to write to a field file.

◆ _usrwrk_output_prefix

const std::vector<std::string>* NekRSProblem::_usrwrk_output_prefix = nullptr
protected

Filename prefix to use for naming the field files containing the nrs->o_usrwrk array slots.

◆ _usrwrk_slots

std::set<unsigned int> NekRSProblem::_usrwrk_slots
protected

Usrwrk slots managed by Cardinal.

◆ _write_fld_files

const bool& NekRSProblem::_write_fld_files
protected

Whether to disable output file writing by NekRS and replace it by output file writing in Cardinal. Suppose the case name is 'channel'. If this parameter is false, then NekRS will write output files as usual, with names like channel0.f00001 for write step 1, channel0.f00002 for write step 2, and so on. If true, then NekRS itself does not output any files like this, and instead output files are written with names a01channel0.f00001, a01channel0.f00002 (for first Nek app), a02channel0.f00001, a02channel0.f00002 (for second Nek app), and so on. This feature should only be used when running repeated Nek sub apps so that the output from each app is retained. Otherwise, if running N repeated Nek sub apps, only a single output file is obtained because each app overwrites the output files of the other apps in the order that the apps reach the nekrs::outfld function.


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