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         429 : PorousFlowDarcyVelocityMaterial::validParams()
      16             : {
      17         429 :   InputParameters params = PorousFlowMaterial::validParams();
      18         858 :   params.addRequiredParam<RealVectorValue>("gravity",
      19             :                                            "Gravitational acceleration vector downwards (m/s^2)");
      20         429 :   params.addClassDescription("This Material calculates the Darcy velocity for all phases");
      21         858 :   params.addPrivateParam<std::string>("pf_material_type", "darcy_velocity");
      22         429 :   params.set<bool>("at_nodes") = false;
      23         429 :   return params;
      24           0 : }
      25             : 
      26         341 : PorousFlowDarcyVelocityMaterial::PorousFlowDarcyVelocityMaterial(const InputParameters & parameters)
      27             :   : PorousFlowMaterial(parameters),
      28         682 :     _num_phases(_dictator.numPhases()),
      29         341 :     _num_var(_dictator.numVariables()),
      30         682 :     _permeability(getMaterialProperty<RealTensorValue>("PorousFlow_permeability_qp")),
      31         341 :     _dpermeability_dvar(
      32         341 :         getMaterialProperty<std::vector<RealTensorValue>>("dPorousFlow_permeability_qp_dvar")),
      33         682 :     _dpermeability_dgradvar(getMaterialProperty<std::vector<std::vector<RealTensorValue>>>(
      34             :         "dPorousFlow_permeability_qp_dgradvar")),
      35         682 :     _fluid_density(getMaterialProperty<std::vector<Real>>("PorousFlow_fluid_phase_density_qp")),
      36         682 :     _dfluid_density_dvar(getMaterialProperty<std::vector<std::vector<Real>>>(
      37             :         "dPorousFlow_fluid_phase_density_qp_dvar")),
      38         682 :     _fluid_viscosity(getMaterialProperty<std::vector<Real>>("PorousFlow_viscosity_qp")),
      39         341 :     _dfluid_viscosity_dvar(
      40         341 :         getMaterialProperty<std::vector<std::vector<Real>>>("dPorousFlow_viscosity_qp_dvar")),
      41         341 :     _relative_permeability(
      42         341 :         getMaterialProperty<std::vector<Real>>("PorousFlow_relative_permeability_qp")),
      43         682 :     _drelative_permeability_dvar(getMaterialProperty<std::vector<std::vector<Real>>>(
      44             :         "dPorousFlow_relative_permeability_qp_dvar")),
      45         682 :     _grad_p(getMaterialProperty<std::vector<RealGradient>>("PorousFlow_grad_porepressure_qp")),
      46         682 :     _dgrad_p_dgradvar(getMaterialProperty<std::vector<std::vector<Real>>>(
      47             :         "dPorousFlow_grad_porepressure_qp_dgradvar")),
      48         682 :     _dgrad_p_dvar(getMaterialProperty<std::vector<std::vector<RealGradient>>>(
      49             :         "dPorousFlow_grad_porepressure_qp_dvar")),
      50         682 :     _gravity(getParam<RealVectorValue>("gravity")),
      51         341 :     _darcy_velocity(declareProperty<std::vector<RealVectorValue>>("PorousFlow_darcy_velocity_qp")),
      52         341 :     _ddarcy_velocity_dvar(declareProperty<std::vector<std::vector<RealVectorValue>>>(
      53             :         "dPorousFlow_darcy_velocity_qp_dvar")),
      54         341 :     _ddarcy_velocity_dgradvar(
      55         341 :         declareProperty<std::vector<std::vector<std::vector<RealVectorValue>>>>(
      56         341 :             "dPorousFlow_darcy_velocity_qp_dgradvar"))
      57             : {
      58         341 :   if (_nodal_material == true)
      59           2 :     mooseError("PorousFlowDarcyVelocityMaterial is only defined for at_nodes = false");
      60         339 : }
      61             : 
      62             : void
      63        6080 : PorousFlowDarcyVelocityMaterial::computeQpProperties()
      64             : {
      65        6080 :   _darcy_velocity[_qp].resize(_num_phases);
      66             : 
      67       15030 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
      68        8950 :     _darcy_velocity[_qp][ph] =
      69        8950 :         -(_permeability[_qp] * (_grad_p[_qp][ph] - _fluid_density[_qp][ph] * _gravity) *
      70        8950 :           _relative_permeability[_qp][ph] / _fluid_viscosity[_qp][ph]);
      71             : 
      72        6080 :   _ddarcy_velocity_dvar[_qp].resize(_num_phases);
      73       15030 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
      74        8950 :     _ddarcy_velocity_dvar[_qp][ph].resize(_num_var);
      75             : 
      76       15030 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
      77        9880 :     for (unsigned v = 0; v < _num_var; ++v)
      78             :     {
      79         930 :       _ddarcy_velocity_dvar[_qp][ph][v] =
      80         930 :           -_dpermeability_dvar[_qp][v] * (_grad_p[_qp][ph] - _fluid_density[_qp][ph] * _gravity) *
      81         930 :           _relative_permeability[_qp][ph] / _fluid_viscosity[_qp][ph];
      82             :       _ddarcy_velocity_dvar[_qp][ph][v] -=
      83        1860 :           _permeability[_qp] *
      84         930 :           (_dgrad_p_dvar[_qp][ph][v] - _dfluid_density_dvar[_qp][ph][v] * _gravity) *
      85         930 :           _relative_permeability[_qp][ph] / _fluid_viscosity[_qp][ph];
      86             :       _ddarcy_velocity_dvar[_qp][ph][v] -=
      87        1860 :           _permeability[_qp] * (_grad_p[_qp][ph] - _fluid_density[_qp][ph] * _gravity) *
      88         930 :           (_drelative_permeability_dvar[_qp][ph][v] / _fluid_viscosity[_qp][ph] -
      89         930 :            _relative_permeability[_qp][ph] * _dfluid_viscosity_dvar[_qp][ph][v] /
      90             :                std::pow(_fluid_viscosity[_qp][ph], 2));
      91             :     }
      92             : 
      93        6080 :   _ddarcy_velocity_dgradvar[_qp].resize(_num_phases);
      94       15030 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
      95             :   {
      96        8950 :     _ddarcy_velocity_dgradvar[_qp][ph].resize(LIBMESH_DIM);
      97       35800 :     for (unsigned i = 0; i < LIBMESH_DIM; ++i)
      98       26850 :       _ddarcy_velocity_dgradvar[_qp][ph][i].resize(_num_var);
      99             :   }
     100             : 
     101       15030 :   for (unsigned ph = 0; ph < _num_phases; ++ph)
     102       35800 :     for (unsigned i = 0; i < LIBMESH_DIM; ++i)
     103       29640 :       for (unsigned v = 0; v < _num_var; ++v)
     104             :       {
     105        2790 :         _ddarcy_velocity_dgradvar[_qp][ph][i][v] =
     106        2790 :             -_dpermeability_dgradvar[_qp][i][v] *
     107        2790 :             (_grad_p[_qp][ph] - _fluid_density[_qp][ph] * _gravity);
     108       11160 :         for (unsigned j = 0; j < LIBMESH_DIM; ++j)
     109        8370 :           _ddarcy_velocity_dgradvar[_qp][ph][i][v](j) -=
     110        8370 :               _permeability[_qp](j, i) * _dgrad_p_dgradvar[_qp][ph][v];
     111             :         _ddarcy_velocity_dgradvar[_qp][ph][i][v] *=
     112        2790 :             _relative_permeability[_qp][ph] / _fluid_viscosity[_qp][ph];
     113             :       }
     114        6080 : }