src/SPH/SPHHydroBaseGSRZ.cc
//------------------------------------------------------------------------------
// Enforce the boundary conditions for hydro state fields.
//------------------------------------------------------------------------------
Enforce the boundary conditions for hydro state fields.
Enforce the boundary conditions for hydro state fields .
void
SPHHydroBaseGSRZ::
enforceBoundaries(State<Dim<2> >& state,
                  StateDerivatives<Dim<2> >& derivs) {
  FieldList<Dimension, Scalar> mass = state.fields(HydroFieldNames::mass, 0.0);
  FieldList<Dimension, Vector> pos = state.fields(HydroFieldNames::position, Vector::zero);
  const unsigned numNodeLists = mass.numFields();
  for (unsigned nodeListi = 0; nodeListi != numNodeLists; ++nodeListi) {
    const unsigned n = mass[nodeListi]->numInternalElements();
    for (unsigned i = 0; i != n; ++i) {
      const Scalar circi = 2.0*M_PI*abs(pos(nodeListi, i).y());
      CHECK(circi > 0.0);
      mass(nodeListi, i) /= circi;
    }
  }
  SPHHydroBase<Dim<2> >::enforceBoundaries(state, derivs);
  FieldList<Dimension, SymTensor> H = state.fields(HydroFieldNames::H, SymTensor::zero);
  for (unsigned nodeListi = 0; nodeListi != numNodeLists; ++nodeListi) {
    const unsigned n = mass[nodeListi]->numInternalElements();
    const Scalar nPerh = mass[nodeListi]->nodeList().nodesPerSmoothingScale();
    for (unsigned i = 0; i != n; ++i) {
      Vector& posi = pos(nodeListi, i);
      const SymTensor& Hi = H(nodeListi, i);
      const Scalar zetai = (Hi*posi).y();
      const Scalar ri = posi.y();
      const Scalar hrInvi = zetai*safeInvVar(ri);
      const Scalar rmin = 0.5/(nPerh*hrInvi);
      if (ri < rmin) posi.y(2.0*rmin - ri);
      const Scalar circi = 2.0*M_PI*abs(posi.y());
      mass(nodeListi, i) *= circi;
    }
  }
}
void SPHHydroBaseGSRZ :: enforceBoundaries ( State < Dim < 2 > > & state , StateDerivatives < Dim < 2 > > & derivs ) { FieldList < Dimension , Scalar > mass = state . fields ( HydroFieldNames :: mass , 0.0 ) ; FieldList < Dimension , Vector > pos = state . fields ( HydroFieldNames :: position , Vector :: zero ) ; const unsigned numNodeLists = mass . numFields ( ) ; for ( unsigned nodeListi = 0 ; nodeListi != numNodeLists ; ++ nodeListi ) { const unsigned n = mass [ nodeListi ] -> numInternalElements ( ) ; for ( unsigned i = 0 ; i != n ; ++ i ) { const Scalar circi = 2.0 * M_PI * abs ( pos ( nodeListi , i ) . y ( ) ) ; CHECK ( circi > 0.0 ) ; mass ( nodeListi , i ) /= circi ; } } SPHHydroBase < Dim < 2 > > :: enforceBoundaries ( state , derivs ) ; FieldList < Dimension , SymTensor > H = state . fields ( HydroFieldNames :: H , SymTensor :: zero ) ; for ( unsigned nodeListi = 0 ; nodeListi != numNodeLists ; ++ nodeListi ) { const unsigned n = mass [ nodeListi ] -> numInternalElements ( ) ; const Scalar nPerh = mass [ nodeListi ] -> nodeList ( ) . nodesPerSmoothingScale ( ) ; for ( unsigned i = 0 ; i != n ; ++ i ) { Vector & posi = pos ( nodeListi , i ) ; const SymTensor & Hi = H ( nodeListi , i ) ; const Scalar zetai = ( Hi * posi ) . y ( ) ; const Scalar ri = posi . y ( ) ; const Scalar hrInvi = zetai * safeInvVar ( ri ) ; const Scalar rmin = 0.5 / ( nPerh * hrInvi ) ; if ( ri < rmin ) posi . y ( 2.0 * rmin - ri ) ; const Scalar circi = 2.0 * M_PI * abs ( posi . y ( ) ) ; mass ( nodeListi , i ) *= circi ; } } }
Enforce the boundary conditions for hydro state fields.
Enforce the boundary conditions for hydro state fields .
// Convert the mass to mass/length before BCs are applied. // Apply ordinary SPH BCs. // Scale back to mass. // We also ensure no point approaches the z-axis too closely.