SamplingData< nDim, SolverType > Class Template Reference

This base class is responsible for the sampling data feature that provides the required general methods needed for sampling data at distinct points. More...

#include <samplingdata.h>

Inheritance diagram for SamplingData< nDim, SolverType >:

Collaboration diagram for SamplingData< nDim, SolverType >:

Public Member Functions
	SamplingData (SolverType &solver)
void	init ()
void	save (MBool finalTimestep)
void	setInputOutputProperties ()
MBool	enabled () const
virtual MInt	loadInputFile (const MString NotUsed(inputFileName), const MInt NotUsed(fileNo), std::vector< MFloat > &NotUsed(coordinates))
virtual MString	getBaseName () const
	Return base name of derived class, e.g. used for naming properties and output file names. More...

Protected Member Functions
void	init (SamplingDataSeries &timeSeries)
void	save (SamplingDataSeries &timeSeries, MBool finalTimeStep)
	Saves sampling data Saves all states of conservative variables at specified points. Also writes sampling data files in specified intervals. More...
virtual MBool	hasInputDataFile (const MInt fileNo) const
	Check if an input file exists, might be overwritten in derived class if points are generated. More...
virtual MBool	generatePoints () const
	Default behaviour for point generation which can be overwritten in derived class. More...
virtual MInt	noSamplingPoints (const MInt NotUsed(fileNo))
	Return the number of sampling points which are generated by the derived class (no input file) More...
virtual void	getSamplingPoint (const MInt NotUsed(fileNo), const MInt NotUsed(pointId), MFloat *const NotUsed(coordinates))
	Return the coordinates of the sampling point with the given id when generating points. More...
void	indexXD (const MInt index, const MInt *const nPoints, MInt *indexXD)
	Convert 1D index into 2D/3D index. More...
void	saveSamplingPointCoordinates (SamplingDataSeries &timeSeries)
	Save point coordinates of time series. More...
virtual void	readAdditionalProperties (const MInt NotUsed(fileNo))
	Read additional properties which are required by the derived class. More...
MInt	solverId () const
MBool	isMpiRoot () const
MPI_Comm	mpiComm () const
MInt	noVars () const
	Return total number of sampling variables. More...
MString	getFileBaseName () const
	Return base name of files for used sampling feature. More...
SolverType &	solver ()
	Return reference to solver. More...
const SolverType &	solver () const

Protected Attributes
MBool	m_enabled = false
MBool	m_interpolatePointData
std::vector< SamplingDataSeries >	m_timeSeries {}
	Holds all properties and buffers for each input file. More...
std::vector< MInt >	m_solverSamplingVarIds {}
	List of sampling variables. More...
std::vector< MInt >	m_noSolverSamplingVars {}
	Number of variables for each sampling variable. More...
std::vector< std::vector< MString > >	m_solverSamplingVarNames {}
	List of variable names for each sampling variable. More...

Private Attributes
SolverType &	m_solver

Detailed Description

template<MInt nDim, class SolverType>
class SamplingData< nDim, SolverType >

Author

Fabian Klemp f.kle.nosp@m.mp@a.nosp@m.ia.rw.nosp@m.th-a.nosp@m.achen.nosp@m..de

Date

2016-09-27

Definition at line 35 of file samplingdata.h.

Constructor & Destructor Documentation

SamplingData()

template<MInt nDim, class SolverType >

SamplingData< nDim, SolverType >::SamplingData ( SolverType &  solver )

inline

Definition at line 37 of file samplingdata.h.

: m_solver(solver) {}

Member Function Documentation

enabled()

template<MInt nDim, class SolverType >

MBool SamplingData< nDim, SolverType >::enabled ( ) const

inline

Definition at line 42 of file samplingdata.h.

{ return m_enabled; }

generatePoints()

template<MInt nDim, class SolverType >

virtual MBool SamplingData< nDim, SolverType >::generatePoints ( ) const

inlineprotectedvirtual

Reimplemented in VolumeData< nDim, SolverType >, VolumeData< nDim, DgCartesianSolver< nDim, SysEqn > >, VolumeData< nDim, FvCartesianSolverXD< nDim, SysEqn > >, and VolumeData< nDim, LbSolver >.

Definition at line 75 of file samplingdata.h.

{ return false; }

getBaseName()

template<MInt nDim, class SolverType >

virtual MString SamplingData< nDim, SolverType >::getBaseName ( ) const

inlinevirtual

Reimplemented in PointData< nDim, SolverType >, PointData< nDim, DgCartesianSolver< nDim, SysEqn > >, PointData< nDim, FvCartesianSolverXD< nDim, SysEqn > >, PointData< nDim, LbSolver >, SurfaceData< nDim, SolverType >, SurfaceData< nDim, DgCartesianSolver< nDim, SysEqn > >, SurfaceData< nDim, FvCartesianSolverXD< nDim, SysEqn > >, SurfaceData< nDim, LbSolver >, VolumeData< nDim, SolverType >, VolumeData< nDim, DgCartesianSolver< nDim, SysEqn > >, VolumeData< nDim, FvCartesianSolverXD< nDim, SysEqn > >, and VolumeData< nDim, LbSolver >.

Definition at line 53 of file samplingdata.h.

                                      {
    TERMM(1, "ERROR: not implemented in base class");
    return "";
  }

getFileBaseName()

template<MInt nDim, class SolverType >

MString SamplingData< nDim, SolverType >::getFileBaseName ( ) const

inlineprotected

Definition at line 126 of file samplingdata.h.

{ return getBaseName() + "_data_"; }

getSamplingPoint()

template<MInt nDim, class SolverType >

virtual void SamplingData< nDim, SolverType >::getSamplingPoint ( const MInt   NotUsedfileNo, const MInt   NotUsedpointId, MFloat *const   NotUsedcoordinates  )

inlineprotectedvirtual

Definition at line 103 of file samplingdata.h.

                                                                    {
    TERMM(1, "Not implemented in derived class!");
  }

hasInputDataFile()

template<MInt nDim, class SolverType >

virtual MBool SamplingData< nDim, SolverType >::hasInputDataFile ( const MInt  fileNo ) const

inlineprotectedvirtual

Reimplemented in VolumeData< nDim, SolverType >, VolumeData< nDim, DgCartesianSolver< nDim, SysEqn > >, VolumeData< nDim, FvCartesianSolverXD< nDim, SysEqn > >, and VolumeData< nDim, LbSolver >.

Definition at line 64 of file samplingdata.h.

                                                          {
    const MString propName = getBaseName() + "DataFileName_" + std::to_string(fileNo);
    if(Context::propertyExists(propName, solverId())) {
      const MString fileName = Context::getSolverProperty<MString>(propName, solverId(), AT_);
      return fileName != "";
    } else {
      return false;
    }
  }

indexXD()

template<MInt nDim, class SolverType >

void SamplingData< nDim, SolverType >::indexXD ( const MInt  index, const MInt *const  nPoints, MInt *  indexXD  )

protected

Convert a one dimensional index into a 2D/3D index with the given number of points per direction (index of last dimension varies the fastest)

Definition at line 1023 of file samplingdata.h.

                                                                                                       {
  const MInt noPoints2D = nPoints[0] * nPoints[1];
  const MInt maxNoPoints = (nDim == 2) ? noPoints2D : noPoints2D * nPoints[2];
  if(index < 0 || index >= maxNoPoints) {
    TERMM(1, "invalid index: " + std::to_string(index));
  }
 
  // Last index (varies fastest)
  indexXD[nDim - 1] = index % nPoints[nDim - 1];
  // First/second index in 2D/3D
  indexXD[nDim - 2] = ((index - indexXD[nDim - 1]) / nPoints[nDim - 1]) % nPoints[nDim - 2];
  IF_CONSTEXPR(nDim == 3) {
    // First index in 3D
    indexXD[0] = (index - indexXD[1] * nPoints[1] - indexXD[2]) / (nPoints[2] * nPoints[1]);
  }
}

init() [1/2]

template<MInt nDim, class SolverType >

void SamplingData< nDim, SolverType >::init

Definition at line 571 of file samplingdata.h.

                                          {
  if(!m_enabled || !solver().isActive()) {
    return;
  }
 
  for(MUint i = 0; i < m_timeSeries.size(); i++) {
    init(m_timeSeries[i]);
  }
}

init() [2/2]

template<MInt nDim, class SolverType >

void SamplingData< nDim, SolverType >::init ( SamplingDataSeries &  timeSeries )

protected

Definition at line 590 of file samplingdata.h.

                                                                        {
  using namespace maia;
  using namespace std;
 
 
  // save points and corresponding data in this container for sorting
  // tuple contains following data in that order:
  // coordinates of point, element id of point, line index in input file
  std::vector<std::tuple<std::array<MFloat, nDim>, MInt, MInt>> localTmpPoints;
 
  m_log << "Solver #" << solverId() << ": ";
  // Check for an input file
  if(timeSeries.m_inputFileName != "") {
    m_log << "Loading points for " << getBaseName() << " data class...";
 
    std::vector<MFloat> tmpCoordinates{};
 
    // read in file on root using the I/O method of the derived class
    if(isMpiRoot()) {
      timeSeries.m_noPoints = loadInputFile(timeSeries.m_inputFileName, timeSeries.m_fileNo, tmpCoordinates);
    }
 
    MPI_Bcast(&timeSeries.m_noPoints, 1, type_traits<MInt>::mpiType(), 0, mpiComm(), AT_, "timeSeries.m_noPoints");
 
    if(timeSeries.m_noPoints <= 0) {
      TERMM(1, "ERROR: no sampling points found for input file '" + timeSeries.m_inputFileName + "'");
      m_log << "ERROR: no sampling points found for input file '" << timeSeries.m_inputFileName << "'" << std::endl;
      return;
    }
 
    ScratchSpace<MFloat> coordinatesScratch(timeSeries.m_noPoints * nDim, FUN_, "coordinatesScratch");
    if(isMpiRoot()) {
      copy(tmpCoordinates.begin(), tmpCoordinates.end(), coordinatesScratch.begin());
    }
    MPI_Bcast(&coordinatesScratch[0], timeSeries.m_noPoints * nDim, type_traits<MFloat>::mpiType(), 0, mpiComm(), AT_,
              "coordinatesScratch[0]");
 
    localTmpPoints.resize(timeSeries.m_noPoints);
 
    MFloatTensor tmpCoordinatesTensor(&coordinatesScratch[0], timeSeries.m_noPoints, nDim);
    std::array<MFloat, nDim> curCoordinates{};
    MInt curElementId = -1;
    for(MInt i = 0; i < timeSeries.m_noPoints; i++) {
      for(MInt j = 0; j < nDim; j++) {
        curCoordinates[j] = tmpCoordinatesTensor(i, j);
      }
      // check if the point is relevant to this rank, each point should be globally unique, i.e.,
      // appear only once
      curElementId = solver().getIdAtPoint(&curCoordinates[0], true);
 
      if(curElementId != -1) {
        localTmpPoints[timeSeries.m_noLocalPoints] = make_tuple(curCoordinates, curElementId, i);
        timeSeries.m_noLocalPoints++;
 
        // Precompute interpolation information if possible/required
        solver().initInterpolationForCell(curElementId);
      }
    }
 
    // sanity check if each point ended up on exactly one rank:
    ScratchSpace<MInt> loadedPointsVector(timeSeries.m_noPoints, FUN_, "loadedPointsVector");
    loadedPointsVector.fill(0);
    for(MInt i = 0; i < timeSeries.m_noLocalPoints; i++) {
      loadedPointsVector[get<2>(localTmpPoints[i])]++;
    }
 
    if(isMpiRoot()) {
      MPI_Reduce(MPI_IN_PLACE, &loadedPointsVector[0], loadedPointsVector.size(), type_traits<MInt>::mpiType(), MPI_SUM,
                 0, mpiComm(), AT_, "MPI_IN_PLACE", "loadedPointsVector[0]");
      for(MUint i = 0; i < loadedPointsVector.size(); i++) {
        if(loadedPointsVector[i] != 1) {
          std::ostringstream err;
          err << "Point at line " << i + 1 << " of file " // start line count at one
              << timeSeries.m_inputFileName << " was loaded " << loadedPointsVector[i] << " times. Aborting."
              << std::endl;
          TERMM(1, err.str());
        }
      }
    } else {
      MPI_Reduce(&loadedPointsVector[0], nullptr, loadedPointsVector.size(), type_traits<MInt>::mpiType(), MPI_SUM, 0,
                 mpiComm(), AT_, "loadedPointsVector[0]", "nullptr");
    }
  } else if(timeSeries.m_generatePoints) {
    m_log << "Generating points for " << getBaseName() << " data class...";
    // Generate points based on given properties
    // Get number of points to generate
    const MInt noPoints = noSamplingPoints(timeSeries.m_fileNo);
    timeSeries.m_noPoints = noPoints;
 
    std::array<MFloat, nDim> coordinates{};
    MLong checkSum = 0;
 
    for(MInt i = 0; i < noPoints; i++) {
      // Get the coordinates for this sampling point id
      getSamplingPoint(timeSeries.m_fileNo, i, &coordinates[0]);
 
      // Check if point is relevant for local rank
      const MInt curElementId = solver().getIdAtPoint(&coordinates[0], true);
      if(curElementId != -1) {
        localTmpPoints.push_back(make_tuple(coordinates, curElementId, i));
        timeSeries.m_noLocalPoints++;
        checkSum += (i + 1);
 
        // Precompute interpolation information if possible/required
        solver().initInterpolationForCell(curElementId);
      }
    }
 
    MInt globalNoPoints = timeSeries.m_noLocalPoints;
    MPI_Allreduce(MPI_IN_PLACE, &globalNoPoints, 1, type_traits<MInt>::mpiType(), MPI_SUM, mpiComm(), AT_,
                  "MPI_IN_PLACE", "globalNoPoints");
    if(globalNoPoints != noPoints) {
      TERMM(1, "Global number of points does not match " + std::to_string(globalNoPoints)
                   + " != " + std::to_string(noPoints));
    }
 
    // Use Gauss sum equation to check if all points ended up on exactly one rank
    // Note: this is a heuristic check but it should be very unlikely to get a false positive; this
    // will work up to 2 billion sampling points
    MPI_Allreduce(MPI_IN_PLACE, &checkSum, 1, type_traits<MLong>::mpiType(), MPI_SUM, mpiComm(), AT_, "MPI_IN_PLACE",
                  "checkSum");
    const MLong noPointsTmp = noPoints;
    if(2 * checkSum != (noPointsTmp * noPointsTmp + noPointsTmp)) {
      TERMM(1, "Incorrect check sum for " + std::to_string(noPoints) + " points: " + std::to_string(checkSum));
    }
  } else {
    m_log << "Warning: no sampling data input file specified and generation of points disabled." << std::endl;
    return;
  }
 
  // Create a new MPI Communicator
  ScratchSpace<MInt> pointDataScratch(globalNoDomains(), FUN_, "pointDataScratch");
  std::fill_n(&pointDataScratch[0], globalNoDomains(), -1);
 
  // Get Rank if points are relevant
  MInt rank = -1;
  if(timeSeries.m_noLocalPoints > 0) {
    MPI_Comm_rank(mpiComm(), &rank);
  }
 
  // Combine all ranks to get relevant ranks
  MPI_Allgather(&rank, 1, type_traits<MInt>::mpiType(), &pointDataScratch[0], 1, type_traits<MInt>::mpiType(),
                mpiComm(), AT_, "rank", "pointDataScratch[0]");
 
  // Check for relevant ranks and save them to create the new communicator
  const MInt noPdDomains =
      count_if(pointDataScratch.begin(), pointDataScratch.end(), [](const MInt a) { return a != -1; });
  MIntScratchSpace pdDomains(noPdDomains, AT_, "pdDomains");
  MInt position = 0;
  for(MUint i = 0; i < pointDataScratch.size(); i++) {
    if(pointDataScratch[i] != -1) {
      pdDomains[position] = pointDataScratch[i];
      position++;
    }
  }
 
  // Create new point data mpi group
  MPI_Group globalGroup, localGroup;
  MPI_Comm_group(mpiComm(), &globalGroup, AT_, "globalGroup");
  MPI_Group_incl(globalGroup, noPdDomains, &pdDomains[0], &localGroup, AT_);
 
  // Create new communicator and clean up
  MPI_Comm_create(mpiComm(), localGroup, &timeSeries.m_mpiComm, AT_, "timeSeries.m_mpiComm");
  MPI_Group_free(&globalGroup, AT_);
  MPI_Group_free(&localGroup, AT_);
 
  // Set root for new communicator
  if(timeSeries.m_noLocalPoints > 0) {
    MInt subRank = -1;
    MPI_Comm_rank(timeSeries.m_mpiComm, &subRank);
    if(subRank == 0) {
      timeSeries.m_isMpiRoot = true;
    }
  }
 
 
  // sort points by global id
  stable_sort(localTmpPoints.begin(),
              localTmpPoints.begin() + timeSeries.m_noLocalPoints,
              // Read all input files and save all points into one array
              [this](const std::tuple<std::array<MFloat, nDim>, MInt, MInt>& a,
                     const std::tuple<std::array<MFloat, nDim>, MInt, MInt>& b) {
                return solver().grid().tree().globalId(solver().getCellIdByIndex(get<1>(a)))
                       < solver().grid().tree().globalId(solver().getCellIdByIndex(get<1>(b)));
              });
 
  // split localTmpPoints for it to be easier to handle later on
  if(timeSeries.m_noLocalPoints > 0) {
    timeSeries.m_coordinates.resize(timeSeries.m_noLocalPoints * nDim);
    MFloatTensor coordinatesTensor(&timeSeries.m_coordinates[0], timeSeries.m_noLocalPoints, nDim);
    timeSeries.m_elementIds.resize(timeSeries.m_noLocalPoints);
    timeSeries.m_sortIndex.resize(timeSeries.m_noLocalPoints);
    for(MInt i = 0; i < timeSeries.m_noLocalPoints; i++) {
      for(MInt j = 0; j < nDim; j++) {
        coordinatesTensor(i, j) = get<0>(localTmpPoints[i])[j];
      }
      timeSeries.m_elementIds[i] = get<1>(localTmpPoints[i]);
      timeSeries.m_sortIndex[i] = get<2>(localTmpPoints[i]);
    }
  }
 
  // Save sampling point coordinates if they were not loaded from a file
  if(timeSeries.m_inputFileName == "" && timeSeries.m_noLocalPoints > 0) {
    saveSamplingPointCoordinates(timeSeries);
  }
 
  // Allocate memory to save/buffer point data inbetween writes
  timeSeries.m_maxNoSample = std::max(1, timeSeries.m_writeInterval / timeSeries.m_sampleInterval);
  if(timeSeries.m_noLocalPoints > 0) {
    // Note: +1 for initial data file containing also timestep 0
    timeSeries.m_stateBuffer.resize(timeSeries.m_noLocalPoints * (timeSeries.m_maxNoSample + 1) * noVars());
  }
 
  // Allocate storage for the time step count and the simulation time
  if(timeSeries.m_isMpiRoot) {
    timeSeries.m_timeStepBuffer.resize(timeSeries.m_maxNoSample + 1);
    timeSeries.m_timeBuffer.resize(timeSeries.m_maxNoSample + 1);
  } else {
    timeSeries.m_timeStepBuffer.resize(1);
    timeSeries.m_timeBuffer.resize(1);
  }
 
  m_log << " done" << std::endl << "#" << timeSeries.m_fileNo;
  if(timeSeries.m_inputFileName != "") {
    m_log << " Loaded " << timeSeries.m_noPoints << " points from " << timeSeries.m_inputFileName << " for "
          << getBaseName() << " sampling feature successfully." << std::endl;
  } else if(timeSeries.m_generatePoints) {
    m_log << " Generated " << timeSeries.m_noPoints << " points for " << getBaseName()
          << " sampling feature successfully." << std::endl;
  }
  m_log << "sampleInterval = " << timeSeries.m_sampleInterval << "; startTimeStep = " << timeSeries.m_startTimeStep
        << "; endTimeStep = " << timeSeries.m_endTimeStep << "; writeInterval = " << timeSeries.m_writeInterval
        << std::endl;
}

isMpiRoot()

template<MInt nDim, class SolverType >

MBool SamplingData< nDim, SolverType >::isMpiRoot ( ) const

inlineprotected

Definition at line 119 of file samplingdata.h.

{ return m_solver.domainId() == 0; }

loadInputFile()

template<MInt nDim, class SolverType >

virtual MInt SamplingData< nDim, SolverType >::loadInputFile ( const MString   NotUsedinputFileName, const MInt   NotUsedfileNo, std::vector< MFloat > &  NotUsedcoordinates  )

inlinevirtual

*** Virtual functions that need to be implemented in derived class *** Load input file and determine number of sampling points and point coordinates

Definition at line 46 of file samplingdata.h.

                                                                      {
    TERMM(1, "ERROR: not implemented in base class");
    return -1;
  }

mpiComm()

template<MInt nDim, class SolverType >

MPI_Comm SamplingData< nDim, SolverType >::mpiComm ( ) const

inlineprotected

Definition at line 120 of file samplingdata.h.

{ return m_solver.mpiComm(); }

noSamplingPoints()

template<MInt nDim, class SolverType >

virtual MInt SamplingData< nDim, SolverType >::noSamplingPoints ( const MInt   NotUsedfileNo )

inlineprotectedvirtual

Definition at line 97 of file samplingdata.h.

                                                            {
    TERMM(1, "Not implemented in derived class!");
    return -1;
  }

noVars()

template<MInt nDim, class SolverType >

MInt SamplingData< nDim, SolverType >::noVars ( ) const

inlineprotected

Definition at line 124 of file samplingdata.h.

{ return std::accumulate(m_noSolverSamplingVars.begin(), m_noSolverSamplingVars.end(), 0); }

readAdditionalProperties()

template<MInt nDim, class SolverType >

virtual void SamplingData< nDim, SolverType >::readAdditionalProperties ( const MInt   NotUsedfileNo )

inlineprotectedvirtual

Definition at line 115 of file samplingdata.h.

{};

save() [1/2]

template<MInt nDim, class SolverType >

void SamplingData< nDim, SolverType >::save

(

MBool

finalTimestep

)

Definition at line 832 of file samplingdata.h.

                                                             {
  if(!m_enabled) {
    return;
  }
 
  for(MUint i = 0; i < m_timeSeries.size(); i++) {
    save(m_timeSeries[i], finalTimeStep);
  }
}

save() [2/2]

template<MInt nDim, class SolverType >

void SamplingData< nDim, SolverType >::save ( SamplingDataSeries &  timeSeries, MBool  finalTimeStep  )

protected

Author

Fabian Klemp f.kle.nosp@m.mp@a.nosp@m.ia.rw.nosp@m.th-a.nosp@m.achen.nosp@m..de

Date

2016-07-22

Assumes beeing called every timestep including timeStep 0 and the finalTimeStep.

Definition at line 855 of file samplingdata.h.

                                                                                             {
  using namespace std;
 
  // Check if time step inverval is reached
  if(solver().getCurrentTimeStep() < timeSeries.m_startTimeStep
     || solver().getCurrentTimeStep() > timeSeries.m_endTimeStep) {
    return;
  }
 
  // Check for correct timestep interval
  if(!(solver().getCurrentTimeStep() % timeSeries.m_sampleInterval == 0 || finalTimeStep)) {
    return;
  }
 
  // Compute sampling variables and exchange data for interpolation (if not done already in this time step)
  solver().calcSamplingVariables(m_solverSamplingVarIds, true);
 
  // Check if rank belongs to an sampling data MPI communicator
  // Note: sampling variables are computed by all ranks as the interpolation might require halo cell
  // values
  if(timeSeries.m_noLocalPoints < 1 || timeSeries.m_noPoints == 0) {
    return;
  }
 
  // buffer time and timestep
  if(timeSeries.m_isMpiRoot) {
    timeSeries.m_timeStepBuffer[timeSeries.m_sampleIndex] = solver().getCurrentTimeStep();
    timeSeries.m_timeBuffer[timeSeries.m_sampleIndex] = solver().time();
  }
 
  // buffer point data states
  MFloatTensor coordinatesTensorBuf(&timeSeries.m_coordinates[0], timeSeries.m_noLocalPoints, nDim);
  MFloatTensor stateTensorBuf(&timeSeries.m_stateBuffer[0], timeSeries.m_noLocalPoints, timeSeries.m_maxNoSample + 1,
                              noVars());
 
  MBool interpolate = m_interpolatePointData;
  for(MInt i = 0; i < timeSeries.m_noLocalPoints; i++) {
    MInt offset = 0;
    for(MUint var = 0; var < m_solverSamplingVarIds.size(); var++) {
      solver().calcSamplingVarAtPoint(&coordinatesTensorBuf(i, 0), timeSeries.m_elementIds[i],
                                      m_solverSamplingVarIds[var], &stateTensorBuf(i, timeSeries.m_sampleIndex, offset),
                                      interpolate);
      offset += m_noSolverSamplingVars[var];
    }
  }
 
  // increment sample index
  timeSeries.m_sampleIndex++;
 
  MBool isFirstDataFile = false;
  // Note: first data file for simulation starting at timestep 0 contains one value more than
  // all other data files; at a restart the restart time step is not considered here even if it is
  // equal to m_startTimeStep
  // if (solver().getCurrentTimeStep() >= timeSeries.m_startTimeStep
  //   && solver().getCurrentTimeStep() <= timeSeries.m_startTimeStep + timeSeries.m_writeInterval) {
  if(solver().getCurrentTimeStep() >= 0 && solver().getCurrentTimeStep() <= timeSeries.m_writeInterval) {
    isFirstDataFile = true;
  }
 
  const MBool isRestartStep =
      (solver().restartInterval() > 0 && solver().getCurrentTimeStep() % solver().restartInterval() == 0);
  const MBool isFirstSample = (isFirstDataFile && timeSeries.m_sampleIndex == 1);
  const MBool writeAtRestartStep = (isRestartStep && !isFirstSample);
  const MBool writeFirstDataFile = (timeSeries.m_sampleIndex > timeSeries.m_maxNoSample && isFirstDataFile);
  const MBool writeDataFile = ((timeSeries.m_sampleIndex + 1 > timeSeries.m_maxNoSample) && !isFirstDataFile);
 
  // write file if it's time to do so otherwise return
  if(!(writeDataFile || finalTimeStep || writeAtRestartStep || writeFirstDataFile)) {
    return;
  }
 
  using namespace maia::parallel_io;
 
  // calculate how many samples have been recorded in total
  MInt sampleCount = (solver().getCurrentTimeStep() / timeSeries.m_sampleInterval) + 1; // also saving at timeStep 0
 
  // Correct number of samples if finalTimeStep is not in the writeInterval
  if(finalTimeStep && !(solver().getCurrentTimeStep() % timeSeries.m_writeInterval == 0)) {
    sampleCount++;
  }
 
  MInt fileNo = 0;
  if(!finalTimeStep || (solver().getCurrentTimeStep() % timeSeries.m_writeInterval == 0)) {
    fileNo = (sampleCount - 1) / (timeSeries.m_writeInterval / timeSeries.m_sampleInterval);
  } else {
    fileNo = ceil((static_cast<MFloat>(sampleCount - 1)
                   / static_cast<MFloat>(timeSeries.m_writeInterval / timeSeries.m_sampleInterval)));
  }
 
  // Output file name
  std::ostringstream fileName;
  MString baseFileName = getFileBaseName();
  fileName << solver().outputDir() << baseFileName << std::to_string(timeSeries.m_fileNo) << "_" << setw(8)
           << setfill('0') << fileNo << ParallelIo::fileExt();
 
  ParallelIo file(fileName.str(), PIO_REPLACE, timeSeries.m_mpiComm);
 
  file.defineArray(PIO_INT, "timeStep", timeSeries.m_sampleIndex);
  file.defineArray(PIO_FLOAT, "time", timeSeries.m_sampleIndex);
  file.defineArray(PIO_INT, "sortIndex", timeSeries.m_noPoints);
 
  file.setAttribute("si[i] = index in input", "description", "sortIndex");
 
  ParallelIo::size_type dimSizes[] = {timeSeries.m_noPoints, timeSeries.m_sampleIndex, noVars()};
  file.defineArray(PIO_FLOAT, "pointStates", 3, &dimSizes[0]);
 
  std::ostringstream descPs;
  // Please be aware that these attributes get parsed in the corresponding post
  // processing script. So if you change them please check if it breaks the
  // script
  file.setAttribute(descPs.str(), "description", "pointStates");
  file.setAttribute("point index", "dim_0", "pointStates");
  file.setAttribute("sample index", "dim_1", "pointStates");
  file.setAttribute("var index", "dim_2", "pointStates");
  file.setAttribute(timeSeries.m_inputFileName, "inputFileName", "pointStates");
 
  MInt varOffset = 0;
  for(MUint var = 0; var < m_solverSamplingVarIds.size(); var++) {
    for(MInt i = 0; i < m_noSolverSamplingVars[var]; i++) {
      file.setAttribute(m_solverSamplingVarNames[var][i], "var_" + std::to_string(varOffset + i), "pointStates");
    }
    varOffset += m_noSolverSamplingVars[var];
  }
 
  // ensure data is consecutive in memory
  if(timeSeries.m_sampleIndex < timeSeries.m_maxNoSample + 1) {
    MFloatTensor stateTensor(&timeSeries.m_stateBuffer[0], timeSeries.m_noLocalPoints, timeSeries.m_maxNoSample + 1,
                             noVars());
    MFloatMatrix finalStateTensor(&timeSeries.m_stateBuffer[0], timeSeries.m_noLocalPoints, timeSeries.m_sampleIndex,
                                  noVars());
    for(MInt i = 0; i < timeSeries.m_noLocalPoints; i++) {
      for(MInt j = 0; j < timeSeries.m_sampleIndex; j++) {
        for(MInt k = 0; k < noVars(); k++) {
          finalStateTensor(i, j, k) = stateTensor(i, j, k);
        }
      }
    }
  }
 
  // root writes timestep and time data
  if(timeSeries.m_isMpiRoot) {
    file.setOffset(timeSeries.m_sampleIndex, 0);
  } else {
    file.setOffset(0, 0);
  }
  file.writeArray(&timeSeries.m_timeStepBuffer[0], "timeStep");
  file.writeArray(&timeSeries.m_timeBuffer[0], "time");
 
  // calculate offset for reverse index and pointStates
  ParallelIo::size_type offset, totalCount;
  file.calcOffset(timeSeries.m_noLocalPoints, &offset, &totalCount, timeSeries.m_mpiComm);
 
  // write sort index
  file.setOffset(timeSeries.m_noLocalPoints, offset);
  file.writeArray(&timeSeries.m_sortIndex[0], "sortIndex");
 
  // write state data to file
  file.setOffset(timeSeries.m_noLocalPoints, offset, 3);
  file.writeArray(&timeSeries.m_stateBuffer[0], "pointStates");
 
  // reset sample index
  timeSeries.m_sampleIndex = 0;
}

saveSamplingPointCoordinates()

template<MInt nDim, class SolverType >

void SamplingData< nDim, SolverType >::saveSamplingPointCoordinates ( SamplingDataSeries &  timeSeries )

protected

Save the point coordinates of the given time series, needs to be called only by the active communicator

Definition at line 1044 of file samplingdata.h.

                                                                                                {
  TRACE();
  using namespace maia::parallel_io;
 
  std::ostringstream fileName;
  fileName << solver().outputDir() << getFileBaseName() << "coordinates_" << timeSeries.m_fileNo
           << ParallelIo::fileExt();
  ParallelIo file(fileName.str(), PIO_REPLACE, timeSeries.m_mpiComm);
 
  // Add spatial dimensions
  file.setAttribute(nDim, "nDim");
 
  file.defineArray(PIO_INT, "sortIndex", timeSeries.m_noPoints);
  file.setAttribute("si[i] = index in input", "description", "sortIndex");
 
  ParallelIo::size_type dimSizes[] = {timeSeries.m_noPoints, nDim};
  file.defineArray(PIO_FLOAT, "coordinates", 2, &dimSizes[0]);
 
  ParallelIo::size_type offset, totalCount;
  file.calcOffset(timeSeries.m_noLocalPoints, &offset, &totalCount, timeSeries.m_mpiComm);
 
  // Write sort index
  file.setOffset(timeSeries.m_noLocalPoints, offset);
  file.writeArray(&timeSeries.m_sortIndex[0], "sortIndex");
 
  // Write coordinates
  file.setOffset(timeSeries.m_noLocalPoints, offset, 2);
  file.writeArray(&timeSeries.m_coordinates[0], "coordinates");
}

setInputOutputProperties()

template<MInt nDim, class SolverType >

void SamplingData< nDim, SolverType >::setInputOutputProperties

(

)

solver() [1/2]

template<MInt nDim, class SolverType >

SolverType & SamplingData< nDim, SolverType >::solver ( )

inlineprotected

Definition at line 129 of file samplingdata.h.

{ return m_solver; }

solver() [2/2]

template<MInt nDim, class SolverType >

const SolverType & SamplingData< nDim, SolverType >::solver ( ) const

inlineprotected

Definition at line 130 of file samplingdata.h.

{ return m_solver; }

solverId()

template<MInt nDim, class SolverType >

MInt SamplingData< nDim, SolverType >::solverId ( ) const

inlineprotected

Definition at line 118 of file samplingdata.h.

{ return m_solver.solverId(); }

Member Data Documentation

m_enabled

template<MInt nDim, class SolverType >

MBool SamplingData< nDim, SolverType >::m_enabled = false

protected

Definition at line 134 of file samplingdata.h.

m_interpolatePointData

template<MInt nDim, class SolverType >

MBool SamplingData< nDim, SolverType >::m_interpolatePointData

protected

Definition at line 135 of file samplingdata.h.

m_noSolverSamplingVars

template<MInt nDim, class SolverType >

std::vector<MInt> SamplingData< nDim, SolverType >::m_noSolverSamplingVars {}

protected

Definition at line 143 of file samplingdata.h.

m_solver

template<MInt nDim, class SolverType >

SolverType& SamplingData< nDim, SolverType >::m_solver

private

Definition at line 150 of file samplingdata.h.

m_solverSamplingVarIds

template<MInt nDim, class SolverType >

std::vector<MInt> SamplingData< nDim, SolverType >::m_solverSamplingVarIds {}

protected

Definition at line 141 of file samplingdata.h.

m_solverSamplingVarNames

template<MInt nDim, class SolverType >

std::vector<std::vector<MString> > SamplingData< nDim, SolverType >::m_solverSamplingVarNames {}

protected

Definition at line 145 of file samplingdata.h.

m_timeSeries

template<MInt nDim, class SolverType >

std::vector<SamplingDataSeries> SamplingData< nDim, SolverType >::m_timeSeries {}

protected

Definition at line 138 of file samplingdata.h.

---

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

• /home/gitlab-runner/scratch/builds/oxpnswJ6/1/aia/m-AIA/m-AIA/src/POST/samplingdata.h