Line data    Source code

       1             : //* This file is part of the MOOSE framework
       2             : //* https://mooseframework.inl.gov
       3             : //*
       4             : //* All rights reserved, see COPYRIGHT for full restrictions
       5             : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
       6             : //*
       7             : //* Licensed under LGPL 2.1, please see LICENSE for details
       8             : //* https://www.gnu.org/licenses/lgpl-2.1.html
       9             : 
      10             : #include "PorousFlowDarcyVelocityMaterial.h"
      11             : 
      12             : registerMooseObject("PorousFlowApp", PorousFlowDarcyVelocityMaterial);
      13             : 
      14             : InputParameters
      15         253 : PorousFlowDarcyVelocityMaterial::validParams()
      16             : {
      17         253 :   InputParameters params = PorousFlowMaterial::validParams();
      18         506 :   params.addRequiredParam<RealVectorValue>("gravity",
      19             :                                            "Gravitational acceleration vector downwards (m/s^2)");
      20         253 :   params.addClassDescription("This Material calculates the Darcy velocity for all phases");
      21         506 :   params.addPrivateParam<std::string>("pf_material_type", "darcy_velocity");
      22         253 :   params.set<bool>("at_nodes") = false;
      23         253 :   return params;
      24           0 : }
      25             : 
      26         197 : PorousFlowDarcyVelocityMaterial::PorousFlowDarcyVelocityMaterial(const InputParameters & parameters)
      27             :   : PorousFlowMaterial(parameters),
      28         394 :     _num_phases(_dictator.numPhases()),
      29         197 :     _num_var(_dictator.numVariables()),
      30         394 :     _permeability(getMaterialProperty<RealTensorValue>("PorousFlow_permeability_qp")),
      31         197 :     _dpermeability_dvar(
      32         197 :         getMaterialProperty<std::vector<RealTensorValue>>("dPorousFlow_permeability_qp_dvar")),
      33         394 :     _dpermeability_dgradvar(getMaterialProperty<std::vector<std::vector<RealTensorValue>>>(
      34             :         "dPorousFlow_permeability_qp_dgradvar")),
      35         394 :     _fluid_density(getMaterialProperty<std::vector<Real>>("PorousFlow_fluid_phase_density_qp")),
      36         394 :     _dfluid_density_dvar(getMaterialProperty<std::vector<std::vector<Real>>>(
      37             :         "dPorousFlow_fluid_phase_density_qp_dvar")),
      38         394 :     _fluid_viscosity(getMaterialProperty<std::vector<Real>>("PorousFlow_viscosity_qp")),
      39         197 :     _dfluid_viscosity_dvar(
      40         197 :         getMaterialProperty<std::vector<std::vector<Real>>>("dPorousFlow_viscosity_qp_dvar")),
      41         197 :     _relative_permeability(
      42         197 :         getMaterialProperty<std::vector<Real>>("PorousFlow_relative_permeability_qp")),
      43         394 :     _drelative_permeability_dvar(getMaterialProperty<std::vector<std::vector<Real>>>(
      44             :         "dPorousFlow_relative_permeability_qp_dvar")),
      45         394 :     _grad_p(getMaterialProperty<std::vector<RealGradient>>("PorousFlow_grad_porepressure_qp")),
      46         394 :     _dgrad_p_dgradvar(getMaterialProperty<std::vector<std::vector<Real>>>(
      47             :         "dPorousFlow_grad_porepressure_qp_dgradvar")),
      48         394 :     _dgrad_p_dvar(getMaterialProperty<std::vector<std::vector<RealGradient>>>(
      49             :         "dPorousFlow_grad_porepressure_qp_dvar")),
      50         394 :     _gravity(getParam<RealVectorValue>("gravity")),
      51         197 :     _darcy_velocity(declareProperty<std::vector<RealVectorValue>>("PorousFlow_darcy_velocity_qp")),
      52         197 :     _ddarcy_velocity_dvar(declareProperty<std::vector<std::vector<RealVectorValue>>>(
      53             :         "dPorousFlow_darcy_velocity_qp_dvar")),
      54         197 :     _ddarcy_velocity_dgradvar(
      55         197 :         declareProperty<std::vector<std::vector<std::vector<RealVectorValue>>>>(
      56         197 :             "dPorousFlow_darcy_velocity_qp_dgradvar"))
      57             : {
      58         197 :   if (_nodal_material == true)
      59           2 :     mooseError("PorousFlowDarcyVelocityMaterial is only defined for at_nodes = false");
      60         195 : }
      61             : 
      62             : void
      63        4168 : PorousFlowDarcyVelocityMaterial::computeQpProperties()
      64             : {
      65        4168 :   _darcy_velocity[_qp].resize(_num_phases);
      66             : 
      67       10284 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
      68        6116 :     _darcy_velocity[_qp][ph] =
      69        6116 :         -(_permeability[_qp] * (_grad_p[_qp][ph] - _fluid_density[_qp][ph] * _gravity) *
      70        6116 :           _relative_permeability[_qp][ph] / _fluid_viscosity[_qp][ph]);
      71             : 
      72        4168 :   _ddarcy_velocity_dvar[_qp].resize(_num_phases);
      73       10284 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
      74        6116 :     _ddarcy_velocity_dvar[_qp][ph].resize(_num_var);
      75             : 
      76       10284 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
      77        6860 :     for (unsigned v = 0; v < _num_var; ++v)
      78             :     {
      79         744 :       _ddarcy_velocity_dvar[_qp][ph][v] =
      80         744 :           -_dpermeability_dvar[_qp][v] * (_grad_p[_qp][ph] - _fluid_density[_qp][ph] * _gravity) *
      81         744 :           _relative_permeability[_qp][ph] / _fluid_viscosity[_qp][ph];
      82             :       _ddarcy_velocity_dvar[_qp][ph][v] -=
      83        1488 :           _permeability[_qp] *
      84         744 :           (_dgrad_p_dvar[_qp][ph][v] - _dfluid_density_dvar[_qp][ph][v] * _gravity) *
      85         744 :           _relative_permeability[_qp][ph] / _fluid_viscosity[_qp][ph];
      86             :       _ddarcy_velocity_dvar[_qp][ph][v] -=
      87        1488 :           _permeability[_qp] * (_grad_p[_qp][ph] - _fluid_density[_qp][ph] * _gravity) *
      88         744 :           (_drelative_permeability_dvar[_qp][ph][v] / _fluid_viscosity[_qp][ph] -
      89         744 :            _relative_permeability[_qp][ph] * _dfluid_viscosity_dvar[_qp][ph][v] /
      90             :                std::pow(_fluid_viscosity[_qp][ph], 2));
      91             :     }
      92             : 
      93        4168 :   _ddarcy_velocity_dgradvar[_qp].resize(_num_phases);
      94       10284 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
      95             :   {
      96        6116 :     _ddarcy_velocity_dgradvar[_qp][ph].resize(LIBMESH_DIM);
      97       24464 :     for (unsigned i = 0; i < LIBMESH_DIM; ++i)
      98       18348 :       _ddarcy_velocity_dgradvar[_qp][ph][i].resize(_num_var);
      99             :   }
     100             : 
     101       10284 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
     102       24464 :     for (unsigned i = 0; i < LIBMESH_DIM; ++i)
     103       20580 :       for (unsigned v = 0; v < _num_var; ++v)
     104             :       {
     105        2232 :         _ddarcy_velocity_dgradvar[_qp][ph][i][v] =
     106        2232 :             -_dpermeability_dgradvar[_qp][i][v] *
     107        2232 :             (_grad_p[_qp][ph] - _fluid_density[_qp][ph] * _gravity);
     108        8928 :         for (unsigned j = 0; j < LIBMESH_DIM; ++j)
     109        6696 :           _ddarcy_velocity_dgradvar[_qp][ph][i][v](j) -=
     110        6696 :               _permeability[_qp](j, i) * _dgrad_p_dgradvar[_qp][ph][v];
     111             :         _ddarcy_velocity_dgradvar[_qp][ph][i][v] *=
     112        2232 :             _relative_permeability[_qp][ph] / _fluid_viscosity[_qp][ph];
     113             :       }
     114        4168 : }