OneSweepPerRungeKuttaStepAbstract.template.cpp

Go to the documentation of this file.

// **********************************************************************************************
// ***                                     !!!WARNING!!!                                      ***
// *** WARNING: AUTO GENERATED FILE! DO NOT MODIFY BY HAND! YOUR CHANGES WILL BE OVERWRITTEN! ***
// ***                                     !!!WARNING!!!                                      ***
// ***                  Generated by Peano's Python API: www.peano-framework.org              ***
// **********************************************************************************************
#include "{{CLASSNAME}}.h"
#include "Constants.h"
 
#include <algorithm>
 
#include "tarch/multicore/Lock.h"
#include "tarch/NonCriticalAssertions.h"
 
 
tarch::logging::Log   {{NAMESPACE | join("::")}}::{{CLASSNAME}}::{{CLASSNAME}}::_log( "{{NAMESPACE | join("::")}}::{{CLASSNAME}}::{{CLASSNAME}}" );
 
std::bitset<Dimensions> {{NAMESPACE | join("::")}}::{{CLASSNAME}}::{{CLASSNAME}}::PeriodicBC = {{NAMESPACE | join("::")}}::PeriodicBC;
 
 
{{NAMESPACE | join("::")}}::{{CLASSNAME}}::{{CLASSNAME}}():
  _solverState(SolverState::GridConstruction),
  _minTimeStamp(0.0),
  _maxTimeStamp(0.0),
  _minTimeStampThisTimeStep(0.0),
  _maxTimeStampThisTimeStep(0.0),
  // This way, we still compute the patch size without an overflow
  _minGridCellH(std::numeric_limits<double>::max()/NumberOfGridCellsPerPatchPerAxis/2.0),
  _maxGridCellH(0.0),
  _minGridCellHThisTimeStep(std::numeric_limits<double>::max()/NumberOfGridCellsPerPatchPerAxis/2.0),
  _maxGridCellHThisTimeStep(0.0),
  _minTimeStepSize(std::numeric_limits<double>::max()),
  _maxTimeStepSize(0.0),
  _minTimeStepSizeThisTimeStep(std::numeric_limits<double>::max()),
  _maxTimeStepSizeThisTimeStep(0.0),
  _patchUpdates(0) {
  {{CONSTRUCTOR_IMPLEMENTATION}}
}
 
 
double {{NAMESPACE | join("::")}}::{{CLASSNAME}}::getMaxMeshSize() const {
  return getMaxPatchSize();
}
 
 
double {{NAMESPACE | join("::")}}::{{CLASSNAME}}::getMinMeshSize() const {
  return getMinPatchSize();
}
 
 
double {{NAMESPACE | join("::")}}::{{CLASSNAME}}::getMinPatchSize(bool ofCurrentlyRunningGridSweep) const {
  return getMinGridCellSize(ofCurrentlyRunningGridSweep) * NumberOfGridCellsPerPatchPerAxis;
}
 
 
double {{NAMESPACE | join("::")}}::{{CLASSNAME}}::getMaxPatchSize(bool ofCurrentlyRunningGridSweep) const {
  return getMaxGridCellSize(ofCurrentlyRunningGridSweep) * NumberOfGridCellsPerPatchPerAxis;
}
 
 
double {{NAMESPACE | join("::")}}::{{CLASSNAME}}::getMinGridCellSize(bool ofCurrentlyRunningGridSweep) const {
  return ofCurrentlyRunningGridSweep ? _minGridCellHThisTimeStep : _minGridCellH;
}
 
 
double {{NAMESPACE | join("::")}}::{{CLASSNAME}}::getMaxGridCellSize(bool ofCurrentlyRunningGridSweep) const {
  return ofCurrentlyRunningGridSweep ? _maxGridCellHThisTimeStep : _maxGridCellH;
}
 
 
double {{NAMESPACE | join("::")}}::{{CLASSNAME}}::getMinTimeStamp(bool ofCurrentlyRunningGridSweep) const {
  return ofCurrentlyRunningGridSweep ? _minTimeStampThisTimeStep : _minTimeStamp;
}
 
 
double {{NAMESPACE | join("::")}}::{{CLASSNAME}}::getMaxTimeStamp(bool ofCurrentlyRunningGridSweep) const {
  return ofCurrentlyRunningGridSweep ? _maxTimeStampThisTimeStep : _maxTimeStamp;
}
 
 
double {{NAMESPACE | join("::")}}::{{CLASSNAME}}::getMinTimeStepSize() const {
  return _minTimeStepSize;
}
 
 
double {{NAMESPACE | join("::")}}::{{CLASSNAME}}::getMaxTimeStepSize() const {
  return _maxTimeStepSize;
}
 
 
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::suspendSolversForOneGridSweep() {
  _solverState = SolverState::Suspended;
}
 
 
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::update(double timeStepSize, double timeStamp, double patchSize) {
  tarch::multicore::Lock lock(_semaphore);
 
  if ( tarch::la::greater(timeStepSize,0.0) ) {
    assertion1(timeStepSize<std::numeric_limits<double>::max()/10.0,timeStepSize);
    assertion1(timeStepSize>=0.0,timeStepSize);
    _minTimeStepSizeThisTimeStep = std::min(timeStepSize,_minTimeStepSizeThisTimeStep);
    _maxTimeStepSizeThisTimeStep = std::max(timeStepSize,_maxTimeStepSizeThisTimeStep);
    assertion2(_minTimeStepSizeThisTimeStep<=_maxTimeStepSizeThisTimeStep, _minTimeStepSizeThisTimeStep, _maxTimeStepSizeThisTimeStep );
 
    _patchUpdates++;
 
    _localMinTimeStampThisTimeStep = std::min(_localMinTimeStampThisTimeStep, timeStamp);
    _localMaxTimeStampThisTimeStep = std::max(_localMaxTimeStampThisTimeStep, timeStamp);
  }
 
  assertion1(timeStamp<std::numeric_limits<double>::max()/10.0,timeStamp);
  assertion1(timeStamp>=0.0,timeStamp);
  _maxTimeStampThisTimeStep = std::max(timeStamp,_maxTimeStampThisTimeStep);
  _minTimeStampThisTimeStep = std::min(timeStamp,_minTimeStampThisTimeStep);
  // This assertion makes sense, but if we reset the time stamp size (in couple systems) then
  // it is valid.
  // assertion6( tarch::la::smallerEquals(_minTimeStamp,_maxTimeStampThisTimeStep), _minTimeStamp, _minTimeStampThisTimeStep, _maxTimeStampThisTimeStep, timeStepSize, timeStamp, patchSize );
 
  assertion1(patchSize<std::numeric_limits<double>::max()/NumberOfGridCellsPerPatchPerAxis,patchSize);
  assertion1(patchSize>0.0,patchSize);
  _maxGridCellHThisTimeStep = std::max(_maxGridCellHThisTimeStep,patchSize / NumberOfGridCellsPerPatchPerAxis);
  _minGridCellHThisTimeStep = std::min(_minGridCellHThisTimeStep,patchSize / NumberOfGridCellsPerPatchPerAxis);
 
  logDebug( "update(double,double,double)", "updated with dt=" << timeStepSize << ", t=" << timeStamp << ", h=" << patchSize << ": t(new)=" << _minTimeStampThisTimeStep );
}
 
 
{{NAMESPACE | join("::")}}::{{CLASSNAME}}::SolverState {{NAMESPACE | join("::")}}::{{CLASSNAME}}::getSolverState() const {
  return _solverState;
}
 
 
{% if REFINEMENT_CRITERION_IMPLEMENTATION!="<user-defined>" %}
::exahype2::RefinementCommand {{NAMESPACE | join("::")}}::{{CLASSNAME}}::refinementCriterion(
  const double * __restrict__ Q, // Q[{{NUMBER_OF_UNKNOWNS}}+{{NUMBER_OF_AUXILIARY_VARIABLES}}],
  const tarch::la::Vector<Dimensions,double>&  meshCellCentre,
  const tarch::la::Vector<Dimensions,double>&  meshCellH,
  double                                       t
) {
  {% if REFINEMENT_CRITERION_IMPLEMENTATION=="<empty>" %}
  ::exahype2::RefinementCommand result = ::exahype2::RefinementCommand::Keep;
 
 
  return result;
  {% else %}
  {{REFINEMENT_CRITERION_IMPLEMENTATION}}
  {% endif %}
}
{% endif %}
 
 
{% if INITIAL_CONDITIONS_IMPLEMENTATION!="<user-defined>" %}
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::initialCondition(
  double * __restrict__                        Q,
  const tarch::la::Vector<Dimensions,double>&  meshCellCentre,
  const tarch::la::Vector<Dimensions,double>&  meshCellH,
  bool                                         gridIsConstructed
) {
  {{INITIAL_CONDITIONS_IMPLEMENTATION}}
}
{% endif %}
 
 
{% if BOUNDARY_CONDITIONS_IMPLEMENTATION!="<user-defined>"  and BOUNDARY_CONDITIONS_IMPLEMENTATION!="<none>" %}
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::boundaryConditions(
  const double * __restrict__                  Qinside,   // Qinside[{{NUMBER_OF_UNKNOWNS}}+{{NUMBER_OF_AUXILIARY_VARIABLES}}]
  double * __restrict__                        Qoutside,  // Qoutside[{{NUMBER_OF_UNKNOWNS}}+{{NUMBER_OF_AUXILIARY_VARIABLES}}]
  const tarch::la::Vector<Dimensions,double>&  faceCentre,
  const tarch::la::Vector<Dimensions,double>&  meshCellH,
  double                                       t,
  int                                          normal
) {
  {% if BOUNDARY_CONDITIONS_IMPLEMENTATION!="<empty>" %}
  {{BOUNDARY_CONDITIONS_IMPLEMENTATION}}
  {% endif %}
}
{% endif %}
 
 
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::startGridConstructionStep() {
  assertion( _solverState == SolverState::GridConstruction );
}
 
 
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::finishGridConstructionStep() {
}
 
 
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::startGridInitialisationStep() {
  _solverState = SolverState::GridInitialisation;
}
 
 
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::finishGridInitialisationStep() {
}
 
 
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::startTimeStep(
  double globalMinTimeStamp,
  double globalMaxTimeStamp,
  double globalMinTimeStepSize,
  double globalMaxTimeStepSize
) {
  switch (_solverState) {
    case SolverState::GridConstruction:
      _solverState = SolverState::GridInitialisation;
      break;
  {% for item in range(RK_STEPS-1) %}
    case SolverState::RungeKuttaSubStep{{item}}:
      _solverState = SolverState::RungeKuttaSubStep{{item+1}};
      break;
  {% endfor %}
    case SolverState::Suspended:
    case SolverState::Plotting:
    case SolverState::PlottingAfterGridInitialisation:
    case SolverState::GridInitialisation:
    case SolverState::RungeKuttaSubStep{{RK_STEPS-1}}:
      _solverState = SolverState::RungeKuttaSubStep0;
      break;
    {% if RK_STEPS>1 %}
    case SolverState::RungeKuttaSubStep0AfterGridInitialisation:
      _solverState = SolverState::RungeKuttaSubStep1;
      break;
    {% else %}
    case SolverState::RungeKuttaSubStep0AfterGridInitialisation:
      _solverState = SolverState::RungeKuttaSubStep0;
      break;
    {% endif %}
  }
  
  {{START_TIME_STEP_IMPLEMENTATION}}
 
  if ( isFirstGridSweepOfTimeStep() ) {
    _minTimeStamp = _minTimeStampThisTimeStep;
    _maxTimeStamp = _maxTimeStampThisTimeStep;
 
    _minTimeStampThisTimeStep = std::numeric_limits<double>::max();
    _maxTimeStampThisTimeStep = std::numeric_limits<double>::min();
 
    _minTimeStepSize = _minTimeStepSizeThisTimeStep;
    _maxTimeStepSize = _maxTimeStepSizeThisTimeStep;
 
    _minTimeStepSizeThisTimeStep = std::numeric_limits<double>::max();
    _maxTimeStepSizeThisTimeStep = std::numeric_limits<double>::min();
 
    _minGridCellH = _minGridCellHThisTimeStep;
    _maxGridCellH = _maxGridCellHThisTimeStep;
 
    _maxGridCellHThisTimeStep      = 0.0;
    _minGridCellHThisTimeStep      = std::numeric_limits<double>::max()/NumberOfGridCellsPerPatchPerAxis/2.0;
 
    _localMinTimeStampThisTimeStep = std::numeric_limits<double>::max();
    _localMaxTimeStampThisTimeStep = std::numeric_limits<double>::min();
 
    _patchUpdates    = 0;
  }
}
 
 
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::finishTimeStep() {
  if ( isLastGridSweepOfTimeStep() ) {
 
  #ifdef Parallel
  double newMinTimeStamp = _minTimeStampThisTimeStep;
  double newMaxTimeStamp = _maxTimeStampThisTimeStep;
  double newMinGridCellH   = _minGridCellHThisTimeStep;
  double newMaxGridCellH   = _maxGridCellHThisTimeStep;
  double newMinTimeStepSize = _minTimeStepSizeThisTimeStep;
  double newMaxTimeStepSize = _maxTimeStepSizeThisTimeStep;
 
  double newLocalMinTimeStamp = _localMinTimeStampThisTimeStep;
  double newLocalMaxTimeStamp = _localMaxTimeStampThisTimeStep;
 
  int newPatchUpdates = _patchUpdates;
 
  tarch::mpi::Rank::getInstance().allReduce(
      &newMinTimeStamp,
      &_minTimeStampThisTimeStep,
      1, MPI_DOUBLE,
      MPI_MIN,
      [&]() -> void { tarch::services::ServiceRepository::getInstance().receiveDanglingMessages(); }
      );
  tarch::mpi::Rank::getInstance().allReduce(
      &newMaxTimeStamp,
      &_maxTimeStampThisTimeStep,
      1, MPI_DOUBLE,
      MPI_MAX,
      [&]() -> void { tarch::services::ServiceRepository::getInstance().receiveDanglingMessages(); }
      );
 
  tarch::mpi::Rank::getInstance().allReduce(
      &newMinGridCellH,
      &_minGridCellHThisTimeStep,
      1, MPI_DOUBLE,
      MPI_MIN,
      [&]() -> void { tarch::services::ServiceRepository::getInstance().receiveDanglingMessages(); }
      );
  tarch::mpi::Rank::getInstance().allReduce(
      &newMaxGridCellH,
      &_maxGridCellHThisTimeStep,
      1, MPI_DOUBLE,
      MPI_MAX,
      [&]() -> void { tarch::services::ServiceRepository::getInstance().receiveDanglingMessages(); }
      );
 
  tarch::mpi::Rank::getInstance().allReduce(
      &newMinTimeStepSize,
      &_minTimeStepSizeThisTimeStep,
      1, MPI_DOUBLE,
      MPI_MIN,
      [&]() -> void { tarch::services::ServiceRepository::getInstance().receiveDanglingMessages(); }
      );
  tarch::mpi::Rank::getInstance().allReduce(
      &newMaxTimeStepSize,
      &_maxTimeStepSizeThisTimeStep,
      1, MPI_DOUBLE,
      MPI_MAX,
      [&]() -> void { tarch::services::ServiceRepository::getInstance().receiveDanglingMessages(); }
      );
 
  tarch::mpi::Rank::getInstance().allReduce(
      &newLocalMinTimeStamp,
      &_localMinTimeStampThisTimeStep,
      1, MPI_DOUBLE,
      MPI_MIN,
      [&]() -> void { tarch::services::ServiceRepository::getInstance().receiveDanglingMessages(); }
      );
  tarch::mpi::Rank::getInstance().allReduce(
      &newLocalMaxTimeStamp,
      &_localMaxTimeStampThisTimeStep,
      1, MPI_DOUBLE,
      MPI_MAX,
      [&]() -> void { tarch::services::ServiceRepository::getInstance().receiveDanglingMessages(); }
      );
 
  tarch::mpi::Rank::getInstance().allReduce(
      &newPatchUpdates,
      &_patchUpdates,
      1, MPI_INT,
      MPI_SUM,
      [&]() -> void { tarch::services::ServiceRepository::getInstance().receiveDanglingMessages(); }
      );
  #endif
 
    if ( _minTimeStampThisTimeStep>=std::numeric_limits<double>::max()/10.0) {
      _minTimeStampThisTimeStep = 0.0;
    }
    if ( _maxTimeStampThisTimeStep>=std::numeric_limits<double>::max()/10.0 ) {
      _maxTimeStampThisTimeStep = 0.0;
    }
  }
 
  {{FINISH_TIME_STEP_IMPLEMENTATION}}
}
 
 
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::startPlottingStep(
  double globalMinTimeStamp,
  double globalMaxTimeStamp,
  double globalMinTimeStepSize,
  double globalMaxTimeStepSize
) {
  if ( _solverState==SolverState::GridInitialisation ) {
    _solverState = SolverState::PlottingAfterGridInitialisation;
  }
  else {
    _solverState = SolverState::Plotting;
  }
}
 
 
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::finishPlottingStep() {
}
 
 
std::string {{NAMESPACE | join("::")}}::{{CLASSNAME}}::toString(SolverState state) {
  switch (state) {
    case SolverState::GridConstruction:
      return "grid-construction";
    case SolverState::Suspended:
      return "suspended";
    case SolverState::GridInitialisation:
      return "grid-initialisation";
    case SolverState::Plotting:
      return "plotting";
    case SolverState::PlottingAfterGridInitialisation:
      return "plotting-after-grid-initialisation";
    case SolverState::RungeKuttaSubStep0AfterGridInitialisation:
      return "runge-kutta-step-0-after-grid-initialisation";
    {% for item in range(RK_STEPS) %}
    case SolverState::RungeKuttaSubStep{{item}}:
      return "runge-kutta-step-{{item}}";
      break;
    {%- endfor %}
  }
  return "<undef>";
}
 
 
bool {{NAMESPACE | join("::")}}::{{CLASSNAME}}::mayPlot() const {
  return isLastGridSweepOfTimeStep()
      or _solverState==SolverState::GridInitialisation
      or _solverState==SolverState::Plotting;
}
 
 
bool {{NAMESPACE | join("::")}}::{{CLASSNAME}}::isFirstGridSweepOfTimeStep() const {
  return _solverState==SolverState::RungeKuttaSubStep0
      or _solverState==SolverState::RungeKuttaSubStep0AfterGridInitialisation;
}
 
 
bool {{NAMESPACE | join("::")}}::{{CLASSNAME}}::isLastGridSweepOfTimeStep() const {
  return _solverState==SolverState::RungeKuttaSubStep{{RK_STEPS-1}};
}
 
 
{{ABSTRACT_SOLVER_USER_DEFINITIONS}}
 
 
 
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::startSimulation() {
 
}
 
 
void {{NAMESPACE | join("::")}}::{{CLASSNAME}}::finishSimulation() {
 
}