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 "IntegratedBC.h"
      11             : 
      12             : // MOOSE includes
      13             : #include "Assembly.h"
      14             : #include "SubProblem.h"
      15             : #include "SystemBase.h"
      16             : #include "MooseVariableFE.h"
      17             : #include "MooseVariableScalar.h"
      18             : #include "FEProblemBase.h"
      19             : 
      20             : #include "libmesh/quadrature.h"
      21             : 
      22             : InputParameters
      23      540683 : IntegratedBC::validParams()
      24             : {
      25      540683 :   InputParameters params = IntegratedBCBase::validParams();
      26      540683 :   return params;
      27             : }
      28             : 
      29        6682 : IntegratedBC::IntegratedBC(const InputParameters & parameters)
      30             :   : IntegratedBCBase(parameters),
      31             :     MooseVariableInterface<Real>(this,
      32             :                                  false,
      33             :                                  "variable",
      34             :                                  Moose::VarKindType::VAR_SOLVER,
      35             :                                  Moose::VarFieldType::VAR_FIELD_STANDARD),
      36       13354 :     _var(*mooseVariable()),
      37        6677 :     _normals(_assembly.normals()),
      38        6677 :     _phi(_assembly.phiFace(_var)),
      39        6677 :     _grad_phi(_assembly.gradPhiFace(_var)),
      40        6677 :     _test(_var.phiFace()),
      41        6677 :     _grad_test(_var.gradPhiFace()),
      42        6677 :     _u(_is_implicit ? _var.sln() : _var.slnOld()),
      43       13359 :     _grad_u(_is_implicit ? _var.gradSln() : _var.gradSlnOld())
      44             : {
      45        6677 :   addMooseVariableDependency(mooseVariable());
      46             : 
      47        6677 :   _save_in.resize(_save_in_strings.size());
      48        6677 :   _diag_save_in.resize(_diag_save_in_strings.size());
      49             : 
      50        6775 :   for (unsigned int i = 0; i < _save_in_strings.size(); i++)
      51             :   {
      52          98 :     MooseVariable * var = &_subproblem.getStandardVariable(_tid, _save_in_strings[i]);
      53             : 
      54          98 :     if (var->feType() != _var.feType())
      55           0 :       paramError(
      56             :           "save_in",
      57             :           "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
      58           0 :           moose::internal::incompatVarMsg(*var, _var));
      59          98 :     _save_in[i] = var;
      60          98 :     var->sys().addVariableToZeroOnResidual(_save_in_strings[i]);
      61          98 :     addMooseVariableDependency(var);
      62             :   }
      63             : 
      64        6677 :   _has_save_in = _save_in.size() > 0;
      65             : 
      66        6726 :   for (unsigned int i = 0; i < _diag_save_in_strings.size(); i++)
      67             :   {
      68          49 :     MooseVariable * var = &_subproblem.getStandardVariable(_tid, _diag_save_in_strings[i]);
      69             : 
      70          49 :     if (var->feType() != _var.feType())
      71           0 :       paramError(
      72             :           "diag_save_in",
      73             :           "saved-in auxiliary variable is incompatible with the object's nonlinear variable: ",
      74           0 :           moose::internal::incompatVarMsg(*var, _var));
      75             : 
      76          49 :     _diag_save_in[i] = var;
      77          49 :     var->sys().addVariableToZeroOnJacobian(_diag_save_in_strings[i]);
      78          49 :     addMooseVariableDependency(var);
      79             :   }
      80             : 
      81        6677 :   _has_diag_save_in = _diag_save_in.size() > 0;
      82        6677 : }
      83             : 
      84             : void
      85     1932922 : IntegratedBC::computeResidual()
      86             : {
      87     1932922 :   prepareVectorTag(_assembly, _var.number());
      88             : 
      89     1932922 :   precalculateResidual();
      90             : 
      91     7331447 :   for (_qp = 0; _qp < _qrule->n_points(); _qp++)
      92             :   {
      93     5398525 :     precalculateQpResidual();
      94    52183169 :     for (_i = 0; _i < _test.size(); _i++)
      95    46784644 :       _local_re(_i) += _JxW[_qp] * _coord[_qp] * computeQpResidual();
      96             :   }
      97             : 
      98     1932922 :   accumulateTaggedLocalResidual();
      99             : 
     100     1932922 :   if (_has_save_in)
     101             :   {
     102        1254 :     libMesh::Threads::spin_mutex::scoped_lock lock(libMesh::Threads::spin_mtx);
     103        2668 :     for (unsigned int i = 0; i < _save_in.size(); i++)
     104        1414 :       _save_in[i]->sys().solution().add_vector(_local_re, _save_in[i]->dofIndices());
     105        1254 :   }
     106     1932922 : }
     107             : 
     108             : void
     109      215670 : IntegratedBC::computeJacobian()
     110             : {
     111      215670 :   prepareMatrixTag(_assembly, _var.number(), _var.number());
     112             : 
     113      215670 :   precalculateJacobian();
     114             : 
     115      889638 :   for (_qp = 0; _qp < _qrule->n_points(); _qp++)
     116             :   {
     117      673968 :     precalculateQpJacobian();
     118     6474553 :     for (_i = 0; _i < _test.size(); _i++)
     119    86764464 :       for (_j = 0; _j < _phi.size(); _j++)
     120    80963879 :         _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpJacobian();
     121             :   }
     122             : 
     123      215670 :   accumulateTaggedLocalMatrix();
     124             : 
     125      215670 :   if (_has_diag_save_in)
     126             :   {
     127         158 :     unsigned int rows = _local_ke.m();
     128         158 :     DenseVector<Number> diag(rows);
     129         790 :     for (unsigned int i = 0; i < rows; i++)
     130         632 :       diag(i) = _local_ke(i, i);
     131             : 
     132         158 :     libMesh::Threads::spin_mutex::scoped_lock lock(libMesh::Threads::spin_mtx);
     133         316 :     for (unsigned int i = 0; i < _diag_save_in.size(); i++)
     134         158 :       _diag_save_in[i]->sys().solution().add_vector(diag, _diag_save_in[i]->dofIndices());
     135         158 :   }
     136      215670 : }
     137             : 
     138             : void
     139      108916 : IntegratedBC::computeOffDiagJacobian(const unsigned int jvar_num)
     140             : {
     141      108916 :   const auto & jvar = getVariable(jvar_num);
     142             : 
     143      108916 :   if (jvar_num == _var.number())
     144             :   {
     145       97106 :     computeJacobian();
     146       97106 :     return;
     147             :   }
     148             : 
     149       11810 :   prepareMatrixTag(_assembly, _var.number(), jvar_num);
     150             : 
     151       11810 :   precalculateOffDiagJacobian(jvar_num);
     152             : 
     153             :   // This (undisplaced) jvar could potentially yield the wrong phi size if this object is acting
     154             :   // on the displaced mesh, so we obtain the variable on the proper system
     155       11810 :   auto phi_size = jvar.dofIndices().size();
     156             :   mooseAssert(phi_size * jvar.count() == _local_ke.n(),
     157             :               "The size of the phi container does not match the number of local Jacobian columns");
     158             : 
     159       46108 :   for (_qp = 0; _qp < _qrule->n_points(); _qp++)
     160             :   {
     161       34298 :     precalculateQpOffDiagJacobian(jvar);
     162      360446 :     for (_i = 0; _i < _test.size(); _i++)
     163      877472 :       for (_j = 0; _j < phi_size; _j++)
     164      551324 :         _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar_num);
     165             :   }
     166             : 
     167       11810 :   accumulateTaggedLocalMatrix();
     168             : }
     169             : 
     170             : void
     171          64 : IntegratedBC::computeOffDiagJacobianScalar(const unsigned int jvar)
     172             : {
     173          64 :   prepareMatrixTag(_assembly, _var.number(), jvar);
     174             : 
     175          64 :   MooseVariableScalar & jv = _sys.getScalarVariable(_tid, jvar);
     176         176 :   for (_qp = 0; _qp < _qrule->n_points(); _qp++)
     177         528 :     for (_i = 0; _i < _test.size(); _i++)
     178         832 :       for (_j = 0; _j < jv.order(); _j++)
     179         416 :         _local_ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobianScalar(jvar);
     180             : 
     181          64 :   accumulateTaggedLocalMatrix();
     182          64 : }
     183             : 
     184             : void
     185         184 : IntegratedBC::computeResidualAndJacobian()
     186             : {
     187         184 :   computeResidual();
     188             : 
     189         824 :   for (const auto & [ivariable, jvariable] : _fe_problem.couplingEntries(_tid, _sys.number()))
     190             :   {
     191         640 :     const unsigned int ivar = ivariable->number();
     192         640 :     const unsigned int jvar = jvariable->number();
     193             : 
     194         640 :     if (ivar != _var.number())
     195         304 :       continue;
     196             : 
     197         336 :     if (_is_implicit)
     198             :     {
     199         336 :       prepareShapes(jvar);
     200         336 :       computeOffDiagJacobian(jvar);
     201             :     }
     202             :   }
     203             : 
     204             :   /// TODO: add nonlocal Jacobians and scalar Jacobians
     205         184 : }