feat: SinglePhaseReactiveTransport implementation

inputFiles/singlePhaseFlow/ReactiveCompressible_1d.xml

@@ -0,0 +1,283 @@
+<?xml version="1.0" ?>
+
+<Problem>
+  <Solvers
+    gravityVector="{ 0.0, 0.0, 0.0 }">
+    <SinglePhaseReactiveTransport
+      name="SinglePhaseReactiveFlow"
+      discretization="singlePhaseTPFA"
+      targetRegions="{ Region1 }"
+      writeLinearSystem="0"
+      logLevel="1">
+      <NonlinearSolverParameters
+        newtonTol="1.0e-6"
+        lineSearchAction="None"
+        maxTimeStepCuts="1"
+        newtonMaxIter="8"
+        maxAllowedResidualNorm="1e40"/>
+      <LinearSolverParameters
+        directParallel="0"/>
+    </SinglePhaseReactiveTransport>
+  </Solvers>
+
+  <Mesh>
+    <InternalMesh
+      name="mesh1"
+      elementTypes="{ C3D8 }"
+      xCoords="{ 0, 10 }"
+      yCoords="{ 0, 1 }"
+      zCoords="{ 0, 1 }"
+      nx="{ 10 }"
+      ny="{ 1 }"
+      nz="{ 1 }"
+      cellBlockNames="{ block1 }"/>
+  </Mesh>
+
+  <Geometry>
+    <Box
+      name="source"
+      xMin="{ -0.01, -0.01, -0.01 }"
+      xMax="{ 1.01, 1.01, 1.01 }"/>
+
+    <Box
+      name="sink"
+      xMin="{ 8.99, -0.01, -0.01 }"
+      xMax="{ 10.01, 1.01, 1.01 }"/>
+  </Geometry>
+
+  <Events
+    maxTime="500">
+    <PeriodicEvent
+      name="outputs"
+      cycleFrequency="1"
+      target="/Outputs/vtkOutput"/>
+
+    <PeriodicEvent
+      name="solverApplications"
+      forceDt="5"
+      target="/Solvers/SinglePhaseReactiveFlow"/>
+  </Events>
+
+  <NumericalMethods>
+    <FiniteVolume>
+      <TwoPointFluxApproximation
+        name="singlePhaseTPFA"/>
+    </FiniteVolume>
+  </NumericalMethods>
+
+  <ElementRegions>
+    <CellElementRegion
+      name="Region1"
+      cellBlocks="{ * }"
+      materialList="{ water, rock }"/>
+  </ElementRegions>
+
+  <Constitutive>
+    <ReactiveCompressibleSinglePhase
+      name="water"
+      defaultDensity="1000"
+      defaultViscosity="0.001"
+      referencePressure="0.0"
+      compressibility="5e-10"
+      viscosibility="0.0"/>
+
+    <CompressibleSolidConstantPermeability
+      name="rock"
+      solidModelName="nullSolid"
+      porosityModelName="rockPorosity"
+      permeabilityModelName="rockPerm"/>
+
+    <NullModel
+      name="nullSolid"/>
+
+    <PressurePorosity
+      name="rockPorosity"
+      defaultReferencePorosity="0.05"
+      referencePressure="0.0"
+      compressibility="1.0e-9"/>
+
+    <ConstantPermeability
+      name="rockPerm"
+      permeabilityComponents="{ 2.0e-14, 2.0e-14, 2.0e-14 }"/>
+  </Constitutive>
+
+  <FieldSpecifications>
+    <FieldSpecification
+      name="Porosity"
+      initialCondition="1"
+      setNames="{ all }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="rockPorosity_referencePorosity"
+      scale="0.05"/>
+
+    <FieldSpecification
+      name="initialPressure"
+      initialCondition="1"
+      setNames="{ all }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="pressure"
+      scale="1e6"/>
+
+    <FieldSpecification
+      name="initialLogSpecies0"
+      initialCondition="1"
+      setNames="{ all }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="0"
+      scale="-8"/>
+
+    <FieldSpecification
+      name="initialLogSpecies1"
+      initialCondition="1"
+      setNames="{ all }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="1"
+      scale="-10"/>
+
+    <FieldSpecification
+      name="initialLogSpecies2"
+      initialCondition="1"
+      setNames="{ all }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="2"
+      scale="-9"/>
+<!-- 
+    <FieldSpecification
+      name="initialLogSpecies3"
+      initialCondition="1"
+      setNames="{ all }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="3"
+      scale="-10"/>
+
+    <FieldSpecification
+      name="initialLogSpecies4"
+      initialCondition="1"
+      setNames="{ all }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="4"
+      scale="-10"/>
+
+    <FieldSpecification
+      name="initialLogSpecies5"
+      initialCondition="1"
+      setNames="{ all }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="5"
+      scale="-10"/>
+
+    <FieldSpecification
+      name="initialLogSpecies6"
+      initialCondition="1"
+      setNames="{ all }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="6"
+      scale="-10"/> -->
+
+    <FieldSpecification
+      name="sourceLogSpecies0"
+      setNames="{ source }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="0"
+      scale="-6"/>
+
+    <FieldSpecification
+      name="sourceLogSpecies1"
+      setNames="{ source }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="1"
+      scale="-8"/>
+
+    <FieldSpecification
+      name="sourceLogSpecies2"
+      setNames="{ source }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="2"
+      scale="-7.2"/>
+
+    <!-- <FieldSpecification
+      name="sourceLogSpecies3"
+      setNames="{ source }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="3"
+      scale="-10"/>
+
+    <FieldSpecification
+      name="sourceLogSpecies4"
+      setNames="{ source }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="4"
+      scale="-10"/>
+
+    <FieldSpecification
+      name="sourceLogSpecies5"
+      setNames="{ source }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="5"
+      scale="-10"/>
+
+    <FieldSpecification
+      name="sourceLogSpecies6"
+      setNames="{ source }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="6"
+      scale="-10"/>  -->
+
+    <FieldSpecification
+      name="sinkLogSpecies0"
+      setNames="{ sink }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="0"
+      scale="-8"/>
+
+    <FieldSpecification
+      name="sinkLogSpecies1"
+      setNames="{ sink }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="1"
+      scale="-10"/>
+
+    <FieldSpecification
+      name="sinkLogSpecies2"
+      setNames="{ sink }"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="logPrimarySpeciesConcentration"
+      component="2"
+      scale="-9"/>
+
+    <FieldSpecification
+      name="sourceTerm"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="pressure"
+      scale="5e6"
+      setNames="{ source }"/>
+
+    <FieldSpecification
+      name="sinkTerm"
+      objectPath="ElementRegions/Region1/block1"
+      fieldName="pressure"
+      scale="1e6"
+      setNames="{ sink }"/>
+  </FieldSpecifications>
+
+  <Outputs>
+    <VTK
+      name="vtkOutput"/>
+  </Outputs>
+</Problem>

src/coreComponents/constitutive/CMakeLists.txt

@@ -135,6 +135,9 @@ set( constitutive_headers
      fluid/singlefluid/SingleFluidSelector.hpp
      fluid/singlefluid/SlurryFluidSelector.hpp
      fluid/singlefluid/ThermalCompressibleSinglePhaseFluid.hpp
+     fluid/singlefluid/reactive/ReactiveCompressibleSinglePhaseFluid.hpp
+     fluid/singlefluid/reactive/ReactiveSingleFluid.hpp
+     fluid/singlefluid/reactive/ThermalReactiveCompressibleSinglePhaseFluid.hpp
      permeability/CarmanKozenyPermeability.hpp
      permeability/ConstantPermeability.hpp
      permeability/DamagePermeability.hpp
@@ -279,6 +282,9 @@ set( constitutive_sources
      fluid/singlefluid/SingleFluidBase.cpp
      fluid/singlefluid/SlurryFluidBase.cpp
      fluid/singlefluid/ThermalCompressibleSinglePhaseFluid.cpp
+     fluid/singlefluid/reactive/ReactiveCompressibleSinglePhaseFluid.cpp
+     fluid/singlefluid/reactive/ReactiveSingleFluid.cpp
+     fluid/singlefluid/reactive/ThermalReactiveCompressibleSinglePhaseFluid.cpp
      permeability/CarmanKozenyPermeability.cpp
      permeability/ConstantPermeability.cpp
      permeability/DamagePermeability.cpp

src/coreComponents/constitutive/fluid/multifluid/reactive/ReactiveFluidSelector.hpp

@@ -21,6 +21,8 @@
 
 #include "constitutive/ConstitutivePassThruHandler.hpp"
 #include "constitutive/fluid/multifluid/reactive/ReactiveBrineFluid.hpp"
+#include "constitutive/fluid/singlefluid/reactive/ReactiveCompressibleSinglePhaseFluid.hpp"
+#include "constitutive/fluid/singlefluid/reactive/ThermalReactiveCompressibleSinglePhaseFluid.hpp"
 
 #include "common/GeosxConfig.hpp"
 
@@ -46,6 +48,22 @@ void constitutiveUpdatePassThru( ReactiveMultiFluid & fluid,
                                ReactiveBrineThermal >::execute( fluid, std::forward< LAMBDA >( lambda ) );
 }
 
+template< typename LAMBDA >
+void constitutiveUpdatePassThru( ReactiveSingleFluid const & fluid,
+                                 LAMBDA && lambda )
+{
+  ConstitutivePassThruHandler< ReactiveCompressibleSinglePhase,
+                               ThermalReactiveCompressibleSinglePhase >::execute( fluid, std::forward< LAMBDA >( lambda ) );
+}
+
+template< typename LAMBDA >
+void constitutiveUpdatePassThru( ReactiveSingleFluid & fluid,
+                                 LAMBDA && lambda )
+{
+  ConstitutivePassThruHandler< ReactiveCompressibleSinglePhase,
+                               ThermalReactiveCompressibleSinglePhase >::execute( fluid, std::forward< LAMBDA >( lambda ) );
+}
+
 } // namespace constitutive
 
 } // namespace geos

src/coreComponents/constitutive/fluid/multifluid/reactive/ReactiveMultiFluid.cpp

@@ -39,7 +39,7 @@ ReactiveMultiFluid::
 
   registerField( fields::reactivefluid::primarySpeciesConcentration{}, &m_primarySpeciesConcentration );
   registerField( fields::reactivefluid::secondarySpeciesConcentration{}, &m_secondarySpeciesConcentration );
-  registerField( fields::reactivefluid::primarySpeciesTotalConcentration{}, &m_primarySpeciesTotalConcentration );
+  registerField( fields::reactivefluid::primarySpeciesAggregateConcentration{}, &m_primarySpeciesTotalConcentration );
   registerField( fields::reactivefluid::kineticReactionRates{}, &m_kineticReactionRates );
 }
 

src/coreComponents/constitutive/fluid/multifluid/reactive/ReactiveMultiFluidFields.hpp

@@ -33,6 +33,7 @@ namespace reactivefluid
 {
 
 using array2dLayoutComp = array2d< real64, compflow::LAYOUT_COMP >;
+using array3dLayoutComp_dC = array3d< real64, compflow::LAYOUT_COMP_DC >;
 
 DECLARE_FIELD( primarySpeciesConcentration,
                "primarySpeciesConcentration",
@@ -42,13 +43,29 @@ DECLARE_FIELD( primarySpeciesConcentration,
                WRITE_AND_READ,
                "primarySpeciesConcentration" );
 
-DECLARE_FIELD( primarySpeciesTotalConcentration,
-               "primarySpeciesTotalConcentration",
+DECLARE_FIELD( primarySpeciesAggregateConcentration,
+               "primarySpeciesAggregateConcentration",
                array2dLayoutComp,
                0,
                LEVEL_0,
                WRITE_AND_READ,
-               "primarySpeciesTotalConcentration" );
+               "primarySpeciesAggregateConcentration" );
+
+DECLARE_FIELD( primarySpeciesAggregateConcentration_n,
+               "primarySpeciesAggregateConcentration_n",
+               array2dLayoutComp,
+               0,
+               LEVEL_0,
+               WRITE_AND_READ,
+               "primarySpeciesAggregateConcentration at the previous timestep" );
+
+DECLARE_FIELD( dPrimarySpeciesAggregateConcentration_dLogPrimaryConc,
+               "dPrimarySpeciesAggregateConcentration_dLogPrimaryConc",
+               array3dLayoutComp_dC,
+               0,
+               LEVEL_0,
+               WRITE_AND_READ,
+               "Deivatives of primarySpeciesAggregateConcentration w.r.t log primary species concentration" );
 
 DECLARE_FIELD( secondarySpeciesConcentration,
                "secondarySpeciesConcentration",