GEOS
Public Types | Static Public Member Functions | List of all members
geos::AcousticTimeSchemeSEM Struct Reference

Public Types

using EXEC_POLICY = parallelDevicePolicy< >
 

Static Public Member Functions

static void LeapFrogWithoutPML (real64 const dt, arrayView1d< real32 > const p_np1, arrayView1d< real32 > const p_n, arrayView1d< real32 > const p_nm1, arrayView1d< real32 const > const mass, arrayView1d< real32 > const stiffnessVector, arrayView1d< real32 const > const damping, arrayView1d< real32 > const rhs, arrayView1d< localIndex const > const freeSurfaceNodeIndicator, SortedArrayView< localIndex const > const solverTargetNodesSet)
 Apply second order Leap-Frog time scheme for isotropic case without PML. More...
 
static void LeapFrogforVTI (localIndex const size, real64 const dt, arrayView1d< real32 > const p_np1, arrayView1d< real32 > const p_n, arrayView1d< real32 > const p_nm1, arrayView1d< real32 > const q_np1, arrayView1d< real32 > const q_n, arrayView1d< real32 > const q_nm1, arrayView1d< real32 const > const mass, arrayView1d< real32 > const stiffnessVector_p, arrayView1d< real32 > const stiffnessVector_q, arrayView1d< real32 const > const damping_p, arrayView1d< real32 const > const damping_pq, arrayView1d< real32 const > const damping_q, arrayView1d< real32 const > const damping_qp, arrayView1d< real32 > const rhs, arrayView1d< localIndex const > const freeSurfaceNodeIndicator, arrayView1d< localIndex const > const lateralSurfaceNodeIndicator, arrayView1d< localIndex const > const bottomSurfaceNodeIndicator)
 Apply second order Leap-Frog time scheme for VTI case without PML. More...
 

Detailed Description

Definition at line 26 of file AcousticTimeSchemeSEMKernel.hpp.

Member Function Documentation

◆ LeapFrogforVTI()

static void geos::AcousticTimeSchemeSEM::LeapFrogforVTI ( localIndex const  size,
real64 const  dt,
arrayView1d< real32 > const  p_np1,
arrayView1d< real32 > const  p_n,
arrayView1d< real32 > const  p_nm1,
arrayView1d< real32 > const  q_np1,
arrayView1d< real32 > const  q_n,
arrayView1d< real32 > const  q_nm1,
arrayView1d< real32 const > const  mass,
arrayView1d< real32 > const  stiffnessVector_p,
arrayView1d< real32 > const  stiffnessVector_q,
arrayView1d< real32 const > const  damping_p,
arrayView1d< real32 const > const  damping_pq,
arrayView1d< real32 const > const  damping_q,
arrayView1d< real32 const > const  damping_qp,
arrayView1d< real32 > const  rhs,
arrayView1d< localIndex const > const  freeSurfaceNodeIndicator,
arrayView1d< localIndex const > const  lateralSurfaceNodeIndicator,
arrayView1d< localIndex const > const  bottomSurfaceNodeIndicator 
)
inlinestatic

Apply second order Leap-Frog time scheme for VTI case without PML.

Parameters
[in]sizeThe number of nodes in the nodeManager
[in]dttime-step
[out]p_np1pressure array at time n+1 (updated here)
[in]p_npressure array at time n
[in]p_nm1pressure array at time n-1
[out]q_np1auxiliary pressure array at time n+1 (updated here)
[in]q_nauxiliary pressure array at time n
[in]q_nm1auxiliary pressure array at time n-1
[in]massthe mass matrix
[in]stiffnessVector_parray containing the product of the stiffness matrix R and the pressure at time n
[in]stiffnessVector_qarray containing the product of the stiffness matrix R and the auxiliary pressure at time n
[in]damping_pthe damping matrix
[in]damping_pqthe damping matrix
[in]damping_qthe damping matrix
[in]damping_qpthe damping matrix
[in]rhsthe right-hand-side
[in]freeSurfaceNodeIndicatorarray which contains indicators to tell if we are on a free-surface boundary or not
[in]lateralSurfaceNodeIndicatorarray which contains indicators to tell if we are on a lateral boundary or not
[in]bottomSurfaceNodeIndicatorarray which contains indicators to telle if we are on the bottom boundary or not

Definition at line 94 of file AcousticTimeSchemeSEMKernel.hpp.

◆ LeapFrogWithoutPML()

static void geos::AcousticTimeSchemeSEM::LeapFrogWithoutPML ( real64 const  dt,
arrayView1d< real32 > const  p_np1,
arrayView1d< real32 > const  p_n,
arrayView1d< real32 > const  p_nm1,
arrayView1d< real32 const > const  mass,
arrayView1d< real32 > const  stiffnessVector,
arrayView1d< real32 const > const  damping,
arrayView1d< real32 > const  rhs,
arrayView1d< localIndex const > const  freeSurfaceNodeIndicator,
SortedArrayView< localIndex const > const  solverTargetNodesSet 
)
inlinestatic

Apply second order Leap-Frog time scheme for isotropic case without PML.

Parameters
[in]dttime-step
[out]p_np1pressure array at time n+1 (updated here)
[in]p_npressure array at time n
[in]p_nm1pressure array at time n-1
[in]massthe mass matrix
[in]stiffnessVectorarray containing the product of the stiffness matrix R and the pressure at time n
[in]dampingthe damping matrix
[in]rhsthe right-hand-side
[in]freeSurfaceNodeIndicatorarray which contains indicators to tell if we are on a free-surface boundary or not
[in]solverTargetNodesSetthe targetted nodeset (useful in particular when we do elasto-acoustic simulation )

Definition at line 45 of file AcousticTimeSchemeSEMKernel.hpp.


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