- T_c0Critical temperature, K
Default:0
C++ Type:double
Description:Critical temperature, K
- allow_imperfect_jacobiansFalsetrue to allow unimplemented property derivative terms to be set to zero for the AD API
Default:False
C++ Type:bool
Description:true to allow unimplemented property derivative terms to be set to zero for the AD API
- e_c0Internal energy at the critical point, J/kg
Default:0
C++ Type:double
Description:Internal energy at the critical point, J/kg
- execute_onTIMESTEP_ENDThe list of flag(s) indicating when this object should be executed, the available options include NONE, INITIAL, LINEAR, NONLINEAR, TIMESTEP_END, TIMESTEP_BEGIN, FINAL, CUSTOM.
Default:TIMESTEP_END
C++ Type:ExecFlagEnum
Description:The list of flag(s) indicating when this object should be executed, the available options include NONE, INITIAL, LINEAR, NONLINEAR, TIMESTEP_END, TIMESTEP_BEGIN, FINAL, CUSTOM.
- fp_typesingle-phase-fpType of the fluid property object
Default:single-phase-fp
C++ Type:FPType
Description:Type of the fluid property object
- gamma1.4gamma value (cp/cv)
Default:1.4
C++ Type:double
Description:gamma value (cp/cv)
- k0.02568Thermal conductivity, W/(m-K)
Default:0.02568
C++ Type:double
Description:Thermal conductivity, W/(m-K)
- molar_mass0.029Constant molar mass of the fluid (kg/mol)
Default:0.029
C++ Type:double
Description:Constant molar mass of the fluid (kg/mol)
- mu1.823e-05Dynamic viscosity, Pa.s
Default:1.823e-05
C++ Type:double
Description:Dynamic viscosity, Pa.s
- rho_c0Critical density, kg/m3
Default:0
C++ Type:double
Description:Critical density, kg/m3
IdealGasFluidProperties
Fluid properties for an ideal gas
A simple formulation that is suitable for ideal gases, where properties are derived from the ideal gas law
Temperature is calculated using the internal energy of an ideal gas
Input Parameters
- allow_duplicate_execution_on_initialFalseIn the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
Default:False
C++ Type:bool
Options:
Description:In the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
- control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Options:
Description:Adds user-defined labels for accessing object parameters via control logic.
- enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Options:
Description:Set the enabled status of the MooseObject.
- force_preauxFalseForces the GeneralUserObject to be executed in PREAUX
Default:False
C++ Type:bool
Options:
Description:Forces the GeneralUserObject to be executed in PREAUX
- force_preicFalseForces the GeneralUserObject to be executed in PREIC during initial setup
Default:False
C++ Type:bool
Options:
Description:Forces the GeneralUserObject to be executed in PREIC during initial setup
- use_displaced_meshFalseWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Default:False
C++ Type:bool
Options:
Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Advanced Parameters
Input Files
- (modules/fluid_properties/test/tests/ideal_gas/test2.i)
- (modules/fluid_properties/test/tests/fp_interrogator/1ph.rho_p.i)
- (modules/navier_stokes/test/tests/bump/bump.i)
- (modules/fluid_properties/fp_interrogator/fp_interrogator.i)
- (modules/fluid_properties/test/tests/ideal_gas/test.i)
- (modules/fluid_properties/test/tests/fp_interrogator/2ph_ncg_p_T.i)
- (modules/fluid_properties/test/tests/fp_interrogator/2ph.p_T.i)
- (modules/porous_flow/test/tests/fluids/ideal_gas.i)
- (modules/fluid_properties/test/tests/fp_interrogator/1ph.p_T.i)
- (modules/fluid_properties/test/tests/fp_interrogator/1ph.rho_e.i)
- (modules/fluid_properties/test/tests/fp_interrogator/2ph.T.i)
- (modules/fluid_properties/test/tests/fp_interrogator/1ph.rho_rhou_rhoE.i)
- (modules/fluid_properties/test/tests/fp_interrogator/err.no_params.i)
- (modules/porous_flow/test/tests/newton_cooling/nc08.i)
- (modules/fluid_properties/test/tests/ics/rho_vapor_mixture_from_pressure_temperature/test.i)
- (modules/fluid_properties/test/tests/fp_interrogator/2ph.p.i)
- (modules/porous_flow/test/tests/jacobian/esbc02.i)
- (modules/fluid_properties/test/tests/functions/saturation_temperature_function/saturation_temperature_function.i)
- (modules/fluid_properties/test/tests/two_phase_fluid_properties_independent/test.i)
- (modules/porous_flow/test/tests/jacobian/esbc01.i)
- (modules/fluid_properties/test/tests/fp_interrogator/vapor_mixture_rho_e.i)
- (modules/navier_stokes/test/tests/step/step.i)
- (modules/fluid_properties/test/tests/functions/saturation_pressure_function/saturation_pressure_function.i)
(modules/fluid_properties/test/tests/ideal_gas/test2.i)
# Test IdealGasFluidPropertiesFluidProperties using pressure and temperature
# Use values for Oxygen at 1 MPa and 350 K from NIST chemistry webook
#
# Input values:
# Cv = 669.8e J/kg/K
# Cp = 938.75 J/kg/K
# M = 31.9988e-3 kg/mol
#
# Expected output:
# density = 10.99591793 kg/m^3
# internal energy = 234.43e3 J/kg
# enthalpy = 328.5625e3 J/kg
# speed of sound = 357.0151605 m/s
[Mesh]
type = GeneratedMesh
dim = 2
nx = 1
ny = 1
[]
[Variables]
[./dummy]
[../]
[]
[AuxVariables]
[./pressure]
family = MONOMIAL
order = CONSTANT
initial_condition = 1e6
[../]
[./temperature]
family = MONOMIAL
order = CONSTANT
initial_condition = 350
[../]
[./density]
family = MONOMIAL
order = CONSTANT
[../]
[./viscosity]
family = MONOMIAL
order = CONSTANT
[../]
[./cp]
family = MONOMIAL
order = CONSTANT
[../]
[./cv]
family = MONOMIAL
order = CONSTANT
[../]
[./internal_energy]
family = MONOMIAL
order = CONSTANT
[../]
[./enthalpy]
family = MONOMIAL
order = CONSTANT
[../]
[./entropy]
family = MONOMIAL
order = CONSTANT
[../]
[./thermal_cond]
family = MONOMIAL
order = CONSTANT
[../]
[./c]
family = MONOMIAL
order = CONSTANT
[../]
[]
[AuxKernels]
[./density]
type = MaterialRealAux
variable = density
property = density
[../]
[./viscosity]
type = MaterialRealAux
variable = viscosity
property = viscosity
[../]
[./cp]
type = MaterialRealAux
variable = cp
property = cp
[../]
[./cv]
type = MaterialRealAux
variable = cv
property = cv
[../]
[./e]
type = MaterialRealAux
variable = internal_energy
property = e
[../]
[./enthalpy]
type = MaterialRealAux
variable = enthalpy
property = h
[../]
[./entropy]
type = MaterialRealAux
variable = entropy
property = s
[../]
[./thermal_cond]
type = MaterialRealAux
variable = thermal_cond
property = k
[../]
[./c]
type = MaterialRealAux
variable = c
property = c
[../]
[]
[Modules]
[./FluidProperties]
[./idealgas]
type = IdealGasFluidProperties
gamma = 1.401537772469394
molar_mass = 0.0319988
[../]
[]
[]
[Materials]
[./fp_mat]
type = FluidPropertiesMaterialPT
pressure = pressure
temperature = temperature
fp = idealgas
[../]
[]
[Kernels]
[./diff]
type = Diffusion
variable = dummy
[../]
[]
[Executioner]
type = Steady
solve_type = NEWTON
[]
[Outputs]
exodus = true
[]
(modules/fluid_properties/test/tests/fp_interrogator/1ph.rho_p.i)
[FluidPropertiesInterrogator]
fp = fp
rho = 1
p = 1e5
[]
[Modules]
[./FluidProperties]
[./fp]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 0.02900055737704918
mu = 1.823e-05
k = 0.02568
[../]
[../]
[]
(modules/navier_stokes/test/tests/bump/bump.i)
# Euler flow of an ideal gas over a Gaussian "bump".
#
# The inlet is a stagnation pressure and temperature BC which
# corresponds to subsonic (M=0.5) flow with a static pressure of 1 atm
# and static temperature of 300K. The outlet consists of a
# weakly-imposed static pressure BC of 1 atm. The top and bottom
# walls of the channel weakly impose the "no normal flow" BC. The
# problem is initialized with freestream flow throughout the domain.
# Although this initial condition is less physically realistic, it
# helps the problem reach steady state more quickly.
#
# There is a sequence of uniformly-refined, geometry-fitted meshes
# from Yidong Xia available for solving this classical subsonic test
# problem (see the Mesh block below). A coarse grid is used for the
# actual regression test, but changing one line in the Mesh block is
# sufficient to run this problem with different meshes. An
# entropy-based error estimate is also provided, and can be used to
# demonstrate convergence of the numerical solution (since the true
# solution should produce zero entropy). The error should converge at
# second-order in this norm.
[Mesh]
# Bi-Linear elements
# file = SmoothBump_quad_ref1_Q1.msh # 84 elems, 65 nodes
# file = SmoothBump_quad_ref2_Q1.msh # 192 elems, 225 nodes
# file = SmoothBump_quad_ref3_Q1.msh # 768 elems, 833 nodes
# file = SmoothBump_quad_ref4_Q1.msh # 3072 elems, 3201 nodes
# file = SmoothBump_quad_ref5_Q1.msh # 12288 elems, 12545 nodes
# Bi-Quadratic elements
# file = SmoothBump_quad_ref0_Q2.msh # 32 elems, 65 nodes
# file = SmoothBump_quad_ref1_Q2.msh # 84 elems, 225 nodes
file = SmoothBump_quad_ref2_Q2.msh # 260 elems, 833 nodes
# file = SmoothBump_quad_ref3_Q2.msh # 900 elems, 3201 nodes
# file = SmoothBump_quad_ref4_Q2.msh # 3332 elems, 12545 nodes
# file = SmoothBump_quad_ref5_Q2.msh # 12804 elems, 49665 nodes
[]
[Modules]
[FluidProperties]
[ideal_gas]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 0.02897024320557491
[]
[]
[CompressibleNavierStokes]
# steady-state or transient
equation_type = transient
# fluid
fluid_properties = ideal_gas
# boundary conditions
stagnation_boundary = 1
stagnation_pressure = 120192.995549849 # Pa, Mach=0.5 at 1 atm
stagnation_temperature = 315 # K, Mach=0.5 at 1 atm
stagnation_flow_direction = '1 0'
no_penetration_boundary = '3 4'
static_pressure_boundary = 2
static_pressure = 101325 # Pa
# variable types, scalings and initial conditions
family = LAGRANGE
order = FIRST
total_energy_scaling = 9.869232667160121e-6
initial_pressure = 101325.
initial_temperature = 300.
initial_velocity = '173.594354746921 0 0' # Mach 0.5: = 0.5*sqrt(gamma*R*T)
[]
[]
[Materials]
[fluid]
type = Air
block = 0 # 'MeshInterior'
rho = rho
rhou = rhou
rhov = rhov
rhoE = rhoE
vel_x = vel_x
vel_y = vel_y
temperature = temperature
enthalpy = enthalpy
# This value is not used in the Euler equations, but it *is* used
# by the stabilization parameter computation, which it decreases
# the amount of artificial viscosity added, so it's best to use a
# realistic value.
dynamic_viscosity = 0.0
fluid_properties = ideal_gas
[]
[]
[Postprocessors]
[entropy_error]
type = NSEntropyError
execute_on = 'initial timestep_end'
block = 0
rho_infty = 1.1768292682926829
p_infty = 101325
rho = rho
pressure = p
fluid_properties = ideal_gas
[]
[]
[Preconditioning]
[SMP_PJFNK]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
dt = 5.e-5
dtmin = 1.e-5
start_time = 0.0
num_steps = 10
nl_rel_tol = 1e-9
nl_max_its = 5
l_tol = 1e-4
l_max_its = 100
# We use trapezoidal quadrature. This improves stability by
# mimicking the "group variable" discretization approach.
[Quadrature]
type = TRAP
order = FIRST
[]
[]
[Outputs]
interval = 1
exodus = true
[]
(modules/fluid_properties/fp_interrogator/fp_interrogator.i)
# The parameters in this block are used to specify the thermodynamic state
# at which to query the fluid properties package
[FluidPropertiesInterrogator]
fp = fp
p = 1e5
T = 300
vel = 10
[]
# The fluid properties (equation of state) to query is defined here
[Modules]
[./FluidProperties]
[./fp]
type = IdealGasFluidProperties
[../]
[../]
[]
(modules/fluid_properties/test/tests/ideal_gas/test.i)
[Mesh]
type = GeneratedMesh
dim = 2
nx = 2
ny = 2
elem_type = QUAD4
[]
[Functions]
[./f_fn]
type = ParsedFunction
value = -4
[../]
[./bc_fn]
type = ParsedFunction
value = 'x*x+y*y'
[../]
[]
[Variables]
[./u]
[../]
[]
[AuxVariables]
[./e]
initial_condition = 6232.5
[../]
[./v]
initial_condition = 0.02493
[../]
[./p]
family = MONOMIAL
order = CONSTANT
[../]
[./T]
family = MONOMIAL
order = CONSTANT
[../]
[./cp]
family = MONOMIAL
order = CONSTANT
[../]
[./cv]
family = MONOMIAL
order = CONSTANT
[../]
[./c]
family = MONOMIAL
order = CONSTANT
[../]
[./mu]
family = MONOMIAL
order = CONSTANT
[../]
[./k]
family = MONOMIAL
order = CONSTANT
[../]
[./g]
family = MONOMIAL
order = CONSTANT
[../]
[]
[AuxKernels]
[./p]
type = MaterialRealAux
variable = p
property = pressure
[../]
[./T]
type = MaterialRealAux
variable = T
property = temperature
[../]
[./cp]
type = MaterialRealAux
variable = cp
property = cp
[../]
[./cv]
type = MaterialRealAux
variable = cv
property = cv
[../]
[./c]
type = MaterialRealAux
variable = c
property = c
[../]
[./mu]
type = MaterialRealAux
variable = mu
property = mu
[../]
[./k]
type = MaterialRealAux
variable = k
property = k
[../]
[./g]
type = MaterialRealAux
variable = g
property = g
[../]
[]
[Modules]
[./FluidProperties]
[./ideal_gas]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 1.000536678700361
[../]
[]
[]
[Materials]
[./fp_mat]
type = FluidPropertiesMaterial
e = e
v = v
fp = ideal_gas
[../]
[]
[Kernels]
[./diff]
type = Diffusion
variable = u
[../]
[./ffn]
type = BodyForce
variable = u
function = f_fn
[../]
[]
[BCs]
[./all]
type = FunctionDirichletBC
variable = u
boundary = 'left right top bottom'
function = bc_fn
[../]
[]
[Executioner]
type = Steady
solve_type = NEWTON
[]
[Outputs]
exodus = true
[]
(modules/fluid_properties/test/tests/fp_interrogator/2ph_ncg_p_T.i)
[FluidPropertiesInterrogator]
fp = fp_2phase_ncg
p = 1e5
T = 372.7559289
x_ncg = '0.1'
[]
[Modules]
[./FluidProperties]
[./fp_nitrogen]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 0.02867055103448276
[../]
[./fp_2phase_ncg]
type = TestTwoPhaseNCGFluidProperties
fp_ncgs = 'fp_nitrogen'
[../]
[../]
[]
(modules/fluid_properties/test/tests/fp_interrogator/2ph.p_T.i)
[FluidPropertiesInterrogator]
fp = fp
p = 1e5
T = 300
[]
[Modules]
[./FluidProperties]
[./fp_liquid]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 0.02900055737704918
mu = 1.823e-05
k = 0.02568
[../]
[./fp_vapor]
type = IdealGasFluidProperties
gamma = 1.1
molar_mass = 0.027714866
mu = 1.7e-05
k = 0.05
[../]
[./fp]
type = TestTwoPhaseFluidProperties
fp_liquid = fp_liquid
fp_vapor = fp_vapor
[../]
[../]
[]
(modules/porous_flow/test/tests/fluids/ideal_gas.i)
# Example of using the IdealGasFluidProperties userobject to provide fluid
# properties for an ideal gas. Use values for hydrogen (H2) at 1 MPa and 50 C.
#
# Input values:
# M = 2.01588e-3 kg/mol
# gamma = 1.4
# viscosity = 9.4393e-6 Pa.s
#
# Expected output:
# density = 750.2854 kg/m^3
# internal energy = 3.33 MJ/kg
# enthalpy = 4.66 MJ/kg
[Mesh]
type = GeneratedMesh
dim = 1
nx = 1
[]
[GlobalParams]
PorousFlowDictator = dictator
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pp'
number_fluid_phases = 1
number_fluid_components = 1
[]
[]
[Variables]
[pp]
initial_condition = 1e6
[]
[]
[Kernels]
[dummy]
type = Diffusion
variable = pp
[]
[]
[AuxVariables]
[temp]
initial_condition = 50.0
[]
[]
[Modules]
[FluidProperties]
[idealgas]
type = IdealGasFluidProperties
molar_mass = 2.01588e-3
gamma = 1.4
mu = 9.4393e-6
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = temp
[]
[ppss]
type = PorousFlow1PhaseFullySaturated
porepressure = pp
[]
[idealgass]
type = PorousFlowSingleComponentFluid
temperature_unit = Celsius
fp = idealgas
phase = 0
[]
[]
[Executioner]
type = Steady
solve_type = Newton
[]
[Postprocessors]
[pressure]
type = ElementIntegralVariablePostprocessor
variable = pp
[]
[temperature]
type = ElementIntegralVariablePostprocessor
variable = temp
[]
[density]
type = ElementIntegralMaterialProperty
mat_prop = 'PorousFlow_fluid_phase_density_qp0'
[]
[viscosity]
type = ElementIntegralMaterialProperty
mat_prop = 'PorousFlow_viscosity_qp0'
[]
[internal_energy]
type = ElementIntegralMaterialProperty
mat_prop = 'PorousFlow_fluid_phase_internal_energy_qp0'
[]
[enthalpy]
type = ElementIntegralMaterialProperty
mat_prop = 'PorousFlow_fluid_phase_enthalpy_qp0'
[]
[]
[Outputs]
execute_on = 'timestep_end'
file_base = ideal_gas
csv = true
[]
(modules/fluid_properties/test/tests/fp_interrogator/1ph.p_T.i)
[FluidPropertiesInterrogator]
fp = fp
p = 1e5
T = 300
[]
[Modules]
[./FluidProperties]
[./fp]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 0.02900055737704918
mu = 1.823e-05
k = 0.02568
[../]
[../]
[]
(modules/fluid_properties/test/tests/fp_interrogator/1ph.rho_e.i)
[FluidPropertiesInterrogator]
fp = fp
rho = 1
e = 2.1502500000e+05
[]
[Modules]
[./FluidProperties]
[./fp]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 0.02900055737704918
mu = 1.823e-05
k = 0.02568
[../]
[../]
[]
(modules/fluid_properties/test/tests/fp_interrogator/2ph.T.i)
[FluidPropertiesInterrogator]
fp = fp
T = 300
[]
[Modules]
[./FluidProperties]
[./fp_liquid]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 0.02900055737704918
mu = 1.823e-05
k = 0.02568
[../]
[./fp_vapor]
type = IdealGasFluidProperties
gamma = 1.1
molar_mass = 0.027714866
mu = 1.7e-05
k = 0.05
[../]
[./fp]
type = TestTwoPhaseFluidProperties
fp_liquid = fp_liquid
fp_vapor = fp_vapor
[../]
[../]
[]
(modules/fluid_properties/test/tests/fp_interrogator/1ph.rho_rhou_rhoE.i)
[FluidPropertiesInterrogator]
fp = fp
rho = 0.5
rhou = 0.5
rhoE = 2.75
[]
[Modules]
[./FluidProperties]
[./fp]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 11.640243719999999
[../]
[../]
[]
(modules/fluid_properties/test/tests/fp_interrogator/err.no_params.i)
[FluidPropertiesInterrogator]
fp = fp
[]
[Modules]
[./FluidProperties]
[./fp]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 0.02900055737704918
mu = 1.823e-05
k = 0.02568
[../]
[../]
[]
(modules/porous_flow/test/tests/newton_cooling/nc08.i)
# Newton cooling from a bar. 1-phase ideal fluid and heat, steady
[Mesh]
type = GeneratedMesh
dim = 2
nx = 100
ny = 1
xmin = 0
xmax = 100
ymin = 0
ymax = 1
[]
[GlobalParams]
PorousFlowDictator = dictator
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pressure temp'
number_fluid_phases = 1
number_fluid_components = 1
[]
[pc]
type = PorousFlowCapillaryPressureVG
m = 0.8
alpha = 1e-5
[]
[]
[Variables]
[pressure]
[]
[temp]
[]
[]
[ICs]
# have to start these reasonably close to their steady-state values
[pressure]
type = FunctionIC
variable = pressure
function = '200-0.5*x'
[]
[temperature]
type = FunctionIC
variable = temp
function = 180+0.1*x
[]
[]
[Kernels]
[flux]
type = PorousFlowAdvectiveFlux
fluid_component = 0
gravity = '0 0 0'
variable = pressure
[]
[heat_advection]
type = PorousFlowHeatAdvection
gravity = '0 0 0'
variable = temp
[]
[]
[Modules]
[FluidProperties]
[idealgas]
type = IdealGasFluidProperties
molar_mass = 1.4
gamma = 1.2
mu = 1.2
[]
[]
[]
[Materials]
[temperature]
type = PorousFlowTemperature
temperature = temp
[]
[ppss]
type = PorousFlow1PhaseP
porepressure = pressure
capillary_pressure = pc
[]
[massfrac]
type = PorousFlowMassFraction
[]
[dens0]
type = PorousFlowSingleComponentFluid
fp = idealgas
phase = 0
[]
[permeability]
type = PorousFlowPermeabilityConst
permeability = '1.1 0 0 0 1.1 0 0 0 1.1'
[]
[relperm]
type = PorousFlowRelativePermeabilityCorey # irrelevant in this fully-saturated situation
n = 2
phase = 0
[]
[]
[BCs]
[leftp]
type = DirichletBC
variable = pressure
boundary = left
value = 200
[]
[leftt]
type = DirichletBC
variable = temp
boundary = left
value = 180
[]
[newtonp]
type = PorousFlowPiecewiseLinearSink
variable = pressure
boundary = right
pt_vals = '-200 0 200'
multipliers = '-200 0 200'
use_mobility = true
use_relperm = true
fluid_phase = 0
flux_function = 0.005 # 1/2/L
[]
[newtont]
type = PorousFlowPiecewiseLinearSink
variable = temp
boundary = right
pt_vals = '-200 0 200'
multipliers = '-200 0 200'
use_mobility = true
use_relperm = true
use_enthalpy = true
fluid_phase = 0
flux_function = 0.005 # 1/2/L
[]
[]
[VectorPostprocessors]
[porepressure]
type = LineValueSampler
variable = pressure
start_point = '0 0.5 0'
end_point = '100 0.5 0'
sort_by = x
num_points = 11
execute_on = timestep_end
[]
[temperature]
type = LineValueSampler
variable = temp
start_point = '0 0.5 0'
end_point = '100 0.5 0'
sort_by = x
num_points = 11
execute_on = timestep_end
[]
[]
[Preconditioning]
[andy]
type = SMP
full = true
[]
[]
[Executioner]
type = Steady
solve_type = Newton
nl_rel_tol = 1E-10
nl_abs_tol = 1E-15
[]
[Outputs]
file_base = nc08
execute_on = timestep_end
exodus = true
[along_line]
type = CSV
execute_vector_postprocessors_on = timestep_end
[]
[]
(modules/fluid_properties/test/tests/ics/rho_vapor_mixture_from_pressure_temperature/test.i)
# Tests the initial condition for mixture density from pressure and temperature.
# This test uses the general vapor mixture fluid properties with steam, air,
# and helium with mass fractions 0.5, 0.3, and 0.2, respectively. The individual
# specific volumes (in m^3/kg) at p = 100 kPa, T = 500 K are:
# steam: 2.298113001
# air: 1.43525
# helium: 10.3855
# For the general vapor mixture, the mixture specific volume is computed as
# v = \sum\limits_i x_i v_i ,
# where x_i is the mass fraction of component i, and v_i is the specific volume
# of component i. Therefore, the correct value for specific volume of the mixture is
# v = 3.65673150050 m^3/kg
# and thus density is
# rho = 0.27346825980066236 kg/m^3
[Mesh]
type = GeneratedMesh
dim = 1
nx = 1
allow_renumbering = false
[]
[Modules]
[FluidProperties]
[fp_steam]
type = StiffenedGasFluidProperties
gamma = 1.43
cv = 1040.0
q = 2.03e6
p_inf = 0.0
q_prime = -2.3e4
k = 0.026
mu = 134.4e-7
M = 0.01801488
rho_c = 322.0
[]
[fp_air]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 28.965197004e-3
[]
[fp_helium]
type = IdealGasFluidProperties
gamma = 1.66
molar_mass = 4.002917432959e-3
[]
[fp_vapor_mixture]
type = IdealRealGasMixtureFluidProperties
fp_primary = fp_steam
fp_secondary = 'fp_air fp_helium'
[]
[]
[]
[AuxVariables]
[rho]
[]
[p]
[]
[T]
[]
[x_air]
[]
[x_helium]
[]
[]
[ICs]
[rho_ic]
type = RhoVaporMixtureFromPressureTemperatureIC
variable = rho
p = p
T = T
x_secondary_vapors = 'x_air x_helium'
fp_vapor_mixture = fp_vapor_mixture
[]
[p_ic]
type = ConstantIC
variable = p
value = 100e3
[]
[T_ic]
type = ConstantIC
variable = T
value = 500
[]
[x_air_ic]
type = ConstantIC
variable = x_air
value = 0.3
[]
[x_helium_ic]
type = ConstantIC
variable = x_helium
value = 0.2
[]
[]
[Executioner]
type = Steady
[]
[Postprocessors]
[rho_test]
type = ElementalVariableValue
elementid = 0
variable = rho
execute_on = 'INITIAL TIMESTEP_END'
[]
[]
[Outputs]
csv = true
execute_on = 'INITIAL'
[]
[Problem]
solve = false
[]
(modules/fluid_properties/test/tests/fp_interrogator/2ph.p.i)
[FluidPropertiesInterrogator]
fp = fp
p = 1e5
[]
[Modules]
[./FluidProperties]
[./fp_liquid]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 0.02900055737704918
mu = 1.823e-05
k = 0.02568
[../]
[./fp_vapor]
type = IdealGasFluidProperties
gamma = 1.1
molar_mass = 0.027714866
mu = 1.7e-05
k = 0.05
[../]
[./fp]
type = TestTwoPhaseFluidProperties
fp_liquid = fp_liquid
fp_vapor = fp_vapor
[../]
[../]
[]
(modules/porous_flow/test/tests/jacobian/esbc02.i)
# Tests the Jacobian of PorousFlowEnthalpySink when pressure
[Mesh]
type = GeneratedMesh
dim = 2
[]
[GlobalParams]
PorousFlowDictator = dictator
at_nodes = true
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pp temp'
number_fluid_phases = 1
number_fluid_components = 1
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0.1
[]
[]
[Variables]
[pp]
initial_condition = 1
[]
[temp]
initial_condition = 2
[]
[]
[AuxVariables]
[pressure]
[]
[]
[Kernels]
[mass0]
type = TimeDerivative
variable = pp
[]
[heat_conduction]
type = TimeDerivative
variable = temp
[]
[]
[Modules]
[FluidProperties]
[simple_fluid]
type = IdealGasFluidProperties
[]
[]
[]
[Materials]
[ppss]
type = PorousFlow1PhaseFullySaturated
porepressure = pp
[]
[simple_fluid]
type = PorousFlowSingleComponentFluid
fp = simple_fluid
phase = 0
[]
[temperature]
type = PorousFlowTemperature
temperature = temp
[]
[]
[BCs]
[left]
type = PorousFlowEnthalpySink
variable = temp
boundary = left
porepressure_var = pressure
T_in = 300
fp = simple_fluid
flux_function = -23
[]
[]
[Preconditioning]
[andy]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = Newton
dt = 0.1
num_steps = 1
nl_rel_tol = 1E-12
nl_abs_tol = 1E-12
petsc_options_iname = '-snes_test_err'
petsc_options_value = '1e-1'
[]
(modules/fluid_properties/test/tests/functions/saturation_temperature_function/saturation_temperature_function.i)
# TestTwoPhaseFluidProperties has the following saturation temperature function:
# T_sat(p) = 2 p
# Thus for p = 5, T_sat should be 10.
[Mesh]
type = GeneratedMesh
dim = 1
nx = 1
[]
[Modules]
[./FluidProperties]
[./fp_liquid]
type = IdealGasFluidProperties
[../]
[./fp_vapor]
type = IdealGasFluidProperties
[../]
[./fp_2phase]
type = TestTwoPhaseFluidProperties
fp_liquid = fp_liquid
fp_vapor = fp_vapor
[../]
[]
[]
[Functions]
[./p]
type = ConstantFunction
value = 5
[../]
[./T_sat]
type = SaturationTemperatureFunction
p = p
fp_2phase = fp_2phase
[../]
[]
[Postprocessors]
[./T_sat_pp]
type = FunctionValuePostprocessor
function = T_sat
execute_on = 'INITIAL'
[../]
[]
[Problem]
solve = false
[]
[Executioner]
type = Steady
[]
[Outputs]
csv = true
[]
(modules/fluid_properties/test/tests/two_phase_fluid_properties_independent/test.i)
# Tests the TwoPhaseFluidPropertiesIndependent class, which takes the names
# of 2 single-phase fluid properties independently. This test uses a dummy
# aux to make sure that the single-phase fluid properties can be recovered
# from the 2-phase fluid properties. A modification to this test checks that
# an error results if one tries to call a 2-phase fluid properties interface
# using this class, which is designed to ensure that the 2 phases are independent.
[Mesh]
type = GeneratedMesh
dim = 1
nx = 1
# Required for NodalVariableValue on distributed mesh
allow_renumbering = false
[]
[Problem]
solve = false
[]
[AuxVariables]
[./p]
initial_condition = 1e5
[../]
[./T]
initial_condition = 300
[../]
[./rho_avg]
[../]
[]
[Modules]
[./FluidProperties]
# rho1 = 1.149425287 kg/m^3
[./fp1]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 0.02867055103448276
[../]
# rho2 = 0.6666666667 kg/m^3
[./fp2]
type = IdealGasFluidProperties
gamma = 1.2
molar_mass = 0.0166289196
[../]
[./fp_2phase]
type = TwoPhaseFluidPropertiesIndependent
fp_liquid = fp1
fp_vapor = fp2
[../]
[]
[]
[AuxKernels]
# correct value (0.5*(rho1 + rho2)) should be: 0.90804597685 kg/m^3
[./rho_avg_aux]
type = TwoPhaseAverageDensityAux
variable = rho_avg
p = p
T = T
fp_2phase = fp_2phase
execute_on = 'initial'
[../]
[]
[Postprocessors]
[./rho_avg_value]
type = NodalVariableValue
variable = rho_avg
nodeid = 0
[../]
[]
[Executioner]
type = Steady
[]
[Outputs]
execute_on = 'timestep_end'
csv = true
[]
(modules/porous_flow/test/tests/jacobian/esbc01.i)
# Tests the Jacobian of PorousFlowEnthalpySink when pore pressure is specified
[Mesh]
type = GeneratedMesh
dim = 2
[]
[GlobalParams]
PorousFlowDictator = dictator
at_nodes = true
[]
[UserObjects]
[dictator]
type = PorousFlowDictator
porous_flow_vars = 'pp temp'
number_fluid_phases = 1
number_fluid_components = 1
[]
[pc]
type = PorousFlowCapillaryPressureConst
pc = 0.1
[]
[]
[Variables]
[pp]
initial_condition = 1
[]
[temp]
initial_condition = 2
[]
[]
[Kernels]
[mass0]
type = TimeDerivative
variable = pp
[]
[heat_conduction]
type = TimeDerivative
variable = temp
[]
[]
[Modules]
[FluidProperties]
[simple_fluid]
type = IdealGasFluidProperties
[]
[]
[]
[Materials]
[ppss]
type = PorousFlow1PhaseFullySaturated
porepressure = pp
[]
[simple_fluid]
type = PorousFlowSingleComponentFluid
fp = simple_fluid
phase = 0
[]
[temperature]
type = PorousFlowTemperature
temperature = temp
[]
[]
[BCs]
[left]
type = PorousFlowEnthalpySink
variable = temp
boundary = left
fluid_phase = 0
T_in = 300
fp = simple_fluid
flux_function = -23
[]
[]
[Preconditioning]
[andy]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = Newton
dt = 0.1
num_steps = 1
nl_rel_tol = 1E-12
nl_abs_tol = 1E-12
petsc_options_iname = '-snes_test_err'
petsc_options_value = '1e-2'
[]
(modules/fluid_properties/test/tests/fp_interrogator/vapor_mixture_rho_e.i)
[FluidPropertiesInterrogator]
fp = fp_vapor_mix
rho = 1.1870052372064208
e = 2477165.9033225174
x_ncg = '0.1'
[]
[Modules]
[./FluidProperties]
[./fp_nitrogen]
type = IdealGasFluidProperties
gamma = 1.4
molar_mass = 0.02867055103448276
[../]
[./fp_primary]
type = IdealGasFluidProperties
gamma = 1.3
molar_mass = 0.027714866
T_c = 126.19
rho_c = 313.189812
[../]
[./fp_vapor_mix]
type = IdealRealGasMixtureFluidProperties
fp_primary = fp_primary
fp_secondary = 'fp_nitrogen'
[../]
[../]
[]
(modules/navier_stokes/test/tests/step/step.i)
# Navier-Stokes (or Euler) flow of an ideal gas over a step.
#
# Note: this problem is not currently a regression test for the
# Navier-Stokes module since it is in some sense ill-posed. As
# discussed in [0], the sharp corner of the step (both forward and
# backward-facing) introduces a singularity in the first derivative of
# the velocity and pressure fields, and therefore produces large
# numerical errors in the neighborhood of these points. Physically,
# this numerical error can be interpreted as causing an artificial
# "boundary layer" to form just above the step, as well as a spurious
# production of entropy even though the flow remains subsonic.
# Nevertheless, the forward-facing step problem in particular remains
# a challenging and well-document test problem for flow solvers, and
# this input file is included to help facilitate its development and
# employment by users of the module.
#
# [0]: Woodward and Colella, "The numerical simulation of
# two-dimenstional fluid flow with strong shocks," Journal of
# Computational Physics 54(1), pp. 115-173, 1984
[Mesh]
type = FileMesh
file = step.e
dim = 2
# uniform_refine = 3
[]
[Modules]
[FluidProperties]
[ideal_gas]
type = IdealGasFluidProperties
gamma = 1.4
[]
[]
[CompressibleNavierStokes]
# steady-state or transient
equation_type = transient
# fluid
fluid_properties = ideal_gas
# boundary conditions
stagnation_boundary = left
stagnation_pressure = 120192.995549849 # Pa, Mach=0.5 at 1 atm
stagnation_temperature = 315 # K, Mach=0.5 at 1 atm
stagnation_flow_direction = '1 0'
no_penetration_boundary = 'top bottom step_top step_left step_right'
static_pressure_boundary = 'right'
static_pressure = 101325 # Pa
# variable types, scalings and initial conditions
family = LAGRANGE
order = FIRST
total_energy_scaling = 9.869232667160121e-6
initial_pressure = 101325.
initial_temperature = 300.
initial_velocity = '173.594354746921 0 0' # Mach 0.5: = 0.5*sqrt(gamma*R*T)
[]
[]
[Materials]
[fluid]
type = Air
block = 1
rho = rho
rhou = rhou
rhov = rhov
rhoE = rhoE
vel_x = vel_x
vel_y = vel_y
temperature = temperature
enthalpy = enthalpy
# This value is not used in the Euler equations, but it *is* used
# by the stabilization parameter computation, which it decreases
# the amount of artificial viscosity added, so it's best to use a
# realistic value.
dynamic_viscosity = 0.0
fluid_properties = ideal_gas
[]
[]
[Preconditioning]
[SMP_PJFNK]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
dt = 5.e-5
dtmin = 1.e-5
start_time = 0.0
num_steps = 10000
nl_rel_tol = 1e-5
nl_abs_tol = 1e-9
# nl_abs_step_tol = 1e-15
nl_max_its = 5
l_tol = 1e-4 # Relative linear tolerance for each Krylov solve
l_max_its = 100 # Number of linear iterations for each Krylov solve
# Specify the order as FIRST, otherwise you will get warnings in DEBUG mode...
[Quadrature]
type = TRAP
order = FIRST
[]
[]
[Outputs]
file_base = step_out
interval = 1
exodus = true
[]
(modules/fluid_properties/test/tests/functions/saturation_pressure_function/saturation_pressure_function.i)
# TestTwoPhaseFluidProperties has the following saturation pressure function:
# p_sat(p) = 3 T
# Thus for T = 5, p_sat should be 15.
[Mesh]
type = GeneratedMesh
dim = 1
nx = 1
[]
[Modules]
[./FluidProperties]
[./fp_liquid]
type = IdealGasFluidProperties
[../]
[./fp_vapor]
type = IdealGasFluidProperties
[../]
[./fp_2phase]
type = TestTwoPhaseFluidProperties
fp_liquid = fp_liquid
fp_vapor = fp_vapor
[../]
[]
[]
[Functions]
[./T]
type = ConstantFunction
value = 5
[../]
[./p_sat]
type = SaturationPressureFunction
T = T
fp_2phase = fp_2phase
[../]
[]
[Postprocessors]
[./p_sat_pp]
type = FunctionValuePostprocessor
function = p_sat
execute_on = 'INITIAL'
[../]
[]
[Problem]
solve = false
[]
[Executioner]
type = Steady
[]
[Outputs]
csv = true
[]