20 "richardsVarNames_UO",
"The UserObject that holds the list of Richards variable names.");
21 params.
addParam<std::vector<UserObjectName>>(
"relperm_UO",
22 "List of names of user objects that " 23 "define relative permeability. Only " 24 "needed if fully_upwind is used");
25 params.
addParam<std::vector<UserObjectName>>(
27 "List of name of user objects that define effective saturation as a function of " 28 "pressure list. Only needed if fully_upwind is used");
29 params.
addParam<std::vector<UserObjectName>>(
"density_UO",
30 "List of names of user objects that " 31 "define the fluid density. Only " 32 "needed if fully_upwind is used");
33 params.
addParam<
bool>(
"fully_upwind",
false,
"Fully upwind the flux");
34 params.
addClassDescription(
"Approximates a borehole in the mesh with given bottomhole pressure, " 35 "and radii using a number of point sinks whose positions are read " 42 _fully_upwind(getParam<bool>(
"fully_upwind")),
44 _num_p(_richards_name_UO.num_v()),
45 _pvar(_richards_name_UO.richards_var_num(_var.number())),
50 getParam<
std::vector<UserObjectName>>(
"density_UO")[_pvar])
52 _seff_UO(_fully_upwind ? &getUserObjectByName<
RichardsSeff>(
53 getParam<
std::vector<UserObjectName>>(
"seff_UO")[_pvar])
56 getParam<
std::vector<UserObjectName>>(
"relperm_UO")[_pvar])
62 _pp(getMaterialProperty<
std::vector<
Real>>(
"porepressure")),
63 _dpp_dv(getMaterialProperty<
std::vector<
std::vector<
Real>>>(
"dporepressure_dv")),
64 _viscosity(getMaterialProperty<
std::vector<
Real>>(
"viscosity")),
66 _dseff_dv(getMaterialProperty<
std::vector<
std::vector<
Real>>>(
"ds_eff_dv")),
67 _rel_perm(getMaterialProperty<
std::vector<
Real>>(
"rel_perm")),
68 _drel_perm_dv(getMaterialProperty<
std::vector<
std::vector<
Real>>>(
"drel_perm_dv")),
69 _density(getMaterialProperty<
std::vector<
Real>>(
"density")),
70 _ddensity_dv(getMaterialProperty<
std::vector<
std::vector<
Real>>>(
"ddensity_dv"))
73 for (
unsigned int pnum = 0; pnum <
_num_p; ++pnum)
79 for (
unsigned int i = 0; i < coupled_vars.size(); i++)
94 std::vector<Real> dseff_dp;
100 for (
unsigned int nodenum = 0; nodenum <
_num_nodes; ++nodenum)
117 for (
unsigned int ph = 0; ph <
_num_p; ++ph)
135 if (character == 0.0)
138 const Real bh_pressure =
167 if (current_dirac_ptid > 0)
175 _rs[current_dirac_ptid]);
176 if ((character < 0.0 && pp < bh_pressure) || (character > 0.0 && pp > bh_pressure))
181 if (current_dirac_ptid + 1 <
_zs.size() ||
_zs.size() == 1)
188 _rs[current_dirac_ptid]);
189 if ((character < 0.0 && pp < bh_pressure) || (character > 0.0 && pp > bh_pressure))
211 if (character == 0.0)
229 if (character == 0.0)
232 const Real bh_pressure =
259 (
_pvar == wrt_num ? 1 : 0);
275 if (current_dirac_ptid > 0)
282 _rs[current_dirac_ptid]);
283 if ((character < 0.0 && pp < bh_pressure) || (character > 0.0 && pp > bh_pressure))
286 (mob * phi * dpp_dv + dmob_dv * phi * (pp - bh_pressure));
289 if (current_dirac_ptid <
_zs.size() - 1 ||
_zs.size() == 1)
296 _rs[current_dirac_ptid]);
297 if ((character < 0.0 && pp < bh_pressure) || (character > 0.0 && pp > bh_pressure))
300 (mob * phi * dpp_dv + dmob_dv * phi * (pp - bh_pressure));
virtual void dseff(std::vector< const VariableValue *> p, unsigned int qp, std::vector< Real > &result) const =0
derivative(s) of effective saturation as a function of porepressure(s) at given quadpoint of the elem...
const MaterialProperty< RealTensorValue > & _permeability
material permeability
const bool _fully_upwind
Whether to use full upwinding.
const MaterialProperty< std::vector< Real > > & _density
fluid density
Real jac(unsigned int wrt_num)
Calculates Jacobian.
const MooseArray< Point > & _q_point
virtual Real drelperm(Real seff) const =0
derivative of relative permeability wrt effective saturation This must be over-ridden in your derived...
virtual Real ddensity(Real p) const =0
derivative of fluid density wrt porepressure This must be over-ridden in derived classes to provide a...
Real wellConstant(const RealTensorValue &perm, const RealTensorValue &rot, const Real &half_len, const Elem *ele, const Real &rad)
Calculates Peaceman's form of the borehole well constant Z Chen, Y Zhang, Well flow models for variou...
const MaterialProperty< std::vector< Real > > & _viscosity
fluid viscosity
Base class for effective saturation as a function of porepressure(s) The functions seff...
Base class for Richards relative permeability classes that provide relative permeability as a functio...
const VariableValue * nodal_var(unsigned int richards_var_num) const
The nodal variable values for the given richards_var_num To extract a the value of pressure variable ...
std::vector< Real > _rs
radii of the borehole
const OutputTools< T >::VariableTestValue & _test
const OutputTools< T >::VariablePhiValue & _phi
virtual Real computeQpJacobian()
Computes the diagonal part of the jacobian.
static const std::string density
unsigned int _num_nodes
number of nodes in this element. Only used if _fully_upwind = true
bool not_richards_var(unsigned int moose_var_num) const
returns true if moose_var_num is not a richards var
unsigned currentPointCachedID()
Approximates a borehole by a sequence of Dirac Points.
This holds maps between pressure_var or pressure_var, sat_var used in RichardsMaterial and kernels...
std::vector< const VariableValue * > _ps_at_nodes
Holds the values of pressures at all the nodes of the element Only used if _fully_upwind = true Eg: _...
static InputParameters validParams()
Creates a new PeacemanBorehole This reads the file containing the lines of the form radius x y z that...
std::vector< Real > _mobility
nodal values of mobility = density*relperm/viscosity These are used if _fully_upwind = true ...
const MaterialProperty< std::vector< Real > > & _pp
fluid porepressure (or porepressures in case of multiphase)
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
virtual Real density(Real p) const =0
fluid density as a function of porepressure This must be over-ridden in derived classes to provide an...
const RichardsRelPerm *const _relperm_UO
user object defining the relative permeability. Only used if _fully_upwind = true ...
TensorValue< Real > RealTensorValue
virtual void computeResidual()
Computes the residual.
const RichardsDensity *const _density_UO
user object defining the density. Only used if _fully_upwind = true
static InputParameters validParams()
Creates a new RichardsBorehole This sets all the variables, but also reads the file containing the li...
const MaterialProperty< std::vector< Real > > & _rel_perm
relative permeability
virtual Real seff(std::vector< const VariableValue *> p, unsigned int qp) const =0
effective saturation as a function of porepressure(s) at given quadpoint of the element ...
const RealVectorValue _unit_weight
unit weight of fluid in borehole (for calculating bottomhole pressure at each Dirac Point) ...
const RichardsSeff *const _seff_UO
user object defining the effective saturation. Only used if _fully_upwind = true
virtual Real relperm(Real seff) const =0
relative permeability as a function of effective saturation This must be over-ridden in your derived ...
virtual Real computeQpOffDiagJacobian(unsigned int jvar)
Computes the off-diagonal part of the jacobian Note: at March2014 this is never called since moose do...
void add(Real contrib)
adds contrib to _total
RichardsSumQuantity & _total_outflow_mass
This is used to hold the total fluid flowing into the borehole Hence, it is positive for production w...
const Elem *const & _current_elem
const std::vector< MooseVariableFieldBase *> & getCoupledMooseVars() const
const MaterialProperty< std::vector< std::vector< Real > > > & _drel_perm_dv
d(relperm_i)/d(variable_j)
const Real _p_bot
bottomhole pressure of borehole
virtual void computeJacobian()
Computes the Jacobian.
void addMooseVariableDependency(MooseVariableFieldBase *var)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
RichardsBorehole(const InputParameters ¶meters)
unsigned int richards_var_num(unsigned int moose_var_num) const
the richards variable number
const unsigned int _pvar
The moose internal variable number of the richards variable of this Dirac Kernel. ...
virtual void computeJacobian() override
std::vector< std::vector< Real > > _dmobility_dv
d(_mobility)/d(variable_ph) (variable_ph is the variable for phase=ph) These are used in the jacobian...
const MaterialProperty< std::vector< std::vector< Real > > > & _dpp_dv
d(porepressure_i)/d(variable_j)
Point _bottom_point
the bottom point of the borehole (where bottom_pressure is defined)
Base class for fluid density as a function of porepressure The functions density, ddensity and d2dens...
std::vector< RealTensorValue > _rot_matrix
rotation matrix used in well_constant calculation
const RichardsVarNames & _richards_name_UO
Defines the richards variables in the simulation.
std::vector< Real > _half_seg_len
0.5*(length of polyline segments between points)
virtual Real value(Real t, const Point &p) const
std::vector< Real > _zs
z points of borehole
virtual Real computeQpResidual()
Computes the Qp residual.
void prepareNodalValues()
calculates the nodal values of pressure, mobility, and derivatives thereof
registerMooseObject("RichardsApp", RichardsBorehole)
const unsigned int _num_p
number of richards variables
const MaterialProperty< std::vector< std::vector< Real > > > & _ddensity_dv
d(density_i)/d(variable_j)
virtual void computeResidual() override
const Function & _character
If positive then the borehole acts as a sink (producion well) for porepressure > borehole pressure...
Approximates a borehole by a sequence of Dirac Points.