[WIP] acoustic - elastic coupling (SEM, 2ⁿᵈ order wave equation)

inputFiles/wavePropagation/acous3D_Q3_small_base.xml

@@ -67,6 +67,14 @@
       scale="1500"
       setNames="{ all }"/>
 
+    <FieldSpecification
+      name="cellDensity"
+      initialCondition="1"
+      objectPath="mesh/FE2/ElementRegions/Region/cb"
+      fieldName="mediumDensityA"
+      scale="1"
+      setNames="{ all }"/>
+
     <!-- 3) The free surface condition the domain -->
     <FieldSpecification
       name="zposFreeSurface"

inputFiles/wavePropagation/acous3D_Q5_small_base.xml

@@ -67,6 +67,14 @@
       scale="1500"
       setNames="{ all }"/>
 
+    <FieldSpecification
+      name="cellDensity"
+      initialCondition="1"
+      objectPath="mesh/FE2/ElementRegions/Region/cb"
+      fieldName="mediumDensityA"
+      scale="1500"
+      setNames="{ all }"/>
+
     <!-- 3) The free surface condition the domain -->
     <FieldSpecification
       name="zposFreeSurface"

inputFiles/wavePropagation/acous3D_pml_smoke.xml

@@ -82,7 +82,7 @@
       name="solverApplications"
       forceDt="0.0025"
       target="/Solvers/acousticSolver"/>
-    
+
     <!-- two events to output pressure in an hdf5 file -->
     <PeriodicEvent
       name="timeHistoryCollection"
@@ -99,7 +99,7 @@
       name="restarts"
       timeFrequency="0.05"
       target="/Outputs/restartOutput"/>
-    
+
   </Events>
 
   <NumericalMethods>
@@ -168,7 +168,7 @@
                   right,
                   front,
                   back}"/>
-    
+
     <!-- The PML field specification assigned to the corresponding region -->
     <PML
       name="pml"
@@ -185,6 +185,15 @@
       scale="2000"
       setNames="{ all }"/>
 
+    <FieldSpecification
+      name="cellDensity1"
+      initialCondition="1"
+      objectPath="ElementRegions/interiorDomain"
+      fieldName="mediumDensityA"
+      scale="1"
+      setNames="{ all }"/>
+
+
     <FieldSpecification
       name="cellVelocity2"
       initialCondition="1"
@@ -193,8 +202,16 @@
       scale="2000"
       setNames="{ all }"/>
 
+    <FieldSpecification
+      name="cellDensity2"
+      initialCondition="1"
+      objectPath="ElementRegions/pmlDomain"
+      fieldName="mediumDensityA"
+      scale="1"
+      setNames="{ all }"/>
+
   </FieldSpecifications>
-  
+
   <!-- collect the pressure values at the nodes -->
   <Tasks>
     <PackCollection
@@ -214,7 +231,7 @@
 
     <Restart
       name="restartOutput"/>
-    
+
   </Outputs>
 
 </Problem>

inputFiles/wavePropagation/acous3D_small_base.xml

@@ -67,6 +67,14 @@
       scale="1500"
       setNames="{ all }"/>
 
+    <FieldSpecification
+      name="cellDensity"
+      initialCondition="1"
+      objectPath="ElementRegions/Region/cb"
+      fieldName="mediumDensityA"
+      scale="1"
+      setNames="{ all }"/>
+
     <!-- 3) The free surface condition the domain -->
     <FieldSpecification
       name="zposFreeSurface"

inputFiles/wavePropagation/acouselas3D.xml

@@ -0,0 +1,230 @@
+<?xml version="1.0" ?>
+
+<Problem>
+  <Solvers>
+    <!-- define the solver, source coordinates, time source frequency, receiver coordinates -->
+    <AcousticSEM
+      name="acousticSolver"
+      discretization="FE1"
+      targetRegions="{ Fluid }"
+      sourceCoordinates="{ { 1575, 2400, 2900 } }"
+      timeSourceFrequency="10"
+      timeSourceDelay="0.115"
+      rickerOrder="3"
+      outputSeismoTrace="1"
+      dtSeismoTrace="1e-3"
+      receiverCoordinates="{ { 3752.2936, 2400, 2933.3333 } }"/>
+    <ElasticSEM
+      name="elasticSolver"
+      discretization="FE1"
+      targetRegions="{ Solid }"
+      sourceCoordinates="{ }"
+      outputSeismoTrace="1"
+      dtSeismoTrace="1e-3"
+      receiverCoordinates="{ { 3752.2936, 2400, 1866.6667 } }"/>
+    <!-- NOTES:
+      - coupled solver must be declared at the end (for interface computations)
+      - discretization is also needed for the correct mesh to be passed through from the coupled solver to the elastic / acoustic solvers
+    -->
+    <AcousticElasticSEM
+      name="acousticelasticSolver"
+      discretization="FE1"
+      targetRegions="{ Solid, Fluid }"
+      acousticSolverName="acousticSolver"
+      elasticSolverName="elasticSolver"/>
+  </Solvers>
+
+  <!--
+  hexahedral mesh generated internally by GEOSX
+  6400x1600x4800m domain
+  # R1-Q5 (pseudo-2D)
+  sourceCoordinates="{ { 1500, 20, 2900  } }"
+  yCoords="{ 0, 40 }"
+  nx="{ 160 }"
+  ny="{ 1 }"
+  nz="{ 60, 60 }"
+  # R2-Q2 (pseudo-2D)
+  sourceCoordinates="{ { 1500, 10, 2900  } }"
+  yCoords="{ 0, 20 }"
+  nx="{ 320 }"
+  ny="{ 1 }"
+  nz="{ 120, 120 }"
+  # R5-Q1 (pseudo-2D)
+  sourceCoordinates="{ { 1500, 4, 2900  } }"
+  yCoords="{ 0, 8 }"
+  nx="{ 800 }"
+  ny="{ 1 }"
+  nz="{ 300, 300 }"
+  NOTES:
+    nx="{ 320 }"
+    ny="{ 240 }"
+    nz="{ 120, 120 }
+    => local threshold for local builds 3D Q1 Δt=1e-4 (18e6 elems, 10Gb per 1 GPU)
+    => Q4 180 procs / GPU 16h
+    => Q1 180 procs / GPU 1e40min
+    nx="{ 160 }"
+    ny="{ 120 }"
+    nz="{ 60, 60 }
+    => local threshold for local builds 3D Q2 Δt=1e-4 (2e6 elems, 4Gb per 1 GPU)
+    nx="{ 80 }"
+    ny="{ 60 }"
+    nz="{ 30, 30 }
+    => coarse 10 min frontal run 3D Q5 Δt=1e-3 (288e3 elems, 8Gb per 6 GPU)
+  -->
+  <Mesh>
+    <InternalMesh
+      name="mesh"
+      elementTypes="{ C3D8 }"
+      xCoords="{ 0, 6400 }"
+      yCoords="{ 0, 4800.0001 }"
+      zCoords="{ 0, 2400, 4800 }"
+      nx="{ 320 }"
+      ny="{ 240 }"
+      nz="{ 120, 120 }"
+      cellBlockNames="{ solid, fluid }"/>
+  </Mesh>
+
+  <NumericalMethods>
+    <FiniteElements>
+      <FiniteElementSpace
+        name="FE1"
+        order="3"
+        formulation="SEM"/>
+    </FiniteElements>
+  </NumericalMethods>
+
+  <Events
+    maxTime="3">
+    <PeriodicEvent
+      name="outputs"
+      timeFrequency="5e-2"
+      targetExactTimestep="0">
+      <PeriodicEvent name="outputs_vtk" target="/Outputs/vtkOutput"/>
+      <!--
+      <PeriodicEvent name="tsk_ux_np1" target="/Tasks/displacementxCollection"/>
+      <PeriodicEvent name="out_ux_np1" target="/Outputs/timeHistoryOutputx"/>
+      <PeriodicEvent name="tsk_uy_np1" target="/Tasks/displacementyCollection"/>
+      <PeriodicEvent name="out_uy_np1" target="/Outputs/timeHistoryOutputy"/>
+      <PeriodicEvent name="tsk_uz_np1" target="/Tasks/displacementzCollection"/>
+      <PeriodicEvent name="out_uz_np1" target="/Outputs/timeHistoryOutputz"/>
+      <PeriodicEvent name="tsk_p_np1" target="/Tasks/pressureCollection"/>
+      <PeriodicEvent name="out_p_np1" target="/Outputs/timeHistoryOutput"/>
+      -->
+    </PeriodicEvent>
+    <!-- control the timestepping here with forceDt -->
+    <PeriodicEvent
+      name="solidfluidSolverApplications"
+      forceDt="1e-4"
+      target="/Solvers/acousticelasticSolver"/>
+    <!-- FIXME: forceDt hack to satisfy pre-timeloop computation -->
+    <SoloEvent
+      name="solidSolverApplications"
+      forceDt="1e-4"
+      endTime="0"
+      target="/Solvers/elasticSolver"/>
+    <SoloEvent
+      name="fluidSolverApplications"
+      forceDt="1e-4"
+      endTime="0"
+      target="/Solvers/acousticSolver"/>
+  </Events>
+
+  <ElementRegions>
+    <CellElementRegion
+      name="Fluid"
+      cellBlocks="{ fluid }"
+      materialList="{ nullModel }"/>
+    <CellElementRegion
+      name="Solid"
+      cellBlocks="{ solid }"
+      materialList="{ nullModel }"/>
+  </ElementRegions>
+
+  <Constitutive>
+    <NullModel
+      name="nullModel"/>
+  </Constitutive>
+
+  <FieldSpecifications>
+    <!-- initial condition for ρ, 𝓋ₚ and 𝓋ₛ -->
+    <FieldSpecification
+      name="acousticVelocity"
+      initialCondition="1"
+      objectPath="mesh/FE1/ElementRegions/Fluid/fluid"
+      fieldName="mediumVelocity"
+      scale="1500.0"
+      setNames="{ all }"/>
+    <FieldSpecification
+      name="acousticDensity"
+      initialCondition="1"
+      objectPath="mesh/FE1/ElementRegions/Fluid/fluid"
+      fieldName="mediumDensityA"
+      scale="1020.0"
+      setNames="{ all }"/>
+    <FieldSpecification
+      name="elasticVelocityVp"
+      initialCondition="1"
+      objectPath="mesh/FE1/ElementRegions/Solid/solid"
+      fieldName="mediumVelocityVp"
+      scale="3400.0"
+      setNames="{ all }"/>
+    <FieldSpecification
+      name="elasticVelocityVs"
+      initialCondition="1"
+      objectPath="mesh/FE1/ElementRegions/Solid/solid"
+      fieldName="mediumVelocityVs"
+      scale="1963.0"
+      setNames="{ all }"/>
+    <FieldSpecification
+      name="elasticDensity"
+      initialCondition="1"
+      objectPath="mesh/FE1/ElementRegions/Solid/solid"
+      fieldName="mediumDensityE"
+      scale="2500.0"
+      setNames="{ all }"/>
+  </FieldSpecifications>
+
+  <Tasks>
+    <PackCollection
+      name="pressureCollection"
+      objectPath="mesh/FE1/nodeManager"
+      fieldName="pressure_np1"/>
+    <PackCollection
+      name="displacementxCollection"
+      objectPath="mesh/FE1/nodeManager"
+      fieldName="displacementx_np1"/>
+    <PackCollection
+      disableCoordCollection="1"
+      name="displacementyCollection"
+      objectPath="mesh/FE1/nodeManager"
+      fieldName="displacementy_np1"/>
+    <PackCollection
+      disableCoordCollection="1"
+      name="displacementzCollection"
+      objectPath="mesh/FE1/nodeManager"
+      fieldName="displacementz_np1"/>
+  </Tasks>
+
+  <Outputs>
+    <VTK
+      name="vtkOutput"
+      levelNames="{ FE1 }"
+      plotLevel="0"/>
+    <TimeHistory
+      name="timeHistoryOutput"
+      sources="{ /Tasks/pressureCollection }"
+      filename="pressure_history"/>
+    <TimeHistory
+      name="timeHistoryOutputx"
+      sources="{ /Tasks/displacementxCollection }"
+      filename="displacementx_history"/>
+    <TimeHistory
+      name="timeHistoryOutputy"
+      sources="{ /Tasks/displacementyCollection }"
+      filename="displacementy_history"/>
+    <TimeHistory
+      name="timeHistoryOutputz"
+      sources="{ /Tasks/displacementzCollection }"
+      filename="displacementz_history"/>
+  </Outputs>
+</Problem>

inputFiles/wavePropagation/benchmarks/elas3D_benchmark_base.xml

@@ -184,7 +184,7 @@
       name="cellDensity"
       initialCondition="1"
       objectPath="ElementRegions/Region/cb"
-      fieldName="mediumDensity"
+      fieldName="mediumDensityE"
       scale="1"
       setNames="{ all }"/>
 

inputFiles/wavePropagation/elas3D_DAS_smoke.xml

@@ -149,7 +149,7 @@
       name="cellDensity"
       initialCondition="1"
       objectPath="ElementRegions/Region/cb"
-      fieldName="mediumDensity"
+      fieldName="mediumDensityE"
       scale="2000"
       setNames="{ all }"/>
 

inputFiles/wavePropagation/elas3D_Q3_small_base.xml

@@ -163,7 +163,7 @@
       name="cellDensity"
       initialCondition="1"
       objectPath="mesh/FE1/ElementRegions/Region/cb"
-      fieldName="mediumDensity"
+      fieldName="mediumDensityE"
       scale="1"
       setNames="{ all }"/>
 

inputFiles/wavePropagation/elas3D_Q5_small_base.xml

@@ -163,7 +163,7 @@
       name="cellDensity"
       initialCondition="1"
       objectPath="mesh/FE1/ElementRegions/Region/cb"
-      fieldName="mediumDensity"
+      fieldName="mediumDensityA"
       scale="1"
       setNames="{ all }"/>
 

inputFiles/wavePropagation/elas3D_small_base.xml

@@ -170,7 +170,7 @@
       name="cellDensity"
       initialCondition="1"
       objectPath="ElementRegions/Region/cb"
-      fieldName="mediumDensity"
+      fieldName="mediumDensityE"
       scale="1"
       setNames="{ all }"/>