Skip to content

Commit

Permalink
wip
Browse files Browse the repository at this point in the history
  • Loading branch information
Pavel Tomin authored and Pavel Tomin committed Oct 19, 2024
1 parent 4b46fee commit 188f5db
Show file tree
Hide file tree
Showing 7 changed files with 115 additions and 252 deletions.
1 change: 1 addition & 0 deletions src/coreComponents/physicsSolvers/fluidFlow/CMakeLists.txt
View file
Original file line number Diff line number Diff line change
Expand Up @@ -62,6 +62,7 @@ set( physicsSolvers_headers
fluidFlow/kernels/compositional/GlobalComponentFractionKernel.hpp
fluidFlow/kernels/compositional/HydrostaticPressureKernel.hpp
fluidFlow/kernels/compositional/IHUPhaseFlux.hpp
fluidFlow/kernels/compositional/HU2PhaseFlux.hpp
fluidFlow/kernels/compositional/KernelLaunchSelector.hpp
fluidFlow/kernels/compositional/PhaseComponentFlux.hpp
fluidFlow/kernels/compositional/PhaseMobilityKernel.hpp
Expand Down
16 changes: 1 addition & 15 deletions src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/C1PPUPhaseFlux.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,6 @@
#include "constitutive/capillaryPressure/layouts.hpp"
#include "mesh/ElementRegionManager.hpp"
#include "physicsSolvers/fluidFlow/kernels/compositional/PotGrad.hpp"
#include "physicsSolvers/fluidFlow/kernels/compositional/PhaseComponentFlux.hpp"


namespace geos
Expand Down Expand Up @@ -89,21 +88,15 @@ struct C1PPUPhaseFlux
ElementViewConst< arrayView2d< real64 const, compflow::USD_PHASE > > const & phaseMob,
ElementViewConst< arrayView3d< real64 const, compflow::USD_PHASE_DC > > const & dPhaseMob,
ElementViewConst< arrayView3d< real64 const, compflow::USD_PHASE_DC > > const & dPhaseVolFrac,
ElementViewConst< arrayView4d< real64 const, constitutive::multifluid::USD_PHASE_COMP > > const & phaseCompFrac,
ElementViewConst< arrayView5d< real64 const, constitutive::multifluid::USD_PHASE_COMP_DC > > const & dPhaseCompFrac,
ElementViewConst< arrayView3d< real64 const, compflow::USD_COMP_DC > > const & dCompFrac_dCompDens,
ElementViewConst< arrayView3d< real64 const, constitutive::multifluid::USD_PHASE > > const & phaseMassDens,
ElementViewConst< arrayView4d< real64 const, constitutive::multifluid::USD_PHASE_DC > > const & dPhaseMassDens,
ElementViewConst< arrayView3d< real64 const, constitutive::cappres::USD_CAPPRES > > const & phaseCapPressure,
ElementViewConst< arrayView4d< real64 const, constitutive::cappres::USD_CAPPRES_DS > > const & dPhaseCapPressure_dPhaseVolFrac,
localIndex & k_up,
real64 & potGrad,
real64 ( &phaseFlux ),
real64 ( & dPhaseFlux_dP )[numFluxSupportPoints],
real64 ( & dPhaseFlux_dC )[numFluxSupportPoints][numComp],
real64 ( & compFlux )[numComp],
real64 ( & dCompFlux_dP )[numFluxSupportPoints][numComp],
real64 ( & dCompFlux_dC )[numFluxSupportPoints][numComp][numComp] )
real64 ( & dPhaseFlux_dC )[numFluxSupportPoints][numComp] )
{
real64 dPresGrad_dP[numFluxSupportPoints]{};
real64 dPresGrad_dC[numFluxSupportPoints][numComp]{};
Expand Down Expand Up @@ -209,13 +202,6 @@ struct C1PPUPhaseFlux
}

potGrad = potGrad * Ttrans;

// choose upstream cell for composition upwinding
k_up = (phaseFlux >= 0) ? 0 : 1;

//distribute on phaseComponentFlux here
PhaseComponentFlux::compute( ip, k_up, seri, sesri, sei, phaseCompFrac, dPhaseCompFrac, dCompFrac_dCompDens, phaseFlux
, dPhaseFlux_dP, dPhaseFlux_dC, compFlux, dCompFlux_dP, dCompFlux_dC );
}
};

Expand Down
55 changes: 17 additions & 38 deletions src/coreComponents/physicsSolvers/fluidFlow/kernels/compositional/FluxComputeKernel.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,6 +36,7 @@
#include "physicsSolvers/fluidFlow/kernels/compositional/C1PPUPhaseFlux.hpp"
#include "physicsSolvers/fluidFlow/kernels/compositional/IHUPhaseFlux.hpp"
#include "physicsSolvers/fluidFlow/kernels/compositional/HU2PhaseFlux.hpp"
#include "physicsSolvers/fluidFlow/kernels/compositional/PhaseComponentFlux.hpp"

namespace geos
{
Expand Down Expand Up @@ -229,6 +230,7 @@ class FluxComputeKernel : public FluxComputeKernelBase
StackVariables & stack,
FUNC && compFluxKernelOp = NoOpFunc{} ) const
{
using namespace isothermalCompositionalMultiphaseFVMKernelUtilities;

// first, compute the transmissibilities at this face
m_stencilWrapper.computeWeights( iconn,
Expand Down Expand Up @@ -270,11 +272,9 @@ class FluxComputeKernel : public FluxComputeKernelBase
real64 dPhaseFlux_dP[numFluxSupportPoints]{};
real64 dPhaseFlux_dC[numFluxSupportPoints][numComp]{};

localIndex k_up = -1;

if( m_kernelFlags.isSet( FluxComputeKernelFlags::C1PPU ) )
{
isothermalCompositionalMultiphaseFVMKernelUtilities::C1PPUPhaseFlux::compute< numComp, numFluxSupportPoints >
C1PPUPhaseFlux::compute< numComp, numFluxSupportPoints >
( m_numPhases,
ip,
m_kernelFlags.isSet( FluxComputeKernelFlags::CapPressure ),
Expand All @@ -285,22 +285,17 @@ class FluxComputeKernel : public FluxComputeKernelBase
m_gravCoef,
m_phaseMob, m_dPhaseMob,
m_dPhaseVolFrac,
m_phaseCompFrac, m_dPhaseCompFrac,
m_dCompFrac_dCompDens,
m_phaseMassDens, m_dPhaseMassDens,
m_phaseCapPressure, m_dPhaseCapPressure_dPhaseVolFrac,
k_up,
potGrad,
phaseFlux,
dPhaseFlux_dP,
dPhaseFlux_dC,
compFlux,
dCompFlux_dP,
dCompFlux_dC );
dPhaseFlux_dC );
}
else if( m_kernelFlags.isSet( FluxComputeKernelFlags::IHU ) )
{
isothermalCompositionalMultiphaseFVMKernelUtilities::IHUPhaseFlux::compute< numComp, numFluxSupportPoints >
IHUPhaseFlux::compute< numComp, numFluxSupportPoints >
( m_numPhases,
ip,
m_kernelFlags.isSet( FluxComputeKernelFlags::CapPressure ),
Expand All @@ -311,22 +306,17 @@ class FluxComputeKernel : public FluxComputeKernelBase
m_gravCoef,
m_phaseMob, m_dPhaseMob,
m_dPhaseVolFrac,
m_phaseCompFrac, m_dPhaseCompFrac,
m_dCompFrac_dCompDens,
m_phaseMassDens, m_dPhaseMassDens,
m_phaseCapPressure, m_dPhaseCapPressure_dPhaseVolFrac,
k_up,
potGrad,
phaseFlux,
dPhaseFlux_dP,
dPhaseFlux_dC,
compFlux,
dCompFlux_dP,
dCompFlux_dC );
dPhaseFlux_dC );
}
else if( m_kernelFlags.isSet( FluxComputeKernelFlags::HU2PH ) )
{
isothermalCompositionalMultiphaseFVMKernelUtilities::HU2PhaseFlux::compute< numComp, numFluxSupportPoints >
HU2PhaseFlux::compute< numComp, numFluxSupportPoints >
( m_numPhases,
ip,
m_kernelFlags.isSet( FluxComputeKernelFlags::CapPressure ),
Expand All @@ -337,22 +327,17 @@ class FluxComputeKernel : public FluxComputeKernelBase
m_gravCoef,
m_phaseMob, m_dPhaseMob,
m_dPhaseVolFrac,
m_phaseCompFrac, m_dPhaseCompFrac,
m_dCompFrac_dCompDens,
m_phaseMassDens, m_dPhaseMassDens,
m_phaseCapPressure, m_dPhaseCapPressure_dPhaseVolFrac,
k_up,
potGrad,
phaseFlux,
dPhaseFlux_dP,
dPhaseFlux_dC,
compFlux,
dCompFlux_dP,
dCompFlux_dC );
dPhaseFlux_dC );
}
else
{
isothermalCompositionalMultiphaseFVMKernelUtilities::PPUPhaseFlux::compute< numComp, numFluxSupportPoints >
PPUPhaseFlux::compute< numComp, numFluxSupportPoints >
( m_numPhases,
ip,
m_kernelFlags.isSet( FluxComputeKernelFlags::CapPressure ),
Expand All @@ -363,35 +348,29 @@ class FluxComputeKernel : public FluxComputeKernelBase
m_gravCoef,
m_phaseMob, m_dPhaseMob,
m_dPhaseVolFrac,
m_phaseCompFrac, m_dPhaseCompFrac,
m_dCompFrac_dCompDens,
m_phaseMassDens, m_dPhaseMassDens,
m_phaseCapPressure, m_dPhaseCapPressure_dPhaseVolFrac,
k_up,
potGrad,
phaseFlux,
dPhaseFlux_dP,
dPhaseFlux_dC,
compFlux,
dCompFlux_dP,
dCompFlux_dC );
dPhaseFlux_dC );
}

std::cout << ip << " " << phaseFlux << std::endl;
for( localIndex ke = 0; ke < numFluxSupportPoints; ++ke )
{
std::cout << dPhaseFlux_dP[ke] << std::endl;
for( localIndex ic = 0; ic < numComp; ++ic )
std::cout << dPhaseFlux_dC[ke][ic] << std::endl;
}
// choose upstream cell for composition upwinding
localIndex const k_up = (phaseFlux >= 0) ? 0 : 1;

//distribute on phaseComponentFlux here
PhaseComponentFlux::compute( ip, k_up, seri, sesri, sei, m_phaseCompFrac, m_dPhaseCompFrac, m_dCompFrac_dCompDens,
phaseFlux, dPhaseFlux_dP, dPhaseFlux_dC, compFlux, dCompFlux_dP, dCompFlux_dC );

// call the lambda in the phase loop to allow the reuse of the phase fluxes and their derivatives
// possible use: assemble the derivatives wrt temperature, and the flux term of the energy equation for this phase
compFluxKernelOp( ip, k, seri, sesri, sei, connectionIndex,
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
Loading

0 comments on commit 188f5db

Please sign in to comment.