Skip to content

Commit

Permalink
refactor: clean up kernel flags usage in flow solver to avoid bugs
Browse files Browse the repository at this point in the history
  • Loading branch information
Pavel Tomin authored and Pavel Tomin committed Dec 19, 2024
1 parent 93f0252 commit f562029
Show file tree
Hide file tree
Showing 17 changed files with 103 additions and 159 deletions.
14 changes: 11 additions & 3 deletions src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseBase.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1349,6 +1349,15 @@ void CompositionalMultiphaseBase::assembleAccumulationAndVolumeBalanceTerms( Dom
{
GEOS_MARK_FUNCTION;

using namespace isothermalCompositionalMultiphaseBaseKernels;

BitFlags< KernelFlags > kernelFlags;
if( m_useTotalMassEquation )
kernelFlags.set( KernelFlags::TotalMassEquation );
if( m_useSimpleAccumulation )
kernelFlags.set( KernelFlags::SimpleAccumulation );


forDiscretizationOnMeshTargets( domain.getMeshBodies(), [&]( string const &,
MeshLevel const & mesh,
arrayView1d< string const > const & regionNames )
Expand All @@ -1371,7 +1380,7 @@ void CompositionalMultiphaseBase::assembleAccumulationAndVolumeBalanceTerms( Dom
createAndLaunch< parallelDevicePolicy<> >( m_numComponents,
m_numPhases,
dofManager.rankOffset(),
m_useTotalMassEquation,
kernelFlags,
dofKey,
subRegion,
fluid,
Expand All @@ -1386,8 +1395,7 @@ void CompositionalMultiphaseBase::assembleAccumulationAndVolumeBalanceTerms( Dom
createAndLaunch< parallelDevicePolicy<> >( m_numComponents,
m_numPhases,
dofManager.rankOffset(),
m_useTotalMassEquation,
m_useSimpleAccumulation,
kernelFlags,
dofKey,
subRegion,
fluid,
Expand Down
59 changes: 42 additions & 17 deletions src/coreComponents/physicsSolvers/fluidFlow/CompositionalMultiphaseFVM.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -164,6 +164,20 @@ void CompositionalMultiphaseFVM::assembleFluxTerms( real64 const dt,
{
GEOS_MARK_FUNCTION;

using namespace isothermalCompositionalMultiphaseFVMKernels;

BitFlags< KernelFlags > kernelFlags;
if( m_hasCapPressure )
kernelFlags.set( KernelFlags::CapPressure );
if( m_hasDiffusion )
kernelFlags.set( KernelFlags::Diffusion );
if( m_hasDispersion )
kernelFlags.set( KernelFlags::Dispersion );
if( m_useTotalMassEquation )
kernelFlags.set( KernelFlags::TotalMassEquation );
if( m_useNewGravity )
kernelFlags.set( KernelFlags::NewGravity );

forDiscretizationOnMeshTargets( domain.getMeshBodies(), [&]( string const &,
MeshLevel const & mesh,
arrayView1d< string const > const & )
Expand All @@ -172,6 +186,13 @@ void CompositionalMultiphaseFVM::assembleFluxTerms( real64 const dt,
FiniteVolumeManager const & fvManager = numericalMethodManager.getFiniteVolumeManager();
FluxApproximationBase const & fluxApprox = fvManager.getFluxApproximation( m_discretizationName );

auto const & upwindingParams = fluxApprox.upwindingParams();
if( upwindingParams.upwindingScheme == UpwindingScheme::C1PPU &&
isothermalCompositionalMultiphaseFVMKernelUtilities::epsC1PPU > 0 )
kernelFlags.set( KernelFlags::C1PPU );
else if( upwindingParams.upwindingScheme == UpwindingScheme::IHU )
kernelFlags.set( KernelFlags::IHU );

string const & elemDofKey = dofManager.getKey( viewKeyStruct::elemDofFieldString() );

fluxApprox.forAllStencils( mesh, [&] ( auto & stencil )
Expand All @@ -187,8 +208,7 @@ void CompositionalMultiphaseFVM::assembleFluxTerms( real64 const dt,
m_numPhases,
dofManager.rankOffset(),
elemDofKey,
m_hasCapPressure,
m_useTotalMassEquation,
kernelFlags,
getName(),
mesh.getElemManager(),
stencilWrapper,
Expand All @@ -206,8 +226,7 @@ void CompositionalMultiphaseFVM::assembleFluxTerms( real64 const dt,
m_numPhases,
dofManager.rankOffset(),
elemDofKey,
m_hasCapPressure,
m_useTotalMassEquation,
kernelFlags,
getName(),
mesh.getElemManager(),
stencilWrapper,
Expand All @@ -229,10 +248,7 @@ void CompositionalMultiphaseFVM::assembleFluxTerms( real64 const dt,
m_numPhases,
dofManager.rankOffset(),
elemDofKey,
m_hasCapPressure,
m_useTotalMassEquation,
m_useNewGravity,
fluxApprox.upwindingParams(),
kernelFlags,
getName(),
mesh.getElemManager(),
stencilWrapper,
Expand All @@ -254,9 +270,7 @@ void CompositionalMultiphaseFVM::assembleFluxTerms( real64 const dt,
m_numPhases,
dofManager.rankOffset(),
elemDofKey,
m_hasDiffusion,
m_hasDispersion,
m_useTotalMassEquation,
kernelFlags,
getName(),
mesh.getElemManager(),
stencilWrapper,
Expand All @@ -272,9 +286,7 @@ void CompositionalMultiphaseFVM::assembleFluxTerms( real64 const dt,
m_numPhases,
dofManager.rankOffset(),
elemDofKey,
m_hasDiffusion,
m_hasDispersion,
m_useTotalMassEquation,
kernelFlags,
getName(),
mesh.getElemManager(),
stencilWrapper,
Expand All @@ -296,6 +308,14 @@ void CompositionalMultiphaseFVM::assembleStabilizedFluxTerms( real64 const dt,
{
GEOS_MARK_FUNCTION;

using namespace isothermalCompositionalMultiphaseFVMKernels;

BitFlags< KernelFlags > kernelFlags;
if( m_hasCapPressure )
kernelFlags.set( KernelFlags::CapPressure );
if( m_useTotalMassEquation )
kernelFlags.set( KernelFlags::TotalMassEquation );

forDiscretizationOnMeshTargets( domain.getMeshBodies(), [&]( string const &,
MeshLevel const & mesh,
arrayView1d< string const > const & )
Expand All @@ -321,8 +341,7 @@ void CompositionalMultiphaseFVM::assembleStabilizedFluxTerms( real64 const dt,
m_numPhases,
dofManager.rankOffset(),
elemDofKey,
m_hasCapPressure,
m_useTotalMassEquation,
kernelFlags,
getName(),
mesh.getElemManager(),
stencilWrapper,
Expand Down Expand Up @@ -1008,6 +1027,12 @@ void CompositionalMultiphaseFVM::applyFaceDirichletBC( real64 const time_n,
GEOS_ERROR_IF( !bcConsistent, GEOS_FMT( "CompositionalMultiphaseBase {}: inconsistent boundary conditions", getDataContext() ) );
}

using namespace isothermalCompositionalMultiphaseFVMKernels;

BitFlags< KernelFlags > kernelFlags;
if( m_useTotalMassEquation )
kernelFlags.set( KernelFlags::TotalMassEquation );

FieldSpecificationManager & fsManager = FieldSpecificationManager::getInstance();

NumericalMethodsManager const & numericalMethodManager = domain.getNumericalMethodManager();
Expand Down Expand Up @@ -1070,7 +1095,7 @@ void CompositionalMultiphaseFVM::applyFaceDirichletBC( real64 const time_n,
createAndLaunch< parallelDevicePolicy<> >( m_numComponents,
m_numPhases,
dofManager.rankOffset(),
m_useTotalMassEquation,
kernelFlags,
elemDofKey,
getName(),
faceManager,
Expand Down
9 changes: 1 addition & 8 deletions src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/AccumulationKernel.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -474,8 +474,7 @@ class AccumulationKernelFactory
createAndLaunch( integer const numComps,
integer const numPhases,
globalIndex const rankOffset,
integer const useTotalMassEquation,
integer const useSimpleAccumulation,
BitFlags< KernelFlags > kernelFlags,
string const dofKey,
ElementSubRegionBase const & subRegion,
constitutive::MultiFluidBase const & fluid,
Expand All @@ -488,12 +487,6 @@ class AccumulationKernelFactory
integer constexpr NUM_COMP = NC();
integer constexpr NUM_DOF = NC()+1;

BitFlags< KernelFlags > kernelFlags;
if( useTotalMassEquation )
kernelFlags.set( KernelFlags::TotalMassEquation );
if( useSimpleAccumulation )
kernelFlags.set( KernelFlags::SimpleAccumulation );

AccumulationKernel< NUM_COMP, NUM_DOF > kernel( numPhases, rankOffset, dofKey, subRegion,
fluid, solid, localMatrix, localRhs, kernelFlags );
AccumulationKernel< NUM_COMP, NUM_DOF >::template launch< POLICY >( subRegion.size(), kernel );
Expand Down
11 changes: 3 additions & 8 deletions ...s/physicsSolvers/fluidFlow/kernels/compositional/DiffusionDispersionFluxComputeKernel.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -725,9 +725,7 @@ class DiffusionDispersionFluxComputeKernelFactory
integer const numPhases,
globalIndex const rankOffset,
string const & dofKey,
integer const hasDiffusion,
integer const hasDispersion,
integer const useTotalMassEquation,
BitFlags< KernelFlags > kernelFlags,
string const & solverName,
ElementRegionManager const & elemManager,
STENCILWRAPPER const & stencilWrapper,
Expand All @@ -744,10 +742,6 @@ class DiffusionDispersionFluxComputeKernelFactory
elemManager.constructArrayViewAccessor< globalIndex, 1 >( dofKey );
dofNumberAccessor.setName( solverName + "/accessors/" + dofKey );

BitFlags< KernelFlags > kernelFlags;
if( useTotalMassEquation )
kernelFlags.set( KernelFlags::TotalMassEquation );

using kernelType = DiffusionDispersionFluxComputeKernel< NUM_COMP, NUM_DOF, STENCILWRAPPER >;
typename kernelType::CompFlowAccessors compFlowAccessors( elemManager, solverName );
typename kernelType::MultiFluidAccessors multiFluidAccessors( elemManager, solverName );
Expand All @@ -760,7 +754,8 @@ class DiffusionDispersionFluxComputeKernelFactory
diffusionAccessors, dispersionAccessors, porosityAccessors,
dt, localMatrix, localRhs, kernelFlags );
kernelType::template launch< POLICY >( stencilWrapper.size(),
hasDiffusion, hasDispersion,
kernelFlags.isSet( KernelFlags::Diffusion ),
kernelFlags.isSet( KernelFlags::Dispersion ),
kernel );
} );
}
Expand Down
9 changes: 1 addition & 8 deletions ...omponents/physicsSolvers/fluidFlow/kernels/compositional/DissipationFluxComputeKernel.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -349,8 +349,7 @@ class FluxComputeKernelFactory
integer const numPhases,
globalIndex const rankOffset,
string const & dofKey,
integer const hasCapPressure,
integer const useTotalMassEquation,
BitFlags< isothermalCompositionalMultiphaseFVMKernels::KernelFlags > kernelFlags,
string const & solverName,
ElementRegionManager const & elemManager,
STENCILWRAPPER const & stencilWrapper,
Expand All @@ -374,12 +373,6 @@ class FluxComputeKernelFactory
elemManager.constructArrayViewAccessor< globalIndex, 1 >( dofKey );
dofNumberAccessor.setName( solverName + "/accessors/" + dofKey );

BitFlags< isothermalCompositionalMultiphaseFVMKernels::KernelFlags > kernelFlags;
if( hasCapPressure )
kernelFlags.set( isothermalCompositionalMultiphaseFVMKernels::KernelFlags::CapPressure );
if( useTotalMassEquation )
kernelFlags.set( isothermalCompositionalMultiphaseFVMKernels::KernelFlags::TotalMassEquation );

using KERNEL_TYPE = FluxComputeKernel< NUM_COMP, NUM_DOF, STENCILWRAPPER >;
typename KERNEL_TYPE::CompFlowAccessors compFlowAccessors( elemManager, solverName );
typename KERNEL_TYPE::MultiFluidAccessors multiFluidAccessors( elemManager, solverName );
Expand Down
21 changes: 2 additions & 19 deletions src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/FluxComputeKernel.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -360,7 +360,7 @@ class FluxComputeKernel : public FluxComputeKernelBase
k_up, seri[k_up], sesri[k_up], sei[k_up], potGrad,
phaseFlux, dPhaseFlux_dP, dPhaseFlux_dC );

} // loop over phases
} // loop over phases

/// populate local flux vector and derivatives
for( integer ic = 0; ic < numComp; ++ic )
Expand Down Expand Up @@ -527,10 +527,7 @@ class FluxComputeKernelFactory
integer const numPhases,
globalIndex const rankOffset,
string const & dofKey,
integer const hasCapPressure,
integer const useTotalMassEquation,
integer const useNewGravity,
UpwindingParameters upwindingParams,
BitFlags< KernelFlags > kernelFlags,
string const & solverName,
ElementRegionManager const & elemManager,
STENCILWRAPPER const & stencilWrapper,
Expand All @@ -547,20 +544,6 @@ class FluxComputeKernelFactory
elemManager.constructArrayViewAccessor< globalIndex, 1 >( dofKey );
dofNumberAccessor.setName( solverName + "/accessors/" + dofKey );

BitFlags< KernelFlags > kernelFlags;
if( hasCapPressure )
kernelFlags.set( KernelFlags::CapPressure );
if( useTotalMassEquation )
kernelFlags.set( KernelFlags::TotalMassEquation );
if( useNewGravity )
kernelFlags.set( KernelFlags::NewGravity );
if( upwindingParams.upwindingScheme == UpwindingScheme::C1PPU &&
isothermalCompositionalMultiphaseFVMKernelUtilities::epsC1PPU > 0 )
kernelFlags.set( KernelFlags::C1PPU );
else if( upwindingParams.upwindingScheme == UpwindingScheme::IHU )
kernelFlags.set( KernelFlags::IHU );


using kernelType = FluxComputeKernel< NUM_COMP, NUM_DOF, STENCILWRAPPER >;
typename kernelType::CompFlowAccessors compFlowAccessors( elemManager, solverName );
typename kernelType::MultiFluidAccessors multiFluidAccessors( elemManager, solverName );
Expand Down
14 changes: 8 additions & 6 deletions src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/FluxComputeKernelBase.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -45,17 +45,19 @@ enum class KernelFlags
{
/// Flag to specify whether capillary pressure is used or not
CapPressure = 1 << 0, // 1
/// Flag to specify whether diffusion is used or not
Diffusion = 1 << 1, // 2
/// Flag to specify whether dispersion is used or not
Dispersion = 1 << 2, // 4
/// Flag indicating whether total mass equation is formed or not
TotalMassEquation = 1 << 1, // 2
TotalMassEquation = 1 << 3, // 8
/// Flag indicating whether new gravity treatment is used or not
NewGravity = 1 << 2, // 4
NewGravity = 1 << 4, // 16
/// Flag indicating whether C1-PPU is used or not
C1PPU = 1 << 3, // 8
C1PPU = 1 << 5, // 32
/// Flag indicating whether IHU is used or not
IHU = 1 << 4 // 16
IHU = 1 << 6 // 64
/// Add more flags like that if needed:
// Flag6 = 1 << 5, // 32
// Flag7 = 1 << 6, // 64
// Flag8 = 1 << 7 //128
};

Expand Down
9 changes: 1 addition & 8 deletions ...Components/physicsSolvers/fluidFlow/kernels/compositional/StabilizedFluxComputeKernel.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -327,8 +327,7 @@ class FluxComputeKernelFactory
integer const numPhases,
globalIndex const rankOffset,
string const & dofKey,
integer const hasCapPressure,
integer const useTotalMassEquation,
BitFlags< isothermalCompositionalMultiphaseFVMKernels::KernelFlags > kernelFlags,
string const & solverName,
ElementRegionManager const & elemManager,
STENCILWRAPPER const & stencilWrapper,
Expand All @@ -346,12 +345,6 @@ class FluxComputeKernelFactory
elemManager.constructArrayViewAccessor< globalIndex, 1 >( dofKey );
dofNumberAccessor.setName( solverName + "/accessors/" + dofKey );

BitFlags< isothermalCompositionalMultiphaseFVMKernels::KernelFlags > kernelFlags;
if( hasCapPressure )
kernelFlags.set( isothermalCompositionalMultiphaseFVMKernels::KernelFlags::CapPressure );
if( useTotalMassEquation )
kernelFlags.set( isothermalCompositionalMultiphaseFVMKernels::KernelFlags::TotalMassEquation );

using KERNEL_TYPE = FluxComputeKernel< NUM_COMP, NUM_DOF, STENCILWRAPPER >;
typename KERNEL_TYPE::CompFlowAccessors compFlowAccessors( elemManager, solverName );
typename KERNEL_TYPE::MultiFluidAccessors multiFluidAccessors( elemManager, solverName );
Expand Down
Loading

0 comments on commit f562029

Please sign in to comment.