diff --git a/src/coreComponents/constitutive/CMakeLists.txt b/src/coreComponents/constitutive/CMakeLists.txt index c763fd49e49..b8a3d02c811 100644 --- a/src/coreComponents/constitutive/CMakeLists.txt +++ b/src/coreComponents/constitutive/CMakeLists.txt @@ -97,8 +97,10 @@ set( constitutive_headers fluid/multifluid/compositional/models/CriticalVolume.hpp fluid/multifluid/compositional/models/EquationOfState.hpp fluid/multifluid/compositional/models/FunctionBase.hpp + fluid/multifluid/compositional/models/ImmiscibleWaterDensity.hpp fluid/multifluid/compositional/models/ImmiscibleWaterFlashModel.hpp fluid/multifluid/compositional/models/ImmiscibleWaterParameters.hpp + fluid/multifluid/compositional/models/ImmiscibleWaterViscosity.hpp fluid/multifluid/compositional/models/LohrenzBrayClarkViscosity.hpp fluid/multifluid/compositional/models/LohrenzBrayClarkViscosityImpl.hpp fluid/multifluid/compositional/models/NegativeTwoPhaseFlashModel.hpp @@ -242,8 +244,10 @@ set( constitutive_sources fluid/multifluid/compositional/models/CompositionalDensity.cpp fluid/multifluid/compositional/models/ConstantViscosity.cpp fluid/multifluid/compositional/models/CriticalVolume.cpp + fluid/multifluid/compositional/models/ImmiscibleWaterDensity.cpp fluid/multifluid/compositional/models/ImmiscibleWaterFlashModel.cpp fluid/multifluid/compositional/models/ImmiscibleWaterParameters.cpp + fluid/multifluid/compositional/models/ImmiscibleWaterViscosity.cpp fluid/multifluid/compositional/models/LohrenzBrayClarkViscosity.cpp fluid/multifluid/compositional/models/NegativeTwoPhaseFlashModel.cpp fluid/multifluid/compositional/CompositionalMultiphaseFluid.cpp diff --git a/src/coreComponents/constitutive/fluid/multifluid/MultiFluidSelector.hpp b/src/coreComponents/constitutive/fluid/multifluid/MultiFluidSelector.hpp index 580a965061b..7ca5afe50cc 100644 --- a/src/coreComponents/constitutive/fluid/multifluid/MultiFluidSelector.hpp +++ b/src/coreComponents/constitutive/fluid/multifluid/MultiFluidSelector.hpp @@ -53,6 +53,7 @@ void constitutiveUpdatePassThru( constitutive::MultiFluidBase const & fluid, #if !defined(GEOS_DEVICE_COMPILE) CO2BrineEzrokhiThermalFluid, CompositionalTwoPhaseLohrenzBrayClarkViscosity, + CompositionalThreePhaseLohrenzBrayClarkViscosity, #endif CompositionalTwoPhaseConstantViscosity >::execute( fluid, std::forward< LAMBDA >( lambda ) ); @@ -75,6 +76,7 @@ void constitutiveUpdatePassThru( constitutive::MultiFluidBase & fluid, #if !defined(GEOS_DEVICE_COMPILE) CO2BrineEzrokhiThermalFluid, CompositionalTwoPhaseLohrenzBrayClarkViscosity, + CompositionalThreePhaseLohrenzBrayClarkViscosity, #endif CompositionalTwoPhaseConstantViscosity >::execute( fluid, std::forward< LAMBDA >( lambda ) ); diff --git a/src/coreComponents/constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluid.cpp b/src/coreComponents/constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluid.cpp index aa91a802388..6e95fed3f4c 100644 --- a/src/coreComponents/constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluid.cpp +++ b/src/coreComponents/constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluid.cpp @@ -22,6 +22,7 @@ #include "constitutive/fluid/multifluid/CO2Brine/functions/PVTFunctionHelpers.hpp" #include "constitutive/fluid/multifluid/MultiFluidFields.hpp" #include "codingUtilities/Utilities.hpp" +#include "common/format/StringUtilities.hpp" namespace geos { @@ -51,6 +52,7 @@ CompositionalMultiphaseFluid( string const & name, Group * const parent ) m_parameters( createModelParameters() ) { using InputFlags = dataRepository::InputFlags; + using RestartFlags = dataRepository::RestartFlags; getWrapperBase( viewKeyStruct::componentNamesString() ).setInputFlag( InputFlags::REQUIRED ); getWrapperBase( viewKeyStruct::componentMolarWeightString() ).setInputFlag( InputFlags::REQUIRED ); @@ -80,13 +82,18 @@ CompositionalMultiphaseFluid( string const & name, Group * const parent ) // Link parameters specific to each model m_parameters->registerParameters( this ); + + // Register extra wrappers to enable auto-cloning + registerWrapper( "phaseOrder", &m_phaseOrder ) + .setSizedFromParent( 0 ) + .setRestartFlags( RestartFlags::NO_WRITE ); } template< typename FLASH, typename PHASE1, typename PHASE2, typename PHASE3 > integer CompositionalMultiphaseFluid< FLASH, PHASE1, PHASE2, PHASE3 >::getWaterPhaseIndex() const { - string const expectedWaterPhaseNames[] = { "water" }; - return PVTProps::PVTFunctionHelpers::findName( m_phaseNames, expectedWaterPhaseNames, viewKeyStruct::phaseNamesString() ); + integer const aqueous = static_cast< integer >(PhaseType::AQUEOUS); + return m_phaseOrder.size() > aqueous ? m_phaseOrder[aqueous] : -1; } template< typename FLASH, typename PHASE1, typename PHASE2, typename PHASE3 > @@ -170,6 +177,12 @@ void CompositionalMultiphaseFluid< FLASH, PHASE1, PHASE2, PHASE3 >::postInputIni } } + // Determine the phase ordering + m_phaseOrder.resize( 3 ); + m_phaseOrder[PhaseType::LIQUID] = findPhaseIndex( "oil,liq,liquid" ); + m_phaseOrder[PhaseType::VAPOUR] = findPhaseIndex( "gas,vap,vapor,vapour" ); + m_phaseOrder[PhaseType::AQUEOUS] = findPhaseIndex( "wat,water,aqueous" ); + m_parameters->postInputInitialization( this, *m_componentProperties ); } @@ -210,6 +223,7 @@ CompositionalMultiphaseFluid< FLASH, PHASE1, PHASE2, PHASE3 >::createKernelWrapp *m_phase1, *m_phase2, *m_phase3, + m_phaseOrder.toViewConst(), m_componentMolarWeight, m_useMass, m_phaseFraction.toView(), @@ -247,6 +261,22 @@ void CompositionalMultiphaseFluid< FLASH, PHASE1, PHASE2, PHASE3 >::createModels *m_parameters ); } +template< typename FLASH, typename PHASE1, typename PHASE2, typename PHASE3 > +integer CompositionalMultiphaseFluid< FLASH, PHASE1, PHASE2, PHASE3 >::findPhaseIndex( string names ) const +{ + auto const nameContainer = stringutilities::tokenize( names, ",", true, false ); + + for( integer ip = 0; ip < numFluidPhases(); ++ip ) + { + std::string const phaseName = stringutilities::toLower( m_phaseNames[ip] ); + if( std::find( nameContainer.begin(), nameContainer.end(), phaseName ) != nameContainer.end()) + { + return ip; + } + } + return -1; +} + // Create the fluid models template< typename FLASH, typename PHASE1, typename PHASE2, typename PHASE3 > std::unique_ptr< compositional::ModelParameters > @@ -269,6 +299,11 @@ template class CompositionalMultiphaseFluid< compositional::NegativeTwoPhaseFlashModel, compositional::PhaseModel< compositional::CompositionalDensity, compositional::LohrenzBrayClarkViscosity, compositional::NullModel >, compositional::PhaseModel< compositional::CompositionalDensity, compositional::LohrenzBrayClarkViscosity, compositional::NullModel > >; +template class CompositionalMultiphaseFluid< + compositional::ImmiscibleWaterFlashModel, + compositional::PhaseModel< compositional::CompositionalDensity, compositional::LohrenzBrayClarkViscosity, compositional::NullModel >, + compositional::PhaseModel< compositional::CompositionalDensity, compositional::LohrenzBrayClarkViscosity, compositional::NullModel >, + compositional::PhaseModel< compositional::ImmiscibleWaterDensity, compositional::ImmiscibleWaterViscosity, compositional::NullModel > >; REGISTER_CATALOG_ENTRY( ConstitutiveBase, CompositionalTwoPhaseConstantViscosity, @@ -280,6 +315,11 @@ REGISTER_CATALOG_ENTRY( ConstitutiveBase, string const &, dataRepository::Group * const ) +REGISTER_CATALOG_ENTRY( ConstitutiveBase, + CompositionalThreePhaseLohrenzBrayClarkViscosity, + string const &, + dataRepository::Group * const ) + } // namespace constitutive } // namespace geos diff --git a/src/coreComponents/constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluid.hpp b/src/coreComponents/constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluid.hpp index a8f846aac39..b43047ade63 100644 --- a/src/coreComponents/constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluid.hpp +++ b/src/coreComponents/constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluid.hpp @@ -23,6 +23,9 @@ #include "constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluidUpdates.hpp" #include "constitutive/fluid/multifluid/compositional/models/ConstantViscosity.hpp" #include "constitutive/fluid/multifluid/compositional/models/CompositionalDensity.hpp" +#include "constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterDensity.hpp" +#include "constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterFlashModel.hpp" +#include "constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterViscosity.hpp" #include "constitutive/fluid/multifluid/compositional/models/LohrenzBrayClarkViscosity.hpp" #include "constitutive/fluid/multifluid/compositional/models/NegativeTwoPhaseFlashModel.hpp" #include "constitutive/fluid/multifluid/compositional/models/ModelParameters.hpp" @@ -110,15 +113,27 @@ class CompositionalMultiphaseFluid : public MultiFluidBase virtual void resizeFields( localIndex const size, localIndex const numPts ) override; + enum PhaseType : integer + { + LIQUID = 0, + VAPOUR = 1, + AQUEOUS = 2, + }; + private: // Create the fluid models void createModels(); + integer findPhaseIndex( string names ) const; + static std::unique_ptr< compositional::ModelParameters > createModelParameters(); // Flash model std::unique_ptr< FLASH > m_flash{}; + // Phase ordering + array1d< integer > m_phaseOrder; + // Phase models std::unique_ptr< PHASE1 > m_phase1{}; std::unique_ptr< PHASE2 > m_phase2{}; @@ -142,6 +157,11 @@ using CompositionalTwoPhaseLohrenzBrayClarkViscosity = CompositionalMultiphaseFl compositional::NegativeTwoPhaseFlashModel, compositional::PhaseModel< compositional::CompositionalDensity, compositional::LohrenzBrayClarkViscosity, compositional::NullModel >, compositional::PhaseModel< compositional::CompositionalDensity, compositional::LohrenzBrayClarkViscosity, compositional::NullModel > >; +using CompositionalThreePhaseLohrenzBrayClarkViscosity = CompositionalMultiphaseFluid< + compositional::ImmiscibleWaterFlashModel, + compositional::PhaseModel< compositional::CompositionalDensity, compositional::LohrenzBrayClarkViscosity, compositional::NullModel >, + compositional::PhaseModel< compositional::CompositionalDensity, compositional::LohrenzBrayClarkViscosity, compositional::NullModel >, + compositional::PhaseModel< compositional::ImmiscibleWaterDensity, compositional::ImmiscibleWaterViscosity, compositional::NullModel > >; } /* namespace constitutive */ diff --git a/src/coreComponents/constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluidUpdates.hpp b/src/coreComponents/constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluidUpdates.hpp index 1567f4b4c46..c3a301e98a5 100644 --- a/src/coreComponents/constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluidUpdates.hpp +++ b/src/coreComponents/constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluidUpdates.hpp @@ -46,6 +46,7 @@ class CompositionalMultiphaseFluidUpdates final : public MultiFluidBase::KernelW PHASE1 const & phase1, PHASE2 const & phase2, PHASE3 const & phase3, + arrayView1d< integer const > const & phaseOrder, arrayView1d< real64 const > const & componentMolarWeight, bool const useMass, MultiFluidBase::PhaseProp::ViewType phaseFrac, @@ -122,6 +123,9 @@ class CompositionalMultiphaseFluidUpdates final : public MultiFluidBase::KernelW // Flash kernel wrapper typename FLASH::KernelWrapper m_flash; + // The ordering of phases + arrayView1d< integer const > const m_phaseOrder; + // Phase model kernel wrappers typename PHASE1::KernelWrapper m_phase1; typename PHASE2::KernelWrapper m_phase2; @@ -138,6 +142,7 @@ CompositionalMultiphaseFluidUpdates( compositional::ComponentProperties const & PHASE1 const & phase1, PHASE2 const & phase2, PHASE3 const & phase3, + arrayView1d< integer const > const & phaseOrder, arrayView1d< real64 const > const & componentMolarWeight, bool const useMass, MultiFluidBase::PhaseProp::ViewType phaseFrac, @@ -161,6 +166,7 @@ CompositionalMultiphaseFluidUpdates( compositional::ComponentProperties const & std::move( totalDensity ) ), m_componentProperties( componentProperties.createKernelWrapper() ), m_flash( flash.createKernelWrapper() ), + m_phaseOrder( phaseOrder ), m_phase1( phase1.createKernelWrapper() ), m_phase2( phase2.createKernelWrapper() ), m_phase3( phase3.createKernelWrapper() ), @@ -262,31 +268,31 @@ CompositionalMultiphaseFluidUpdates< FLASH, PHASE1, PHASE2, PHASE3 >::compute( m_phase1.density.compute( m_componentProperties, pressure, temperature, - phaseCompFrac.value[0].toSliceConst(), - phaseDens.value[0], - phaseDens.derivs[0], - phaseMassDensity.value[0], - phaseMassDensity.derivs[0], + phaseCompFrac.value[m_phaseOrder[0]].toSliceConst(), + phaseDens.value[m_phaseOrder[0]], + phaseDens.derivs[m_phaseOrder[0]], + phaseMassDensity.value[m_phaseOrder[0]], + phaseMassDensity.derivs[m_phaseOrder[0]], m_useMass ); m_phase2.density.compute( m_componentProperties, pressure, temperature, - phaseCompFrac.value[1].toSliceConst(), - phaseDens.value[1], - phaseDens.derivs[1], - phaseMassDensity.value[1], - phaseMassDensity.derivs[1], + phaseCompFrac.value[m_phaseOrder[1]].toSliceConst(), + phaseDens.value[m_phaseOrder[1]], + phaseDens.derivs[m_phaseOrder[1]], + phaseMassDensity.value[m_phaseOrder[1]], + phaseMassDensity.derivs[m_phaseOrder[1]], m_useMass ); if constexpr (2 < FLASH::KernelWrapper::getNumberOfPhases()) { m_phase3.density.compute( m_componentProperties, pressure, temperature, - phaseCompFrac.value[2].toSliceConst(), - phaseDens.value[2], - phaseDens.derivs[2], - phaseMassDensity.value[2], - phaseMassDensity.derivs[2], + phaseCompFrac.value[m_phaseOrder[2]].toSliceConst(), + phaseDens.value[m_phaseOrder[2]], + phaseDens.derivs[m_phaseOrder[2]], + phaseMassDensity.value[m_phaseOrder[2]], + phaseMassDensity.derivs[m_phaseOrder[2]], m_useMass ); } @@ -294,31 +300,31 @@ CompositionalMultiphaseFluidUpdates< FLASH, PHASE1, PHASE2, PHASE3 >::compute( m_phase1.viscosity.compute( m_componentProperties, pressure, temperature, - phaseCompFrac.value[0].toSliceConst(), - phaseMassDensity.value[0], - phaseMassDensity.derivs[0].toSliceConst(), - phaseVisc.value[0], - phaseVisc.derivs[0], + phaseCompFrac.value[m_phaseOrder[0]].toSliceConst(), + phaseMassDensity.value[m_phaseOrder[0]], + phaseMassDensity.derivs[m_phaseOrder[0]].toSliceConst(), + phaseVisc.value[m_phaseOrder[0]], + phaseVisc.derivs[m_phaseOrder[0]], m_useMass ); m_phase2.viscosity.compute( m_componentProperties, pressure, temperature, - phaseCompFrac.value[1].toSliceConst(), - phaseMassDensity.value[1], - phaseMassDensity.derivs[1].toSliceConst(), - phaseVisc.value[1], - phaseVisc.derivs[1], + phaseCompFrac.value[m_phaseOrder[1]].toSliceConst(), + phaseMassDensity.value[m_phaseOrder[1]], + phaseMassDensity.derivs[m_phaseOrder[1]].toSliceConst(), + phaseVisc.value[m_phaseOrder[1]], + phaseVisc.derivs[m_phaseOrder[1]], m_useMass ); if constexpr (2 < FLASH::KernelWrapper::getNumberOfPhases()) { m_phase3.viscosity.compute( m_componentProperties, pressure, temperature, - phaseCompFrac.value[2].toSliceConst(), - phaseMassDensity.value[2], - phaseMassDensity.derivs[2].toSliceConst(), - phaseVisc.value[2], - phaseVisc.derivs[2], + phaseCompFrac.value[m_phaseOrder[2]].toSliceConst(), + phaseMassDensity.value[m_phaseOrder[2]], + phaseMassDensity.derivs[m_phaseOrder[2]].toSliceConst(), + phaseVisc.value[m_phaseOrder[2]], + phaseVisc.derivs[m_phaseOrder[2]], m_useMass ); } diff --git a/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterDensity.cpp b/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterDensity.cpp new file mode 100644 index 00000000000..e257b3b984a --- /dev/null +++ b/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterDensity.cpp @@ -0,0 +1,82 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2024 Total, S.A + * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2018-2024 Chevron + * Copyright (c) 2019- GEOS/GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/** + * @file ImmiscibleWaterDensity.cpp + */ + +#include "ImmiscibleWaterDensity.hpp" +#include "ImmiscibleWaterParameters.hpp" +#include "dataRepository/InputFlags.hpp" + +namespace geos +{ + +namespace constitutive +{ + +namespace compositional +{ + +ImmiscibleWaterDensityUpdate::ImmiscibleWaterDensityUpdate( real64 const waterMolecularWeight, + real64 const referencePressure, + real64 const referenceTemperature, + real64 const density, + real64 const compressibility, + real64 const expansionCoefficient ): + m_waterMolecularWeight( waterMolecularWeight ), + m_referencePressure( referencePressure ), + m_referenceTemperature( referenceTemperature ), + m_density( density ), + m_compressibility( compressibility ), + m_expansionCoefficient( expansionCoefficient ) +{} + +ImmiscibleWaterDensity::ImmiscibleWaterDensity( string const & name, + ComponentProperties const & componentProperties, + integer const phaseIndex, + ModelParameters const & modelParameters ): + FunctionBase( name, componentProperties ), + m_parameters( modelParameters ) +{ + GEOS_UNUSED_VAR( phaseIndex ); + integer const h2oIndex = ImmiscibleWaterParameters::getWaterComponentIndex( componentProperties ); + GEOS_THROW_IF_LT_MSG( h2oIndex, 0, "Water component not found", InputError ); + m_waterMolecularWeight = componentProperties.getComponentMolarWeight()[h2oIndex]; +} + +ImmiscibleWaterDensity::KernelWrapper +ImmiscibleWaterDensity::createKernelWrapper() const +{ + ImmiscibleWaterParameters const * waterParameters = m_parameters.get< ImmiscibleWaterParameters >(); + return KernelWrapper( m_waterMolecularWeight, + waterParameters->m_waterReferencePressure, + waterParameters->m_waterReferenceTemperature, + waterParameters->m_waterDensity, + waterParameters->m_waterCompressibility, + waterParameters->m_waterExpansionCoefficient ); +} + +std::unique_ptr< ModelParameters > +ImmiscibleWaterDensity::createParameters( std::unique_ptr< ModelParameters > parameters ) +{ + return ImmiscibleWaterParameters::create( std::move( parameters ) ); +} + +} // namespace compositional + +} // namespace constitutive + +} // namespace geos diff --git a/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterDensity.hpp b/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterDensity.hpp new file mode 100644 index 00000000000..9d258085649 --- /dev/null +++ b/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterDensity.hpp @@ -0,0 +1,140 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2024 Total, S.A + * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2018-2024 Chevron + * Copyright (c) 2019- GEOS/GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/** + * @file ImmiscibleWaterDensity.cpp + */ + +#ifndef GEOS_CONSTITUTIVE_FLUID_MULTIFLUID_COMPOSITIONAL_MODELS_IMMISCIBLEWATERDENSITY_HPP_ +#define GEOS_CONSTITUTIVE_FLUID_MULTIFLUID_COMPOSITIONAL_MODELS_IMMISCIBLEWATERDENSITY_HPP_ + +#include "FunctionBase.hpp" +#include "constitutive/fluid/multifluid/Layouts.hpp" + +namespace geos +{ + +namespace constitutive +{ + +namespace compositional +{ + +class ImmiscibleWaterDensityUpdate final : public FunctionBaseUpdate +{ + using Deriv = geos::constitutive::multifluid::DerivativeOffset; +public: + ImmiscibleWaterDensityUpdate( real64 const waterMolecularWeight, + real64 const referencePressure, + real64 const referenceTemperature, + real64 const density, + real64 const compressibility, + real64 const expansionCoefficient ); + + template< integer USD1, integer USD2 > + GEOS_HOST_DEVICE + void compute( ComponentProperties::KernelWrapper const & componentProperties, + real64 const & pressure, + real64 const & temperature, + arraySlice1d< real64 const, USD1 > const & phaseComposition, + real64 & molarDensity, + arraySlice1d< real64, USD2 > const & dMolarDensity, + real64 & massDensity, + arraySlice1d< real64, USD2 > const & dMassDensity, + bool useMass ) const; + +private: + real64 const m_waterMolecularWeight; + real64 const m_referencePressure; + real64 const m_referenceTemperature; + real64 const m_density; + real64 const m_compressibility; + real64 const m_expansionCoefficient; +}; + +class ImmiscibleWaterDensity : public FunctionBase +{ +public: + ImmiscibleWaterDensity( string const & name, + ComponentProperties const & componentProperties, + integer const phaseIndex, + ModelParameters const & modelParameters ); + + static string catalogName() { return "ImmiscibleWaterDensity"; } + + virtual FunctionType functionType() const override + { + return FunctionType::DENSITY; + } + + /// Type of kernel wrapper for in-kernel update + using KernelWrapper = ImmiscibleWaterDensityUpdate; + + /** + * @brief Create an update kernel wrapper. + * @return the wrapper + */ + KernelWrapper createKernelWrapper() const; + + // Create parameters unique to this model + static std::unique_ptr< ModelParameters > createParameters( std::unique_ptr< ModelParameters > parameters ); + +private: + ModelParameters const & m_parameters; + real64 m_waterMolecularWeight{0.0}; +}; + +template< integer USD1, integer USD2 > +GEOS_HOST_DEVICE +void ImmiscibleWaterDensityUpdate::compute( + ComponentProperties::KernelWrapper const & componentProperties, + real64 const & pressure, + real64 const & temperature, + arraySlice1d< real64 const, USD1 > const & phaseComposition, + real64 & molarDensity, + arraySlice1d< real64, USD2 > const & dMolarDensity, + real64 & massDensity, + arraySlice1d< real64, USD2 > const & dMassDensity, + bool useMass ) const +{ + GEOS_UNUSED_VAR( componentProperties ); + GEOS_UNUSED_VAR( phaseComposition ); + GEOS_UNUSED_VAR( useMass ); + + LvArray::forValuesInSlice( dMolarDensity, setZero ); + LvArray::forValuesInSlice( dMassDensity, setZero ); + + real64 const density = m_density * + LvArray::math::exp( m_compressibility * (pressure - m_referencePressure) ) * + LvArray::math::exp( -m_expansionCoefficient * (temperature - m_referenceTemperature) ); + real64 const dDensity_dp = m_compressibility * density; + real64 const dDensity_dT = -m_expansionCoefficient * density; + + massDensity = density; + dMassDensity[Deriv::dP] = dDensity_dp; + dMassDensity[Deriv::dT] = dDensity_dT; + + molarDensity = density / m_waterMolecularWeight; + dMolarDensity[Deriv::dP] = dDensity_dp / m_waterMolecularWeight; + dMolarDensity[Deriv::dT] = dDensity_dT / m_waterMolecularWeight; +} + +} // end namespace compositional + +} // end namespace constitutive + +} // end namespace geos + +#endif //GEOS_CONSTITUTIVE_FLUID_MULTIFLUID_COMPOSITIONAL_MODELS_IMMISCIBLEWATERDENSITY_HPP_ diff --git a/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterParameters.cpp b/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterParameters.cpp index 87c1d659435..ddb64316270 100644 --- a/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterParameters.cpp +++ b/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterParameters.cpp @@ -63,7 +63,41 @@ integer ImmiscibleWaterParameters::getWaterComponentIndex( ComponentProperties c void ImmiscibleWaterParameters::registerParametersImpl( MultiFluidBase * fluid ) { - GEOS_UNUSED_VAR( fluid ); + fluid->registerWrapper( viewKeyStruct::waterReferencePressureString(), &m_waterReferencePressure ). + setInputFlag( dataRepository::InputFlags::REQUIRED ). + setDescription( "The reference pressure for water density and viscosity" ); + + fluid->registerWrapper( viewKeyStruct::waterReferenceTemperatureString(), &m_waterReferenceTemperature ). + setInputFlag( dataRepository::InputFlags::OPTIONAL ). + setDefaultValue( m_waterReferenceTemperature ). + setDescription( "The reference temperature for water density and viscosity" ); + + fluid->registerWrapper( viewKeyStruct::waterDensityString(), &m_waterDensity ). + setInputFlag( dataRepository::InputFlags::REQUIRED ). + setDescription( "The water density at the reference pressure and temperature" ); + + fluid->registerWrapper( viewKeyStruct::waterViscosityString(), &m_waterViscosity ). + setInputFlag( dataRepository::InputFlags::REQUIRED ). + setDescription( "The water viscosity at the reference pressure and temperature" ); + + fluid->registerWrapper( viewKeyStruct::waterCompressibilityString(), &m_waterCompressibility ). + setInputFlag( dataRepository::InputFlags::REQUIRED ). + setDescription( "The compressibility of water" ); + + fluid->registerWrapper( viewKeyStruct::waterViscosityCompressibilityString(), &m_waterViscosityCompressibility ). + setInputFlag( dataRepository::InputFlags::OPTIONAL ). + setDefaultValue( m_waterViscosityCompressibility ). + setDescription( "The compressibility (normalized derivative with respect to pressure) of the water viscosity" ); + + fluid->registerWrapper( viewKeyStruct::waterExpansionCoefficientString(), &m_waterExpansionCoefficient ). + setInputFlag( dataRepository::InputFlags::OPTIONAL ). + setDefaultValue( m_waterExpansionCoefficient ). + setDescription( "The volumetric coefficient of thermal expansion of water" ); + + fluid->registerWrapper( viewKeyStruct::waterViscosityExpansionCoefficientString(), &m_waterViscosityExpansionCoefficient ). + setInputFlag( dataRepository::InputFlags::OPTIONAL ). + setDefaultValue( m_waterViscosityExpansionCoefficient ). + setDescription( "The coefficient of thermal expansion (normalized derivative with respect to temperature) of water viscosity" ); } void ImmiscibleWaterParameters::postInputInitializationImpl( MultiFluidBase const * fluid, @@ -80,6 +114,24 @@ void ImmiscibleWaterParameters::postInputInitializationImpl( MultiFluidBase cons GEOS_FMT( "{}: water component not found '{}'", fluid->getFullName(), MultiFluidBase::viewKeyStruct::componentNamesString() ), InputError ); + + // Pretty much everything should be positive + auto const checkLowerBound = [&]( real64 const & value, real64 const & bound, string const & attribute ) + { + GEOS_THROW_IF_LT_MSG( value, bound, + GEOS_FMT( "{}: invalid number of value in attribute '{}'. Should be greater than {}", + fluid->getFullName(), bound, attribute ), + InputError ); + }; + + real64 constexpr epsilon = MultiFluidConstants::epsilon; + + checkLowerBound( m_waterDensity, epsilon, viewKeyStruct::waterDensityString()); + checkLowerBound( m_waterViscosity, epsilon, viewKeyStruct::waterViscosityString()); + checkLowerBound( m_waterCompressibility, 0.0, viewKeyStruct::waterCompressibilityString()); + checkLowerBound( m_waterViscosityCompressibility, 0.0, viewKeyStruct::waterViscosityCompressibilityString()); + checkLowerBound( m_waterExpansionCoefficient, 0.0, viewKeyStruct::waterExpansionCoefficientString()); + checkLowerBound( m_waterViscosityExpansionCoefficient, 0.0, viewKeyStruct::waterViscosityExpansionCoefficientString()); } } // end namespace compositional diff --git a/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterParameters.hpp b/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterParameters.hpp index 9404e565f83..da1e0e4ddf4 100644 --- a/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterParameters.hpp +++ b/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterParameters.hpp @@ -44,6 +44,27 @@ class ImmiscibleWaterParameters : public ModelParameters static integer getWaterComponentIndex( ComponentProperties const & componentProperties ); + struct viewKeyStruct + { + static constexpr char const * waterReferencePressureString() { return "waterReferencePressure"; } + static constexpr char const * waterReferenceTemperatureString() { return "waterReferenceTemperature"; } + static constexpr char const * waterDensityString() { return "waterDensity"; } + static constexpr char const * waterViscosityString() { return "waterViscosity"; } + static constexpr char const * waterCompressibilityString() { return "waterCompressibility"; } + static constexpr char const * waterViscosityCompressibilityString() { return "waterViscosityCompressibility"; } + static constexpr char const * waterExpansionCoefficientString() { return "waterExpansionCoefficient"; } + static constexpr char const * waterViscosityExpansionCoefficientString() { return "waterViscosityExpansionCoefficient"; } + }; + + real64 m_waterReferencePressure; + real64 m_waterReferenceTemperature{293.15}; + real64 m_waterDensity; + real64 m_waterViscosity; + real64 m_waterCompressibility; + real64 m_waterViscosityCompressibility{0.0}; + real64 m_waterExpansionCoefficient{0.0}; + real64 m_waterViscosityExpansionCoefficient{0.0}; + protected: void registerParametersImpl( MultiFluidBase * fluid ) override; diff --git a/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterViscosity.cpp b/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterViscosity.cpp new file mode 100644 index 00000000000..d1e5cd30867 --- /dev/null +++ b/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterViscosity.cpp @@ -0,0 +1,75 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2024 Total, S.A + * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2018-2024 Chevron + * Copyright (c) 2019- GEOS/GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/** + * @file ImmiscibleWaterViscosity.cpp + */ + +#include "ImmiscibleWaterViscosity.hpp" +#include "ImmiscibleWaterParameters.hpp" + +namespace geos +{ + +namespace constitutive +{ + +namespace compositional +{ + +ImmiscibleWaterViscosityUpdate::ImmiscibleWaterViscosityUpdate( real64 const referencePressure, + real64 const referenceTemperature, + real64 const viscosity, + real64 const compressibility, + real64 const expansionCoefficient ): + m_referencePressure( referencePressure ), + m_referenceTemperature( referenceTemperature ), + m_viscosity( viscosity ), + m_compressibility( compressibility ), + m_expansionCoefficient( expansionCoefficient ) +{} + +ImmiscibleWaterViscosity::ImmiscibleWaterViscosity( string const & name, + ComponentProperties const & componentProperties, + integer const phaseIndex, + ModelParameters const & modelParameters ): + FunctionBase( name, componentProperties ), + m_parameters( modelParameters ) +{ + GEOS_UNUSED_VAR( phaseIndex ); +} + +ImmiscibleWaterViscosity::KernelWrapper +ImmiscibleWaterViscosity::createKernelWrapper() const +{ + ImmiscibleWaterParameters const * waterParameters = m_parameters.get< ImmiscibleWaterParameters >(); + return KernelWrapper( waterParameters->m_waterReferencePressure, + waterParameters->m_waterReferenceTemperature, + waterParameters->m_waterViscosity, + waterParameters->m_waterViscosityCompressibility, + waterParameters->m_waterViscosityExpansionCoefficient ); +} + +std::unique_ptr< ModelParameters > +ImmiscibleWaterViscosity::createParameters( std::unique_ptr< ModelParameters > parameters ) +{ + return ImmiscibleWaterParameters::create( std::move( parameters ) ); +} + +} // namespace compositional + +} // namespace constitutive + +} // end namespace geos diff --git a/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterViscosity.hpp b/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterViscosity.hpp new file mode 100644 index 00000000000..815769110ab --- /dev/null +++ b/src/coreComponents/constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterViscosity.hpp @@ -0,0 +1,130 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2024 Total, S.A + * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2018-2024 Chevron + * Copyright (c) 2019- GEOS/GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/** + * @file ImmiscibleWaterViscosity.hpp + */ + +#ifndef GEOS_CONSTITUTIVE_FLUID_MULTIFLUID_COMPOSITIONAL_MODELS_IMMISCIBLEWATERVISCOSITY_HPP_ +#define GEOS_CONSTITUTIVE_FLUID_MULTIFLUID_COMPOSITIONAL_MODELS_IMMISCIBLEWATERVISCOSITY_HPP_ + +#include "FunctionBase.hpp" +#include "constitutive/fluid/multifluid/Layouts.hpp" + +namespace geos +{ + +namespace constitutive +{ + +namespace compositional +{ + +class ImmiscibleWaterViscosityUpdate final : public FunctionBaseUpdate +{ + using Deriv = geos::constitutive::multifluid::DerivativeOffset; +public: + ImmiscibleWaterViscosityUpdate( real64 const referencePressure, + real64 const referenceTemperature, + real64 const viscosity, + real64 const compressibility, + real64 const expansionCoefficient ); + + template< integer USD1, integer USD2 > + GEOS_HOST_DEVICE + void compute( ComponentProperties::KernelWrapper const & componentProperties, + real64 const & pressure, + real64 const & temperature, + arraySlice1d< real64 const, USD1 > const & phaseComposition, + real64 const & density, + arraySlice1d< real64 const, USD2 > const & dDensity, + real64 & viscosity, + arraySlice1d< real64, USD2 > const & dViscosity, + bool useMass ) const; + +private: + real64 const m_referencePressure; + real64 const m_referenceTemperature; + real64 const m_viscosity; + real64 const m_compressibility; + real64 const m_expansionCoefficient; +}; + +class ImmiscibleWaterViscosity : public FunctionBase +{ +public: + ImmiscibleWaterViscosity( string const & name, + ComponentProperties const & componentProperties, + integer const phaseIndex, + ModelParameters const & modelParameters ); + + static string catalogName() { return ""; } + + FunctionType functionType() const override + { + return FunctionType::VISCOSITY; + } + + /// Type of kernel wrapper for in-kernel update + using KernelWrapper = ImmiscibleWaterViscosityUpdate; + + /** + * @brief Create an update kernel wrapper. + * @return the wrapper + */ + KernelWrapper createKernelWrapper() const; + + // Create parameters unique to this model + static std::unique_ptr< ModelParameters > createParameters( std::unique_ptr< ModelParameters > parameters ); + +private: + ModelParameters const & m_parameters; +}; + +template< integer USD1, integer USD2 > +GEOS_HOST_DEVICE +GEOS_FORCE_INLINE +void ImmiscibleWaterViscosityUpdate::compute( ComponentProperties::KernelWrapper const & componentProperties, + real64 const & pressure, + real64 const & temperature, + arraySlice1d< real64 const, USD1 > const & phaseComposition, + real64 const & density, + arraySlice1d< real64 const, USD2 > const & dDensity, + real64 & viscosity, + arraySlice1d< real64, USD2 > const & dViscosity, + bool useMass ) const +{ + GEOS_UNUSED_VAR( componentProperties ); + GEOS_UNUSED_VAR( phaseComposition ); + GEOS_UNUSED_VAR( density ); + GEOS_UNUSED_VAR( dDensity ); + GEOS_UNUSED_VAR( useMass ); + + LvArray::forValuesInSlice( dViscosity, setZero ); + + viscosity = m_viscosity * + LvArray::math::exp( m_compressibility * (pressure - m_referencePressure) ) * + LvArray::math::exp( -m_expansionCoefficient * (temperature - m_referenceTemperature) ); + dViscosity[Deriv::dP] = m_compressibility * viscosity; + dViscosity[Deriv::dT] = -m_expansionCoefficient * viscosity; +} + +} // end namespace compositional + +} // end namespace constitutive + +} // end namespace geos + +#endif //GEOS_CONSTITUTIVE_FLUID_MULTIFLUID_COMPOSITIONAL_MODELS_IMMISCIBLEWATERVISCOSITY_HPP_ diff --git a/src/coreComponents/constitutive/unitTests/CMakeLists.txt b/src/coreComponents/constitutive/unitTests/CMakeLists.txt index d357ad9290c..aa095d493c8 100644 --- a/src/coreComponents/constitutive/unitTests/CMakeLists.txt +++ b/src/coreComponents/constitutive/unitTests/CMakeLists.txt @@ -8,6 +8,7 @@ set( gtest_geosx_tests testElasticIsotropic.cpp testKValueInitialization.cpp testImmiscibleWaterFlashModel.cpp + testImmiscibleWaterProperties.cpp testLohrenzBrayClarkViscosity.cpp testModifiedCamClay.cpp testMultiFluidSelector.cpp diff --git a/src/coreComponents/constitutive/unitTests/testImmiscibleWaterProperties.cpp b/src/coreComponents/constitutive/unitTests/testImmiscibleWaterProperties.cpp new file mode 100644 index 00000000000..6b31ceffc9e --- /dev/null +++ b/src/coreComponents/constitutive/unitTests/testImmiscibleWaterProperties.cpp @@ -0,0 +1,285 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2024 Total, S.A + * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2023-2024 Chevron + * Copyright (c) 2019- GEOS/GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +// Source includes +#include "codingUtilities/UnitTestUtilities.hpp" +#include "constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterParameters.hpp" +#include "constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterDensity.hpp" +#include "constitutive/fluid/multifluid/compositional/models/ImmiscibleWaterViscosity.hpp" +#include "TestFluid.hpp" +#include "TestFluidUtilities.hpp" + +using namespace geos::constitutive::compositional; + +namespace geos +{ +namespace testing +{ + +template< int NC > +using TestData = std::tuple< + real64 const, // pressure + real64 const, // temperature + Feed< NC > const, // phase composition + real64 const, // expected molar density + real64 const, // expected mass density + real64 const // expected viscosity + >; + +template< int NC > +struct FluidData {}; + +template<> +struct FluidData< 3 > +{ + static std::unique_ptr< TestFluid< 3 > > createFluid() + { + auto fluid = TestFluid< 3 >::create( {Fluid::C1, Fluid::C10, Fluid::H2O} ); + std::array< real64, 3 > const bics = {0.25, 0.0, 0.0}; + fluid->setBinaryCoefficients( bics ); + return fluid; + } +}; + +template< int NC > +class ImmiscibleWaterPropertiesTestFixture : public ::testing::TestWithParam< TestData< NC > > +{ + static constexpr real64 relTol = 1.0e-5; + static constexpr real64 absTol = 1.0e-7; + static constexpr int numComps = NC; + static constexpr int numDofs = NC + 2; + using Deriv = geos::constitutive::multifluid::DerivativeOffset; +public: + ImmiscibleWaterPropertiesTestFixture() + : m_fluid( FluidData< NC >::createFluid() ) + { + ComponentProperties const & componentProperties = this->m_fluid->getComponentProperties(); + + m_parameters = ImmiscibleWaterDensity::createParameters( std::make_unique< ModelParameters >() ); + m_parameters = ImmiscibleWaterViscosity::createParameters( std::move( m_parameters ) ); + + auto * waterParameters = const_cast< ImmiscibleWaterParameters * >(m_parameters->get< ImmiscibleWaterParameters >()); + waterParameters->m_waterReferencePressure = 215.0e5; + waterParameters->m_waterDensity = 1020.0; + waterParameters->m_waterCompressibility = 4.1483E-10; + waterParameters->m_waterViscosityCompressibility = 2.0E-11; + waterParameters->m_waterViscosity = 0.32929e-3; + + m_density = std::make_unique< ImmiscibleWaterDensity >( "PhaseDensity", componentProperties, 0, *m_parameters ); + m_viscosity = std::make_unique< ImmiscibleWaterViscosity >( "PhaseViscosity", componentProperties, 0, *m_parameters ); + } + + ~ImmiscibleWaterPropertiesTestFixture() = default; + + void testProperties( TestData< NC > const & data ) + { + real64 const pressure = std::get< 0 >( data ); + real64 const temperature = std::get< 1 >( data ); + stackArray1d< real64, numComps > phaseComposition; + TestFluid< NC >::createArray( phaseComposition, std::get< 2 >( data )); + real64 const expectedMolarDensity = std::get< 3 >( data ); + real64 const expectedMassDensity = std::get< 4 >( data ); + real64 const expectedViscosity = std::get< 5 >( data ); + + real64 molarDensity = 0.0; + real64 massDensity = 0.0; + real64 viscosity = 0.0; + stackArray1d< real64, numDofs > tempDerivs( numDofs ); + + auto componentProperties = m_fluid->createKernelWrapper(); + auto densityKernelWrapper = m_density->createKernelWrapper(); + auto viscosityKernelWrapper = m_viscosity->createKernelWrapper(); + + densityKernelWrapper.compute( componentProperties, + pressure, + temperature, + phaseComposition.toSliceConst(), + molarDensity, + tempDerivs.toSlice(), + massDensity, + tempDerivs.toSlice(), + false ); + + viscosityKernelWrapper.compute( componentProperties, + pressure, + temperature, + phaseComposition.toSliceConst(), + massDensity, + tempDerivs.toSliceConst(), + viscosity, + tempDerivs.toSlice(), + false ); + + checkRelativeError( molarDensity, expectedMolarDensity, relTol, absTol ); + checkRelativeError( massDensity, expectedMassDensity, relTol, absTol ); + checkRelativeError( viscosity, expectedViscosity, relTol, absTol ); + } + + void testPropertyDerivatives( TestData< NC > const & data ) + { + real64 const pressure = std::get< 0 >( data ); + real64 const temperature = std::get< 1 >( data ); + stackArray1d< real64, numComps > phaseComposition; + TestFluid< NC >::createArray( phaseComposition, std::get< 2 >( data )); + + auto componentProperties = m_fluid->createKernelWrapper(); + auto densityKernelWrapper = m_density->createKernelWrapper(); + auto viscosityKernelWrapper = m_viscosity->createKernelWrapper(); + + real64 molarDensity = 0.0; + real64 massDensity = 0.0; + real64 viscosity = 0.0; + stackArray1d< real64, numDofs > molarDensityDerivs( numDofs ); + stackArray1d< real64, numDofs > massDensityDerivs( numDofs ); + stackArray1d< real64, numDofs > viscosityDerivs( numDofs ); + stackArray1d< real64, 3 > derivatives( 3 ); + + densityKernelWrapper.compute( componentProperties, + pressure, + temperature, + phaseComposition.toSliceConst(), + molarDensity, + molarDensityDerivs.toSlice(), + massDensity, + massDensityDerivs.toSlice(), + false ); + + viscosityKernelWrapper.compute( componentProperties, + pressure, + temperature, + phaseComposition.toSliceConst(), + massDensity, + massDensityDerivs.toSliceConst(), + viscosity, + viscosityDerivs.toSlice(), + false ); + + // Viscosity values are very small so we will inflate the values to avoid false positives due + // to the absolute value check + real64 constexpr viscosityScale = 1.0e6; + + auto calculateProperties = [&]( real64 const p, real64 const t, auto const & zmf, auto & values ) { + stackArray1d< real64, numDofs > tempDerivs( numDofs ); + densityKernelWrapper.compute( componentProperties, p, t, zmf.toSliceConst(), + values[0], tempDerivs.toSlice(), values[1], tempDerivs.toSlice(), false ); + viscosityKernelWrapper.compute( componentProperties, p, t, zmf.toSliceConst(), + values[1], tempDerivs.toSliceConst(), values[2], tempDerivs.toSlice(), false ); + values[2] *= viscosityScale; + }; + + auto concatDerivatives = [&]( int idof ){ + derivatives[0] = molarDensityDerivs[idof]; + derivatives[1] = massDensityDerivs[idof]; + derivatives[2] = viscosityScale * viscosityDerivs[idof]; + }; + + // Compare against numerical derivatives + // -- Pressure derivative + concatDerivatives( Deriv::dP ); + real64 const dp = 1.0e-4 * pressure; + internal::testNumericalDerivative< 3 >( pressure, dp, derivatives.toSliceConst(), + [&]( real64 const p, auto & values ) { + calculateProperties( p, temperature, phaseComposition, values ); + }, absTol, relTol ); + + // -- Temperature derivative + concatDerivatives( Deriv::dT ); + real64 const dT = 1.0e-6 * temperature; + internal::testNumericalDerivative< 3 >( temperature, dT, derivatives.toSliceConst(), + [&]( real64 const t, auto & values ) { + calculateProperties( pressure, t, phaseComposition, values ); + }, absTol, relTol ); + + // -- Composition derivatives derivative + real64 const dz = 1.0e-7; + for( integer ic = 0; ic < NC; ++ic ) + { + concatDerivatives( Deriv::dC+ic ); + internal::testNumericalDerivative< 3 >( 0.0, dz, derivatives.toSliceConst(), + [&]( real64 const z, auto & values ) { + stackArray1d< real64, numComps > zmf( numComps ); + for( integer jc = 0; jc < numComps; ++jc ) + { + zmf[jc] = phaseComposition[jc]; + } + zmf[ic] += z; + calculateProperties( pressure, temperature, zmf, values ); + }, absTol, relTol ); + } + } + +protected: + std::unique_ptr< TestFluid< NC > > m_fluid{}; + std::unique_ptr< ImmiscibleWaterDensity > m_density{}; + std::unique_ptr< ImmiscibleWaterViscosity > m_viscosity{}; + std::unique_ptr< ModelParameters > m_parameters{}; +}; + +using ImmiscibleWaterProperties3 = ImmiscibleWaterPropertiesTestFixture< 3 >; + +TEST_P( ImmiscibleWaterProperties3, testProperties ) +{ + testProperties( GetParam() ); +} + +TEST_P( ImmiscibleWaterProperties3, testPropertyDerivatives ) +{ + testPropertyDerivatives( GetParam() ); +} + +//------------------------------------------------------------------------------- +// Data +//------------------------------------------------------------------------------- + +/* UNCRUSTIFY-OFF */ + +INSTANTIATE_TEST_SUITE_P( + ImmiscibleWaterProperties, ImmiscibleWaterProperties3, + ::testing::Values( + TestData< 3 >( 2.0e+06, 293.15, { 0.30, 0.30, 0.40 }, 5.616333e+04, 1.011782e+03, 3.291616e-04 ), + TestData< 3 >( 1.5e+07, 293.15, { 0.30, 0.30, 0.40 }, 5.646702e+04, 1.017253e+03, 3.292472e-04 ), + TestData< 3 >( 6.0e+07, 293.15, { 0.30, 0.30, 0.40 }, 5.753101e+04, 1.036421e+03, 3.295437e-04 ), + TestData< 3 >( 2.0e+06, 353.15, { 0.30, 0.30, 0.40 }, 5.616333e+04, 1.011782e+03, 3.291616e-04 ), + TestData< 3 >( 1.5e+07, 353.15, { 0.30, 0.30, 0.40 }, 5.646702e+04, 1.017253e+03, 3.292472e-04 ), + TestData< 3 >( 6.0e+07, 353.15, { 0.30, 0.30, 0.40 }, 5.753101e+04, 1.036421e+03, 3.295437e-04 ), + TestData< 3 >( 2.0e+06, 573.15, { 0.30, 0.30, 0.40 }, 5.616333e+04, 1.011782e+03, 3.291616e-04 ), + TestData< 3 >( 1.5e+07, 573.15, { 0.30, 0.30, 0.40 }, 5.646702e+04, 1.017253e+03, 3.292472e-04 ), + TestData< 3 >( 6.0e+07, 573.15, { 0.30, 0.30, 0.40 }, 5.753101e+04, 1.036421e+03, 3.295437e-04 ), + TestData< 3 >( 2.0e+06, 293.15, { 0.00, 0.00, 1.00 }, 5.616333e+04, 1.011782e+03, 3.291616e-04 ), + TestData< 3 >( 1.5e+07, 293.15, { 0.00, 0.00, 1.00 }, 5.646702e+04, 1.017253e+03, 3.292472e-04 ), + TestData< 3 >( 6.0e+07, 293.15, { 0.00, 0.00, 1.00 }, 5.753101e+04, 1.036421e+03, 3.295437e-04 ), + TestData< 3 >( 2.0e+06, 353.15, { 0.00, 0.00, 1.00 }, 5.616333e+04, 1.011782e+03, 3.291616e-04 ), + TestData< 3 >( 1.5e+07, 353.15, { 0.00, 0.00, 1.00 }, 5.646702e+04, 1.017253e+03, 3.292472e-04 ), + TestData< 3 >( 6.0e+07, 353.15, { 0.00, 0.00, 1.00 }, 5.753101e+04, 1.036421e+03, 3.295437e-04 ), + TestData< 3 >( 2.0e+06, 573.15, { 0.00, 0.00, 1.00 }, 5.616333e+04, 1.011782e+03, 3.291616e-04 ), + TestData< 3 >( 1.5e+07, 573.15, { 0.00, 0.00, 1.00 }, 5.646702e+04, 1.017253e+03, 3.292472e-04 ), + TestData< 3 >( 6.0e+07, 573.15, { 0.00, 0.00, 1.00 }, 5.753101e+04, 1.036421e+03, 3.295437e-04 ), + TestData< 3 >( 2.0e+06, 293.15, { 0.20, 0.80, 0.00 }, 5.616333e+04, 1.011782e+03, 3.291616e-04 ), + TestData< 3 >( 1.5e+07, 293.15, { 0.20, 0.80, 0.00 }, 5.646702e+04, 1.017253e+03, 3.292472e-04 ), + TestData< 3 >( 6.0e+07, 293.15, { 0.20, 0.80, 0.00 }, 5.753101e+04, 1.036421e+03, 3.295437e-04 ), + TestData< 3 >( 2.0e+06, 353.15, { 0.20, 0.80, 0.00 }, 5.616333e+04, 1.011782e+03, 3.291616e-04 ), + TestData< 3 >( 1.5e+07, 353.15, { 0.20, 0.80, 0.00 }, 5.646702e+04, 1.017253e+03, 3.292472e-04 ), + TestData< 3 >( 6.0e+07, 353.15, { 0.20, 0.80, 0.00 }, 5.753101e+04, 1.036421e+03, 3.295437e-04 ), + TestData< 3 >( 2.0e+06, 573.15, { 0.20, 0.80, 0.00 }, 5.616333e+04, 1.011782e+03, 3.291616e-04 ), + TestData< 3 >( 1.5e+07, 573.15, { 0.20, 0.80, 0.00 }, 5.646702e+04, 1.017253e+03, 3.292472e-04 ), + TestData< 3 >( 6.0e+07, 573.15, { 0.20, 0.80, 0.00 }, 5.753101e+04, 1.036421e+03, 3.295437e-04 ) + ) +); + +/* UNCRUSTIFY-ON */ + +} // testing + +} // geos diff --git a/src/coreComponents/constitutiveDrivers/CMakeLists.txt b/src/coreComponents/constitutiveDrivers/CMakeLists.txt index 376d4d34573..4cc69dd4935 100644 --- a/src/coreComponents/constitutiveDrivers/CMakeLists.txt +++ b/src/coreComponents/constitutiveDrivers/CMakeLists.txt @@ -39,6 +39,7 @@ set( constitutiveDrivers_sources fluid/multiFluid/CO2Brine/PVTDriverRunTestCO2BrinePhillipsThermalFluid.cpp fluid/multiFluid/CO2Brine/PVTDriverRunTestCO2BrineEzrokhiFluid.cpp fluid/multiFluid/CO2Brine/PVTDriverRunTestCO2BrineEzrokhiThermalFluid.cpp + fluid/multiFluid/compositional/PVTDriverRunTestCompositionalThreePhaseLohrenzBrayClarkViscosity.cpp fluid/multiFluid/compositional/PVTDriverRunTestCompositionalTwoPhaseConstantViscosity.cpp fluid/multiFluid/compositional/PVTDriverRunTestCompositionalTwoPhaseLohrenzBrayClarkViscosity.cpp fluid/multiFluid/reactive/ReactiveFluidDriver.cpp diff --git a/src/coreComponents/constitutiveDrivers/fluid/multiFluid/compositional/PVTDriverRunTestCompositionalThreePhaseLohrenzBrayClarkViscosity.cpp b/src/coreComponents/constitutiveDrivers/fluid/multiFluid/compositional/PVTDriverRunTestCompositionalThreePhaseLohrenzBrayClarkViscosity.cpp new file mode 100644 index 00000000000..c7ae9568177 --- /dev/null +++ b/src/coreComponents/constitutiveDrivers/fluid/multiFluid/compositional/PVTDriverRunTestCompositionalThreePhaseLohrenzBrayClarkViscosity.cpp @@ -0,0 +1,27 @@ +/* + * ------------------------------------------------------------------------------------------------------------ + * SPDX-License-Identifier: LGPL-2.1-only + * + * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC + * Copyright (c) 2018-2024 Total, S.A + * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University + * Copyright (c) 2018-2024 Chevron + * Copyright (c) 2019- GEOS/GEOSX Contributors + * All rights reserved + * + * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details. + * ------------------------------------------------------------------------------------------------------------ + */ + +/* + * PVTDriverRunTestCompositionalThreePhaseLohrenzBrayClarkViscosity.cpp + */ + +#include "constitutiveDrivers/fluid/multiFluid/PVTDriverRunTest.hpp" +#include "constitutive/fluid/multifluid/compositional/CompositionalMultiphaseFluid.hpp" + +namespace geos +{ +template void PVTDriver::runTest< constitutive::CompositionalThreePhaseLohrenzBrayClarkViscosity >( + constitutive::CompositionalThreePhaseLohrenzBrayClarkViscosity &, arrayView2d< real64 > const & ); +} diff --git a/src/coreComponents/schema/schema.xsd b/src/coreComponents/schema/schema.xsd index 6617fd74792..3235c8a9794 100644 --- a/src/coreComponents/schema/schema.xsd +++ b/src/coreComponents/schema/schema.xsd @@ -609,6 +609,10 @@ + + + + @@ -4817,6 +4821,7 @@ Level 0 outputs no specific information for this solver. Higher levels require m + @@ -5149,6 +5154,55 @@ The expected format is "{ waterMax, oilMax }", in that order--> + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/src/coreComponents/schema/schema.xsd.other b/src/coreComponents/schema/schema.xsd.other index c2a891421d3..5de00ab14f2 100644 --- a/src/coreComponents/schema/schema.xsd.other +++ b/src/coreComponents/schema/schema.xsd.other @@ -1453,6 +1453,7 @@ + @@ -1931,6 +1932,56 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + @@ -1970,6 +2021,8 @@ + + @@ -2018,6 +2071,8 @@ + + diff --git a/src/coreComponents/unitTests/constitutiveTests/CMakeLists.txt b/src/coreComponents/unitTests/constitutiveTests/CMakeLists.txt index 523c6c89500..f59c739a956 100644 --- a/src/coreComponents/unitTests/constitutiveTests/CMakeLists.txt +++ b/src/coreComponents/unitTests/constitutiveTests/CMakeLists.txt @@ -28,6 +28,7 @@ set( gtest_pvt_xmls testPVT_CO2Brine.xml testPVT_CO2BrineTables.xml testPVT_PhaseComposition.xml + testPVT_ThreePhaseCompositional.xml ) set( gtest_reactivefluid_xmls diff --git a/src/coreComponents/unitTests/constitutiveTests/testPVT_ThreePhaseCompositional.xml b/src/coreComponents/unitTests/constitutiveTests/testPVT_ThreePhaseCompositional.xml new file mode 100644 index 00000000000..79c4892078e --- /dev/null +++ b/src/coreComponents/unitTests/constitutiveTests/testPVT_ThreePhaseCompositional.xml @@ -0,0 +1,83 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/src/coreComponents/unitTests/constitutiveTests/testPVT_ThreePhaseCompositional_blackOil.txt b/src/coreComponents/unitTests/constitutiveTests/testPVT_ThreePhaseCompositional_blackOil.txt new file mode 100644 index 00000000000..53e20dce320 --- /dev/null +++ b/src/coreComponents/unitTests/constitutiveTests/testPVT_ThreePhaseCompositional_blackOil.txt @@ -0,0 +1,113 @@ +# column 1 = time +# column 2 = pressure +# column 3 = temperature +# column 4 = density +# column 5 = total compressibility +# columns 6-8 = phase fractions +# columns 9-11 = phase densities +# columns 12-14 = phase mass densities +# columns 15-17 = phase viscosities +# columns 18-20 = oil phase fractions [C1, C7+, H2O] +# columns 21-23 = gas phase fractions [C1, C7+, H2O] +# columns 24-26 = water phase fractions [C1, C7+, H2O] +0.0000e+00 3.5000e+07 5.5315e+02 6.1560e+03 6.0169e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.5295e+03 5.5295e+03 5.6953e+04 5.2265e+02 5.2265e+02 1.0257e+03 1.0240e-04 1.0240e-04 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +2.0000e-02 3.4000e+07 5.5315e+02 6.1183e+03 6.2707e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.4953e+03 5.4953e+03 5.6930e+04 5.1941e+02 5.1941e+02 1.0253e+03 1.0027e-04 1.0027e-04 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +4.0000e-02 3.3000e+07 5.5315e+02 6.0792e+03 6.5439e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.4598e+03 5.4598e+03 5.6906e+04 5.1606e+02 5.1606e+02 1.0249e+03 9.8135e-05 9.8135e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +6.0000e-02 3.2000e+07 5.5315e+02 6.0387e+03 6.8385e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.4230e+03 5.4230e+03 5.6882e+04 5.1258e+02 5.1258e+02 1.0245e+03 9.6007e-05 9.6007e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +8.0000e-02 3.1000e+07 5.5315e+02 5.9966e+03 7.1573e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.3847e+03 5.3847e+03 5.6859e+04 5.0896e+02 5.0896e+02 1.0240e+03 9.3883e-05 9.3883e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.0000e-01 3.0000e+07 5.5315e+02 5.9528e+03 7.5032e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.3450e+03 5.3450e+03 5.6835e+04 5.0521e+02 5.0521e+02 1.0236e+03 9.1761e-05 9.1761e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.2000e-01 2.9000e+07 5.5315e+02 5.9072e+03 7.8796e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.3036e+03 5.3036e+03 5.6812e+04 5.0129e+02 5.0129e+02 1.0232e+03 8.9642e-05 8.9642e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.4000e-01 2.8000e+07 5.5315e+02 5.8597e+03 8.2906e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.2604e+03 5.2604e+03 5.6788e+04 4.9721e+02 4.9721e+02 1.0228e+03 8.7525e-05 8.7525e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.6000e-01 2.7000e+07 5.5315e+02 5.8100e+03 8.7412e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.2153e+03 5.2153e+03 5.6765e+04 4.9295e+02 4.9295e+02 1.0223e+03 8.5409e-05 8.5409e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.8000e-01 2.6000e+07 5.5315e+02 5.7580e+03 9.2370e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.1682e+03 5.1682e+03 5.6741e+04 4.8849e+02 4.8849e+02 1.0219e+03 8.3293e-05 8.3293e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +2.0000e-01 2.5000e+07 5.5315e+02 5.7036e+03 9.7851e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.1187e+03 5.1187e+03 5.6717e+04 4.8382e+02 4.8382e+02 1.0215e+03 8.1175e-05 8.1175e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +2.2000e-01 2.4000e+07 5.5315e+02 5.6174e+03 1.7179e-08 8.8110e-01 6.1877e-03 1.1271e-01 5.0414e+03 4.9180e+03 5.6694e+04 4.8521e+02 1.3086e+02 1.0211e+03 8.1572e-05 2.1495e-05 3.2929e-04 4.8901e-01 5.1099e-01 0.0000e+00 9.3267e-01 6.7327e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +2.4000e-01 2.3000e+07 5.5315e+02 5.5204e+03 1.7647e-08 8.7325e-01 1.4036e-02 1.1271e-01 4.9562e+03 4.7398e+03 5.6670e+04 4.8864e+02 1.2549e+02 1.0206e+03 8.2813e-05 2.1134e-05 3.2929e-04 4.7405e-01 5.2595e-01 0.0000e+00 9.3352e-01 6.6484e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +2.6000e-01 2.2000e+07 5.5315e+02 5.4225e+03 1.8162e-08 8.6581e-01 2.1475e-02 1.1271e-01 4.8719e+03 4.5585e+03 5.6647e+04 4.9205e+02 1.2016e+02 1.0202e+03 8.4061e-05 2.0788e-05 3.2929e-04 4.5873e-01 5.4127e-01 0.0000e+00 9.3426e-01 6.5741e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +2.8000e-01 2.1000e+07 5.5315e+02 5.3234e+03 1.8732e-08 8.5874e-01 2.8551e-02 1.1271e-01 4.7885e+03 4.3743e+03 5.6623e+04 4.9543e+02 1.1486e+02 1.0198e+03 8.5316e-05 2.0456e-05 3.2929e-04 4.4303e-01 5.5697e-01 0.0000e+00 9.3490e-01 6.5103e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +3.0000e-01 2.0000e+07 5.5315e+02 5.2230e+03 1.9369e-08 8.5198e-01 3.5311e-02 1.1271e-01 4.7061e+03 4.1872e+03 5.6600e+04 4.9880e+02 1.0961e+02 1.0194e+03 8.6578e-05 2.0139e-05 3.2929e-04 4.2693e-01 5.7307e-01 0.0000e+00 9.3542e-01 6.4579e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +3.2000e-01 1.9000e+07 5.5315e+02 5.1210e+03 2.0089e-08 8.4549e-01 4.1796e-02 1.1271e-01 4.6246e+03 3.9973e+03 5.6576e+04 5.0214e+02 1.0438e+02 1.0189e+03 8.7844e-05 1.9835e-05 3.2929e-04 4.1042e-01 5.8958e-01 0.0000e+00 9.3582e-01 6.4178e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +3.4000e-01 1.8000e+07 5.5315e+02 5.0171e+03 2.0911e-08 8.3924e-01 4.8049e-02 1.1271e-01 4.5439e+03 3.8047e+03 5.6553e+04 5.0547e+02 9.9197e+01 1.0185e+03 8.9115e-05 1.9544e-05 3.2929e-04 3.9347e-01 6.0653e-01 0.0000e+00 9.3609e-01 6.3914e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +3.6000e-01 1.7000e+07 5.5315e+02 4.9110e+03 2.1861e-08 8.3318e-01 5.4111e-02 1.1271e-01 4.4641e+03 3.6096e+03 5.6530e+04 5.0878e+02 9.4045e+01 1.0181e+03 9.0389e-05 1.9265e-05 3.2929e-04 3.7608e-01 6.2392e-01 0.0000e+00 9.3620e-01 6.3802e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +3.8000e-01 1.6000e+07 5.5315e+02 4.8022e+03 2.2976e-08 8.2726e-01 6.0026e-02 1.1271e-01 4.3852e+03 3.4119e+03 5.6506e+04 5.1208e+02 8.8928e+01 1.0177e+03 9.1664e-05 1.8998e-05 3.2929e-04 3.5822e-01 6.4178e-01 0.0000e+00 9.3614e-01 6.3864e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +4.0000e-01 1.5000e+07 5.5315e+02 4.6901e+03 2.4303e-08 8.2145e-01 6.5840e-02 1.1271e-01 4.3071e+03 3.2118e+03 5.6483e+04 5.1536e+02 8.3846e+01 1.0173e+03 9.2939e-05 1.8741e-05 3.2929e-04 3.3987e-01 6.6013e-01 0.0000e+00 9.3587e-01 6.4127e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +4.2000e-01 1.4000e+07 5.5315e+02 4.5740e+03 2.5910e-08 8.1568e-01 7.1608e-02 1.1271e-01 4.2299e+03 3.0095e+03 5.6459e+04 5.1864e+02 7.8800e+01 1.0168e+03 9.4210e-05 1.8494e-05 3.2929e-04 3.2102e-01 6.7898e-01 0.0000e+00 9.3537e-01 6.4625e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +4.4000e-01 1.3000e+07 5.5315e+02 4.4527e+03 2.7893e-08 8.0990e-01 7.7391e-02 1.1271e-01 4.1535e+03 2.8050e+03 5.6436e+04 5.2190e+02 7.3789e+01 1.0164e+03 9.5475e-05 1.8255e-05 3.2929e-04 3.0163e-01 6.9837e-01 0.0000e+00 9.3459e-01 6.5406e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +4.6000e-01 1.2000e+07 5.5315e+02 4.3251e+03 3.0392e-08 8.0402e-01 8.3265e-02 1.1271e-01 4.0779e+03 2.5984e+03 5.6412e+04 5.2516e+02 6.8814e+01 1.0160e+03 9.6732e-05 1.8022e-05 3.2929e-04 2.8170e-01 7.1830e-01 0.0000e+00 9.3347e-01 6.6532e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +4.8000e-01 1.1000e+07 5.5315e+02 4.1891e+03 3.3618e-08 7.9796e-01 8.9328e-02 1.1271e-01 4.0031e+03 2.3899e+03 5.6389e+04 5.2841e+02 6.3876e+01 1.0156e+03 9.7975e-05 1.7795e-05 3.2929e-04 2.6120e-01 7.3880e-01 0.0000e+00 9.3191e-01 6.8088e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +5.0000e-01 1.0000e+07 5.5315e+02 4.0424e+03 3.7901e-08 7.9158e-01 9.5711e-02 1.1271e-01 3.9291e+03 2.1796e+03 5.6366e+04 5.3165e+02 5.8976e+01 1.0151e+03 9.9202e-05 1.7569e-05 3.2929e-04 2.4011e-01 7.5989e-01 0.0000e+00 9.2981e-01 7.0194e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +5.2000e-01 9.6040e+06 5.5315e+02 3.9806e+03 3.9996e-08 7.8892e-01 9.8365e-02 1.1271e-01 3.9000e+03 2.0959e+03 5.6356e+04 5.3293e+02 5.7046e+01 1.0150e+03 9.9682e-05 1.7479e-05 3.2929e-04 2.3159e-01 7.6841e-01 0.0000e+00 9.2878e-01 7.1215e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +5.4000e-01 9.2080e+06 5.5315e+02 3.9162e+03 4.2382e-08 7.8618e-01 1.0111e-01 1.1271e-01 3.8710e+03 2.0119e+03 5.6347e+04 5.3421e+02 5.5122e+01 1.0148e+03 1.0016e-04 1.7389e-05 3.2929e-04 2.2297e-01 7.7703e-01 0.0000e+00 9.2764e-01 7.2363e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +5.6000e-01 8.8120e+06 5.5315e+02 3.8490e+03 4.5115e-08 7.8332e-01 1.0397e-01 1.1271e-01 3.8422e+03 1.9277e+03 5.6338e+04 5.3549e+02 5.3205e+01 1.0146e+03 1.0063e-04 1.7299e-05 3.2929e-04 2.1425e-01 7.8575e-01 0.0000e+00 9.2635e-01 7.3653e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +5.8000e-01 8.4160e+06 5.5315e+02 3.7785e+03 4.8266e-08 7.8033e-01 1.0696e-01 1.1271e-01 3.8135e+03 1.8432e+03 5.6329e+04 5.3677e+02 5.1293e+01 1.0145e+03 1.0110e-04 1.7207e-05 3.2929e-04 2.0543e-01 7.9457e-01 0.0000e+00 9.2489e-01 7.5105e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +6.0000e-01 8.0200e+06 5.5315e+02 3.7044e+03 5.1924e-08 7.7719e-01 1.1010e-01 1.1271e-01 3.7849e+03 1.7585e+03 5.6319e+04 5.3805e+02 4.9388e+01 1.0143e+03 1.0156e-04 1.7113e-05 3.2929e-04 1.9650e-01 8.0350e-01 0.0000e+00 9.2326e-01 7.6743e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +6.2000e-01 7.6240e+06 5.5315e+02 3.6260e+03 5.6202e-08 7.7386e-01 1.1343e-01 1.1271e-01 3.7565e+03 1.6736e+03 5.6310e+04 5.3933e+02 4.7489e+01 1.0141e+03 1.0202e-04 1.7018e-05 3.2929e-04 1.8748e-01 8.1252e-01 0.0000e+00 9.2141e-01 7.8594e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +6.4000e-01 7.2280e+06 5.5315e+02 3.5428e+03 6.1247e-08 7.7031e-01 1.1698e-01 1.1271e-01 3.7281e+03 1.5884e+03 5.6301e+04 5.4061e+02 4.5597e+01 1.0140e+03 1.0247e-04 1.6920e-05 3.2929e-04 1.7834e-01 8.2166e-01 0.0000e+00 9.1931e-01 8.0692e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +6.6000e-01 6.8320e+06 5.5315e+02 3.4539e+03 6.7247e-08 7.6650e-01 1.2079e-01 1.1271e-01 3.6999e+03 1.5031e+03 5.6292e+04 5.4188e+02 4.3711e+01 1.0138e+03 1.0292e-04 1.6818e-05 3.2929e-04 1.6910e-01 8.3090e-01 0.0000e+00 9.1692e-01 8.3080e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +6.8000e-01 6.4360e+06 5.5315e+02 3.3585e+03 7.4453e-08 7.6236e-01 1.2493e-01 1.1271e-01 3.6718e+03 1.4176e+03 5.6282e+04 5.4316e+02 4.1832e+01 1.0136e+03 1.0335e-04 1.6713e-05 3.2929e-04 1.5975e-01 8.4025e-01 0.0000e+00 9.1419e-01 8.5810e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +7.0000e-01 6.0400e+06 5.5315e+02 3.2554e+03 8.3199e-08 7.5783e-01 1.2946e-01 1.1271e-01 3.6438e+03 1.3319e+03 5.6273e+04 5.4443e+02 3.9959e+01 1.0135e+03 1.0379e-04 1.6602e-05 3.2929e-04 1.5029e-01 8.4971e-01 0.0000e+00 9.1105e-01 8.8948e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +7.2000e-01 5.6440e+06 5.5315e+02 3.1435e+03 9.3938e-08 7.5281e-01 1.3448e-01 1.1271e-01 3.6160e+03 1.2460e+03 5.6264e+04 5.4570e+02 3.8093e+01 1.0133e+03 1.0421e-04 1.6485e-05 3.2929e-04 1.4071e-01 8.5929e-01 0.0000e+00 9.0742e-01 9.2579e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +7.4000e-01 5.2480e+06 5.5315e+02 3.0210e+03 1.0730e-07 7.4718e-01 1.4011e-01 1.1271e-01 3.5883e+03 1.1600e+03 5.6255e+04 5.4698e+02 3.6233e+01 1.0131e+03 1.0463e-04 1.6359e-05 3.2929e-04 1.3102e-01 8.6898e-01 0.0000e+00 9.0319e-01 9.6812e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +7.6000e-01 4.8520e+06 5.5315e+02 2.8861e+03 1.2415e-07 7.4077e-01 1.4652e-01 1.1271e-01 3.5607e+03 1.0738e+03 5.6245e+04 5.4825e+02 3.4381e+01 1.0130e+03 1.0503e-04 1.6224e-05 3.2929e-04 1.2122e-01 8.7878e-01 0.0000e+00 8.9821e-01 1.0179e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +7.8000e-01 4.4560e+06 5.5315e+02 2.7364e+03 1.4578e-07 7.3336e-01 1.5393e-01 1.1271e-01 3.5332e+03 9.8750e+02 5.6236e+04 5.4952e+02 3.2535e+01 1.0128e+03 1.0543e-04 1.6076e-05 3.2929e-04 1.1129e-01 8.8871e-01 0.0000e+00 8.9228e-01 1.0772e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +8.0000e-01 4.0600e+06 5.5315e+02 2.5691e+03 1.7403e-07 7.2461e-01 1.6268e-01 1.1271e-01 3.5058e+03 9.0104e+02 5.6227e+04 5.5079e+02 3.0696e+01 1.0126e+03 1.0582e-04 1.5911e-05 3.2929e-04 1.0125e-01 8.9875e-01 0.0000e+00 8.8514e-01 1.1486e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +8.2000e-01 3.6640e+06 5.5315e+02 2.3811e+03 2.1172e-07 7.1406e-01 1.7323e-01 1.1271e-01 3.4786e+03 8.1446e+02 5.6218e+04 5.5206e+02 2.8865e+01 1.0125e+03 1.0619e-04 1.5725e-05 3.2929e-04 9.1082e-02 9.0892e-01 0.0000e+00 8.7640e-01 1.2360e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +8.4000e-01 3.2680e+06 5.5315e+02 2.1686e+03 2.6326e-07 7.0095e-01 1.8634e-01 1.1271e-01 3.4515e+03 7.2776e+02 5.6208e+04 5.5333e+02 2.7040e+01 1.0123e+03 1.0656e-04 1.5512e-05 3.2929e-04 8.0792e-02 9.1921e-01 0.0000e+00 8.6548e-01 1.3452e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +8.6000e-01 2.8720e+06 5.5315e+02 1.9276e+03 3.3583e-07 6.8411e-01 2.0318e-01 1.1271e-01 3.4245e+03 6.4096e+02 5.6199e+04 5.5460e+02 2.5222e+01 1.0121e+03 1.0691e-04 1.5261e-05 3.2929e-04 7.0376e-02 9.2962e-01 0.0000e+00 8.5148e-01 1.4852e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +8.8000e-01 2.4760e+06 5.5315e+02 1.6549e+03 4.4172e-07 6.6146e-01 2.2583e-01 1.1271e-01 3.3976e+03 5.5407e+02 5.6190e+04 5.5587e+02 2.3412e+01 1.0120e+03 1.0725e-04 1.4958e-05 3.2929e-04 5.9832e-02 9.4017e-01 0.0000e+00 8.3298e-01 1.6702e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +9.0000e-01 2.0800e+06 5.5315e+02 1.3487e+03 6.0350e-07 6.2908e-01 2.5821e-01 1.1271e-01 3.3708e+03 4.6709e+02 5.6181e+04 5.5714e+02 2.1609e+01 1.0118e+03 1.0757e-04 1.4583e-05 3.2929e-04 4.9159e-02 9.5084e-01 0.0000e+00 8.0745e-01 1.9255e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +9.2000e-01 1.6840e+06 5.5315e+02 1.0126e+03 8.6704e-07 5.7850e-01 3.0879e-01 1.1271e-01 3.3442e+03 3.8004e+02 5.6172e+04 5.5841e+02 1.9813e+01 1.0116e+03 1.0788e-04 1.4097e-05 3.2929e-04 3.8354e-02 9.6165e-01 0.0000e+00 7.7004e-01 2.2996e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +9.4000e-01 1.2880e+06 5.5315e+02 6.6052e+02 1.3428e-06 4.8744e-01 3.9985e-01 1.1271e-01 3.3176e+03 2.9292e+02 5.6162e+04 5.5968e+02 1.8025e+01 1.0115e+03 1.0818e-04 1.3437e-05 3.2929e-04 2.7415e-02 9.7259e-01 0.0000e+00 7.1015e-01 2.8985e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +9.6000e-01 8.9200e+05 5.5315e+02 3.2534e+02 2.3985e-06 2.7229e-01 6.1500e-01 1.1271e-01 3.2912e+03 2.0575e+02 5.6153e+04 5.6094e+02 1.6244e+01 1.0113e+03 1.0845e-04 1.2477e-05 3.2929e-04 1.6340e-02 9.8366e-01 0.0000e+00 5.9920e-01 4.0080e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +9.8000e-01 4.9600e+05 5.5315e+02 1.2779e+02 2.1211e-06 8.8729e-01 0.0000e+00 1.1271e-01 1.1342e+02 1.1342e+02 5.6144e+04 1.0720e+01 1.0720e+01 1.0112e+03 1.1716e-05 1.1716e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.0000e+00 1.0000e+05 5.5315e+02 2.4746e+01 1.0099e-05 8.8729e-01 0.0000e+00 1.1271e-01 2.1958e+01 2.1958e+01 5.6135e+04 2.0755e+00 2.0755e+00 1.0110e+03 1.1560e-05 1.1560e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.0200e+00 4.9600e+05 5.8315e+02 1.2012e+02 2.1009e-06 8.8729e-01 0.0000e+00 1.1271e-01 1.0660e+02 1.0660e+02 5.6144e+04 1.0076e+01 1.0076e+01 1.0112e+03 1.2221e-05 1.2221e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.0400e+00 8.9200e+05 5.8315e+02 2.2354e+02 1.2093e-06 8.8729e-01 0.0000e+00 1.1271e-01 1.9843e+02 1.9843e+02 5.6153e+04 1.8756e+01 1.8756e+01 1.0113e+03 1.2405e-05 1.2405e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.0600e+00 1.2880e+06 5.8315e+02 3.4774e+02 1.7667e-06 4.3627e-02 8.4366e-01 1.1271e-01 3.0621e+03 2.9504e+02 5.6162e+04 5.1804e+02 2.7265e+01 1.0115e+03 8.2939e-05 1.2695e-05 3.2929e-04 2.4366e-02 9.7563e-01 0.0000e+00 5.1344e-01 4.8656e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.0800e+00 1.6840e+06 5.8315e+02 6.0704e+02 1.1265e-06 3.0186e-01 5.8543e-01 1.1271e-01 3.0890e+03 3.7828e+02 5.6172e+04 5.1649e+02 2.9209e+01 1.0116e+03 8.2756e-05 1.3409e-05 3.2929e-04 3.6932e-02 9.6307e-01 0.0000e+00 6.1026e-01 3.8974e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.1000e+00 2.0800e+06 5.8315e+02 8.8181e+02 7.8965e-07 4.2838e-01 4.5891e-01 1.1271e-01 3.1161e+03 4.6143e+02 5.6181e+04 5.1495e+02 3.1163e+01 1.0118e+03 8.2556e-05 1.3965e-05 3.2929e-04 4.9339e-02 9.5066e-01 0.0000e+00 6.7192e-01 3.2808e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.1200e+00 2.4760e+06 5.8315e+02 1.1539e+03 5.8364e-07 5.0379e-01 3.8350e-01 1.1271e-01 3.1432e+03 5.4447e+02 5.6190e+04 5.1339e+02 3.3128e+01 1.0120e+03 8.2339e-05 1.4414e-05 3.2929e-04 6.1589e-02 9.3841e-01 0.0000e+00 7.1454e-01 2.8546e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.1400e+00 2.8720e+06 5.8315e+02 1.4127e+03 4.4720e-07 5.5409e-01 3.3320e-01 1.1271e-01 3.1706e+03 6.2740e+02 5.6199e+04 5.1184e+02 3.5105e+01 1.0121e+03 8.2106e-05 1.4788e-05 3.2929e-04 7.3685e-02 9.2632e-01 0.0000e+00 7.4571e-01 2.5429e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.1600e+00 3.2680e+06 5.8315e+02 1.6532e+03 3.5221e-07 5.9019e-01 2.9710e-01 1.1271e-01 3.1980e+03 7.1021e+02 5.6208e+04 5.1028e+02 3.7093e+01 1.0123e+03 8.1860e-05 1.5107e-05 3.2929e-04 8.5630e-02 9.1437e-01 0.0000e+00 7.6944e-01 2.3056e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.1800e+00 3.6640e+06 5.8315e+02 1.8738e+03 2.8372e-07 6.1751e-01 2.6978e-01 1.1271e-01 3.2257e+03 7.9289e+02 5.6218e+04 5.0871e+02 3.9092e+01 1.0125e+03 8.1599e-05 1.5383e-05 3.2929e-04 9.7426e-02 9.0257e-01 0.0000e+00 7.8808e-01 2.1192e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.2000e+00 4.0600e+06 5.8315e+02 2.0746e+03 2.3298e-07 6.3902e-01 2.4827e-01 1.1271e-01 3.2534e+03 8.7543e+02 5.6227e+04 5.0714e+02 4.1103e+01 1.0126e+03 8.1325e-05 1.5628e-05 3.2929e-04 1.0908e-01 8.9092e-01 0.0000e+00 8.0306e-01 1.9694e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.2200e+00 4.4560e+06 5.8315e+02 2.2571e+03 1.9456e-07 6.5648e-01 2.3081e-01 1.1271e-01 3.2813e+03 9.5783e+02 5.6236e+04 5.0557e+02 4.3124e+01 1.0128e+03 8.1039e-05 1.5847e-05 3.2929e-04 1.2058e-01 8.7942e-01 0.0000e+00 8.1535e-01 1.8465e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.2400e+00 4.8520e+06 5.8315e+02 2.4230e+03 1.6492e-07 6.7102e-01 2.1627e-01 1.1271e-01 3.3094e+03 1.0401e+03 5.6245e+04 5.0399e+02 4.5158e+01 1.0130e+03 8.0742e-05 1.6047e-05 3.2929e-04 1.3195e-01 8.6805e-01 0.0000e+00 8.2557e-01 1.7443e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.2600e+00 5.2480e+06 5.8315e+02 2.5742e+03 1.4168e-07 6.8338e-01 2.0390e-01 1.1271e-01 3.3376e+03 1.1221e+03 5.6255e+04 5.0240e+02 4.7203e+01 1.0131e+03 8.0433e-05 1.6230e-05 3.2929e-04 1.4318e-01 8.5682e-01 0.0000e+00 8.3419e-01 1.6581e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.2800e+00 5.6440e+06 5.8315e+02 2.7125e+03 1.2319e-07 6.9409e-01 1.9320e-01 1.1271e-01 3.3660e+03 1.2040e+03 5.6264e+04 5.0081e+02 4.9259e+01 1.0133e+03 8.0115e-05 1.6401e-05 3.2929e-04 1.5427e-01 8.4573e-01 0.0000e+00 8.4154e-01 1.5846e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.3000e+00 6.0400e+06 5.8315e+02 2.8394e+03 1.0830e-07 7.0351e-01 1.8378e-01 1.1271e-01 3.3945e+03 1.2857e+03 5.6273e+04 4.9921e+02 5.1327e+01 1.0135e+03 7.9786e-05 1.6561e-05 3.2929e-04 1.6523e-01 8.3477e-01 0.0000e+00 8.4786e-01 1.5214e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.3200e+00 6.4360e+06 5.8315e+02 2.9565e+03 9.6164e-08 7.1191e-01 1.7538e-01 1.1271e-01 3.4232e+03 1.3673e+03 5.6282e+04 4.9761e+02 5.3407e+01 1.0136e+03 7.9448e-05 1.6713e-05 3.2929e-04 1.7606e-01 8.2394e-01 0.0000e+00 8.5333e-01 1.4667e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.3400e+00 6.8320e+06 5.8315e+02 3.0649e+03 8.6168e-08 7.1948e-01 1.6781e-01 1.1271e-01 3.4520e+03 1.4486e+03 5.6292e+04 4.9600e+02 5.5498e+01 1.0138e+03 7.9102e-05 1.6858e-05 3.2929e-04 1.8676e-01 8.1324e-01 0.0000e+00 8.5811e-01 1.4189e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.3600e+00 7.2280e+06 5.8315e+02 3.1660e+03 7.7854e-08 7.2640e-01 1.6089e-01 1.1271e-01 3.4809e+03 1.5297e+03 5.6301e+04 4.9439e+02 5.7601e+01 1.0140e+03 7.8748e-05 1.6997e-05 3.2929e-04 1.9733e-01 8.0267e-01 0.0000e+00 8.6229e-01 1.3771e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.3800e+00 7.6240e+06 5.8315e+02 3.2605e+03 7.0879e-08 7.3276e-01 1.5453e-01 1.1271e-01 3.5101e+03 1.6106e+03 5.6310e+04 4.9277e+02 5.9717e+01 1.0141e+03 7.8386e-05 1.7132e-05 3.2929e-04 2.0778e-01 7.9222e-01 0.0000e+00 8.6597e-01 1.3403e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.4000e+00 8.0200e+06 5.8315e+02 3.3492e+03 6.4979e-08 7.3868e-01 1.4861e-01 1.1271e-01 3.5394e+03 1.6913e+03 5.6319e+04 4.9114e+02 6.1844e+01 1.0143e+03 7.8016e-05 1.7263e-05 3.2929e-04 2.1811e-01 7.8189e-01 0.0000e+00 8.6923e-01 1.3077e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.4200e+00 8.4160e+06 5.8315e+02 3.4330e+03 5.9949e-08 7.4422e-01 1.4307e-01 1.1271e-01 3.5688e+03 1.7717e+03 5.6329e+04 4.8951e+02 6.3983e+01 1.0145e+03 7.7640e-05 1.7391e-05 3.2929e-04 2.2831e-01 7.7169e-01 0.0000e+00 8.7211e-01 1.2789e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.4400e+00 8.8120e+06 5.8315e+02 3.5124e+03 5.5632e-08 7.4946e-01 1.3783e-01 1.1271e-01 3.5984e+03 1.8519e+03 5.6338e+04 4.8787e+02 6.6135e+01 1.0146e+03 7.7258e-05 1.7517e-05 3.2929e-04 2.3840e-01 7.6160e-01 0.0000e+00 8.7467e-01 1.2533e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.4600e+00 9.2080e+06 5.8315e+02 3.5880e+03 5.1902e-08 7.5444e-01 1.3285e-01 1.1271e-01 3.6282e+03 1.9319e+03 5.6347e+04 4.8623e+02 6.8299e+01 1.0148e+03 7.6869e-05 1.7641e-05 3.2929e-04 2.4838e-01 7.5162e-01 0.0000e+00 8.7695e-01 1.2305e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.4800e+00 9.6040e+06 5.8315e+02 3.6601e+03 4.8659e-08 7.5920e-01 1.2809e-01 1.1271e-01 3.6581e+03 2.0115e+03 5.6356e+04 4.8457e+02 7.0476e+01 1.0150e+03 7.6475e-05 1.7765e-05 3.2929e-04 2.5824e-01 7.4176e-01 0.0000e+00 8.7898e-01 1.2102e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.5000e+00 1.0000e+07 5.8315e+02 3.7291e+03 4.5824e-08 7.6378e-01 1.2351e-01 1.1271e-01 3.6882e+03 2.0910e+03 5.6366e+04 4.8291e+02 7.2665e+01 1.0151e+03 7.6075e-05 1.7888e-05 3.2929e-04 2.6799e-01 7.3201e-01 0.0000e+00 8.8078e-01 1.1922e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.5200e+00 1.1000e+07 5.8315e+02 3.8922e+03 4.0073e-08 7.7474e-01 1.1254e-01 1.1271e-01 3.7648e+03 2.2903e+03 5.6389e+04 4.7869e+02 7.8250e+01 1.0156e+03 7.5043e-05 1.8197e-05 3.2929e-04 2.9214e-01 7.0786e-01 0.0000e+00 8.8451e-01 1.1549e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.5400e+00 1.2000e+07 5.8315e+02 4.0423e+03 3.5793e-08 7.8515e-01 1.0214e-01 1.1271e-01 3.8426e+03 2.4877e+03 5.6412e+04 4.7441e+02 8.3920e+01 1.0160e+03 7.3984e-05 1.8510e-05 3.2929e-04 3.1562e-01 6.8438e-01 0.0000e+00 8.8727e-01 1.1273e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.5600e+00 1.3000e+07 5.8315e+02 4.1824e+03 3.2518e-08 7.9527e-01 9.2016e-02 1.1271e-01 3.9214e+03 2.6831e+03 5.6436e+04 4.7007e+02 8.9678e+01 1.0164e+03 7.2901e-05 1.8830e-05 3.2929e-04 3.3847e-01 6.6153e-01 0.0000e+00 8.8925e-01 1.1075e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.5800e+00 1.4000e+07 5.8315e+02 4.3149e+03 2.9952e-08 8.0535e-01 8.1944e-02 1.1271e-01 4.0014e+03 2.8763e+03 5.6459e+04 4.6568e+02 9.5527e+01 1.0168e+03 7.1797e-05 1.9161e-05 3.2929e-04 3.6073e-01 6.3927e-01 0.0000e+00 8.9060e-01 1.0940e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.6000e+00 1.5000e+07 5.8315e+02 4.4414e+03 2.7899e-08 8.1555e-01 7.1735e-02 1.1271e-01 4.0824e+03 3.0673e+03 5.6483e+04 4.6121e+02 1.0147e+02 1.0173e+03 7.0674e-05 1.9505e-05 3.2929e-04 3.8242e-01 6.1758e-01 0.0000e+00 8.9143e-01 1.0857e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.6200e+00 1.6000e+07 5.8315e+02 4.5631e+03 2.6223e-08 8.2607e-01 6.1222e-02 1.1271e-01 4.1647e+03 3.2559e+03 5.6506e+04 4.5668e+02 1.0752e+02 1.0177e+03 6.9535e-05 1.9866e-05 3.2929e-04 4.0358e-01 5.9642e-01 0.0000e+00 8.9180e-01 1.0820e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.6400e+00 1.7000e+07 5.8315e+02 4.6810e+03 2.4832e-08 8.3705e-01 5.0240e-02 1.1271e-01 4.2481e+03 3.4420e+03 5.6530e+04 4.5206e+02 1.1367e+02 1.0181e+03 6.8382e-05 2.0245e-05 3.2929e-04 4.2423e-01 5.7577e-01 0.0000e+00 8.9179e-01 1.0821e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.6600e+00 1.8000e+07 5.8315e+02 4.7958e+03 2.3660e-08 8.4867e-01 3.8624e-02 1.1271e-01 4.3328e+03 3.6254e+03 5.6553e+04 4.4735e+02 1.1994e+02 1.0185e+03 6.7214e-05 2.0645e-05 3.2929e-04 4.4441e-01 5.5559e-01 0.0000e+00 8.9142e-01 1.0858e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.6800e+00 1.9000e+07 5.8315e+02 4.9081e+03 2.2657e-08 8.6110e-01 2.6192e-02 1.1271e-01 4.4188e+03 3.8061e+03 5.6576e+04 4.4254e+02 1.2633e+02 1.0189e+03 6.6034e-05 2.1068e-05 3.2929e-04 4.6414e-01 5.3586e-01 0.0000e+00 8.9072e-01 1.0928e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.7000e+00 2.0000e+07 5.8315e+02 5.0183e+03 2.1788e-08 8.7455e-01 1.2742e-02 1.1271e-01 4.5061e+03 3.9839e+03 5.6600e+04 4.3761e+02 1.3286e+02 1.0194e+03 6.4840e-05 2.1517e-05 3.2929e-04 4.8347e-01 5.1653e-01 0.0000e+00 8.8972e-01 1.1028e-01 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.7200e+00 2.1000e+07 5.8315e+02 5.1233e+03 1.5507e-08 8.8729e-01 0.0000e+00 1.1271e-01 4.5927e+03 4.5927e+03 5.6623e+04 4.3410e+02 4.3410e+02 1.0198e+03 6.4038e-05 6.4038e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.7400e+00 2.2000e+07 5.8315e+02 5.2003e+03 1.4379e-08 8.8729e-01 0.0000e+00 1.1271e-01 4.6625e+03 4.6625e+03 5.6647e+04 4.4069e+02 4.4069e+02 1.0202e+03 6.5972e-05 6.5972e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.7600e+00 2.3000e+07 5.8315e+02 5.2730e+03 1.3388e-08 8.8729e-01 0.0000e+00 1.1271e-01 4.7283e+03 4.7283e+03 5.6670e+04 4.4691e+02 4.4691e+02 1.0206e+03 6.7891e-05 6.7891e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.7800e+00 2.4000e+07 5.8315e+02 5.3417e+03 1.2511e-08 8.8729e-01 0.0000e+00 1.1271e-01 4.7905e+03 4.7905e+03 5.6694e+04 4.5280e+02 4.5280e+02 1.0211e+03 6.9797e-05 6.9797e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.8000e+00 2.5000e+07 5.8315e+02 5.4068e+03 1.1730e-08 8.8729e-01 0.0000e+00 1.1271e-01 4.8495e+03 4.8495e+03 5.6717e+04 4.5837e+02 4.5837e+02 1.0215e+03 7.1695e-05 7.1695e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.8200e+00 2.6000e+07 5.8315e+02 5.4686e+03 1.1031e-08 8.8729e-01 0.0000e+00 1.1271e-01 4.9055e+03 4.9055e+03 5.6741e+04 4.6367e+02 4.6367e+02 1.0219e+03 7.3587e-05 7.3587e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.8400e+00 2.7000e+07 5.8315e+02 5.5275e+03 1.0402e-08 8.8729e-01 0.0000e+00 1.1271e-01 4.9589e+03 4.9589e+03 5.6765e+04 4.6872e+02 4.6872e+02 1.0223e+03 7.5475e-05 7.5475e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.8600e+00 2.8000e+07 5.8315e+02 5.5837e+03 9.8334e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.0099e+03 5.0099e+03 5.6788e+04 4.7353e+02 4.7353e+02 1.0228e+03 7.7360e-05 7.7360e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.8800e+00 2.9000e+07 5.8315e+02 5.6374e+03 9.3169e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.0586e+03 5.0586e+03 5.6812e+04 4.7814e+02 4.7814e+02 1.0232e+03 7.9245e-05 7.9245e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.9000e+00 3.0000e+07 5.8315e+02 5.6888e+03 8.8460e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.1052e+03 5.1052e+03 5.6835e+04 4.8254e+02 4.8254e+02 1.0236e+03 8.1131e-05 8.1131e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.9200e+00 3.1000e+07 5.8315e+02 5.7381e+03 8.4151e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.1499e+03 5.1499e+03 5.6859e+04 4.8677e+02 4.8677e+02 1.0240e+03 8.3018e-05 8.3018e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.9400e+00 3.2000e+07 5.8315e+02 5.7854e+03 8.0196e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.1929e+03 5.1929e+03 5.6882e+04 4.9083e+02 4.9083e+02 1.0245e+03 8.4907e-05 8.4907e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.9600e+00 3.3000e+07 5.8315e+02 5.8309e+03 7.6552e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.2342e+03 5.2342e+03 5.6906e+04 4.9473e+02 4.9473e+02 1.0249e+03 8.6800e-05 8.6800e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +1.9800e+00 3.4000e+07 5.8315e+02 5.8747e+03 7.3187e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.2739e+03 5.2739e+03 5.6930e+04 4.9849e+02 4.9849e+02 1.0253e+03 8.8696e-05 8.8696e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00 +2.0000e+00 3.5000e+07 5.8315e+02 5.9169e+03 7.0071e-09 8.8729e-01 0.0000e+00 1.1271e-01 5.3122e+03 5.3122e+03 5.6953e+04 5.0211e+02 5.0211e+02 1.0257e+03 9.0597e-05 9.0597e-05 3.2929e-04 5.0000e-01 5.0000e-01 0.0000e+00 9.1515e-01 8.4850e-02 0.0000e+00 0.0000e+00 0.0000e+00 1.0000e+00