Method initializeDynamicRuptureMatrices

src/Initializer/CellLocalMatrices.cpp:459–774

};

void seissol::initializer::initializeDynamicRuptureMatrices( seissol::geometry::MeshReader const&      i_meshReader,
                                                              LTSTree*               io_ltsTree,
                                                              LTS*                   i_lts,
                                                              Lut*                   i_ltsLut,
                                                              LTSTree*               dynRupTree,
                                                              DynamicRupture*        dynRup,
                                                              unsigned*              ltsFaceToMeshFace,
                                                              GlobalData const&      global,
                                                              double etaHack )
{
  real TData[tensor::T::size()];
  real TinvData[tensor::Tinv::size()];
  real APlusData[tensor::star::size(0)];
  real AMinusData[tensor::star::size(0)];

  std::vector<Fault> const& fault = i_meshReader.getFault();
  std::vector<Element> const& elements = i_meshReader.getElements();
  CellDRMapping (*drMapping)[4] = io_ltsTree->var(i_lts->drMapping);
  CellDRMapping (*drMappingDevice)[4] = io_ltsTree->var(i_lts->drMappingDevice);
  CellMaterialData* material = io_ltsTree->var(i_lts->material);
  real** derivatives = io_ltsTree->var(i_lts->derivatives);
  real* (*faceNeighbors)[4] = io_ltsTree->var(i_lts->faceNeighbors);
  real** derivativesDevice = io_ltsTree->var(i_lts->derivativesDevice);
  real* (*faceNeighborsDevice)[4] = io_ltsTree->var(i_lts->faceNeighborsDevice);
  CellLocalInformation* cellInformation = io_ltsTree->var(i_lts->cellInformation);

  unsigned* layerLtsFaceToMeshFace = ltsFaceToMeshFace;

  for (auto& layer : dynRupTree->leaves(Ghost)) {
    real**                                timeDerivativePlus                                        = layer.var(dynRup->timeDerivativePlus);
    real**                                timeDerivativeMinus                                       = layer.var(dynRup->timeDerivativeMinus);
    real**                                timeDerivativePlusDevice                                        = layer.var(dynRup->timeDerivativePlusDevice);
    real**                                timeDerivativeMinusDevice                                       = layer.var(dynRup->timeDerivativeMinusDevice);
    DRGodunovData*                        godunovData                                               = layer.var(dynRup->godunovData);
    real                                (*imposedStatePlus)[tensor::QInterpolated::size()]          = layer.var(dynRup->imposedStatePlus, AllocationPlace::Host);
    real                                (*imposedStateMinus)[tensor::QInterpolated::size()]         = layer.var(dynRup->imposedStateMinus, AllocationPlace::Host);
    real                                (*fluxSolverPlus)[tensor::fluxSolver::size()]               = layer.var(dynRup->fluxSolverPlus, AllocationPlace::Host);
    real                                (*fluxSolverMinus)[tensor::fluxSolver::size()]              = layer.var(dynRup->fluxSolverMinus, AllocationPlace::Host);
    real                                (*imposedStatePlusDevice)[tensor::QInterpolated::size()]          = layer.var(dynRup->imposedStatePlus, AllocationPlace::Device);
    real                                (*imposedStateMinusDevice)[tensor::QInterpolated::size()]         = layer.var(dynRup->imposedStateMinus, AllocationPlace::Device);
    real                                (*fluxSolverPlusDevice)[tensor::fluxSolver::size()]               = layer.var(dynRup->fluxSolverPlus, AllocationPlace::Device);
    real                                (*fluxSolverMinusDevice)[tensor::fluxSolver::size()]              = layer.var(dynRup->fluxSolverMinus, AllocationPlace::Device);
    DRFaceInformation*                    faceInformation                                           = layer.var(dynRup->faceInformation);
    seissol::model::IsotropicWaveSpeeds*  waveSpeedsPlus                                            = layer.var(dynRup->waveSpeedsPlus);
    seissol::model::IsotropicWaveSpeeds*  waveSpeedsMinus                                           = layer.var(dynRup->waveSpeedsMinus);
    seissol::dr::ImpedancesAndEta*        impAndEta                                                 = layer.var(dynRup->impAndEta);
    seissol::dr::ImpedanceMatrices*       impedanceMatrices                                         = layer.var(dynRup->impedanceMatrices);


#ifdef _OPENMP
  #pragma omp parallel for private(TData, TinvData, APlusData, AMinusData) schedule(static)
#endif
    for (unsigned ltsFace = 0; ltsFace < layer.getNumberOfCells(); ++ltsFace) {
      unsigned meshFace = layerLtsFaceToMeshFace[ltsFace];
      assert(fault[meshFace].element >= 0 || fault[meshFace].neighborElement >= 0);

      /// Face information
      faceInformation[ltsFace].meshFace = meshFace;