www.mooseframework.org
MortarConstraint.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
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 "MortarConstraint.h"
11 
12 // MOOSE includes
13 #include "MooseVariable.h"
14 #include "Assembly.h"
15 
16 template <>
19 {
21 }
22 
24  : MortarConstraintBase(parameters),
25  _lambda_dummy(),
26  _lambda(_var ? _var->slnLower() : _lambda_dummy),
27  _u_slave(_slave_var.sln()),
28  _u_master(_master_var.slnNeighbor()),
29  _grad_u_slave(_slave_var.gradSln()),
30  _grad_u_master(_master_var.gradSlnNeighbor()),
31  _phi(nullptr),
32  _grad_phi(nullptr)
33 {
34 }
35 
36 void
38 {
39  unsigned int test_space_size = 0;
40  switch (mortar_type)
41  {
44  test_space_size = _test_slave.size();
45  break;
46 
49  test_space_size = _test_master.size();
50  break;
51 
53  mooseAssert(_var, "LM variable is null");
55  test_space_size = _test.size();
56  break;
57  }
58 
59  for (_qp = 0; _qp < _qrule_msm->n_points(); _qp++)
60  for (_i = 0; _i < test_space_size; _i++)
61  _local_re(_i) += _JxW_msm[_qp] * _coord[_qp] * computeQpResidual(mortar_type);
62 
64 }
65 
66 void
68 {
69  size_t test_space_size = 0;
70  typedef Moose::ConstraintJacobianType JType;
71  typedef Moose::MortarType MType;
72  std::array<JType, 3> jacobian_types;
73 
74  switch (mortar_type)
75  {
76  case MType::Slave:
77  test_space_size = _slave_var.dofIndices().size();
79  break;
80 
81  case MType::Master:
82  test_space_size = _master_var.dofIndicesNeighbor().size();
84  break;
85 
86  case MType::Lower:
87  test_space_size = _var ? _var->dofIndicesLower().size() : 0;
89  break;
90  }
91 
92  std::vector<std::pair<MooseVariableFEBase *, MooseVariableFEBase *>> & ce =
94  for (const auto & it : ce)
95  {
96  MooseVariableFEBase & ivariable = *(it.first);
97  MooseVariableFEBase & jvariable = *(it.second);
98 
99  unsigned int ivar = ivariable.number();
100  unsigned int jvar = jvariable.number();
101 
102  switch (mortar_type)
103  {
104  case MType::Slave:
105  if (ivar != _slave_var.number())
106  continue;
107  break;
108 
109  case MType::Master:
110  if (ivar != _master_var.number())
111  continue;
112  break;
113 
114  case MType::Lower:
115  if (!_var || _var->number() != ivar)
116  continue;
117  break;
118  }
119 
120  std::array<size_t, 3> shape_space_sizes{{jvariable.dofIndices().size(),
121  jvariable.dofIndicesNeighbor().size(),
122  jvariable.dofIndicesLower().size()}};
123  std::array<const VariablePhiValue *, 3> phis;
124  std::array<const VariablePhiGradient *, 3> grad_phis;
125  std::array<const VectorVariablePhiValue *, 3> vector_phis;
126  std::array<const VectorVariablePhiGradient *, 3> vector_grad_phis;
127  if (jvariable.isVector())
128  {
129  const auto & temp_var = static_cast<MooseVariableFE<RealVectorValue> &>(jvariable);
130  vector_phis = {{&temp_var.phiFace(), &temp_var.phiFaceNeighbor(), &temp_var.phiLower()}};
131  vector_grad_phis = {
132  {&temp_var.gradPhiFace(), &temp_var.gradPhiFaceNeighbor(), &temp_var.gradPhiLower()}};
133  }
134  else
135  {
136  const auto & temp_var = static_cast<MooseVariableFE<Real> &>(jvariable);
137  phis = {{&temp_var.phiFace(), &temp_var.phiFaceNeighbor(), &temp_var.phiLower()}};
138  grad_phis = {
139  {&temp_var.gradPhiFace(), &temp_var.gradPhiFaceNeighbor(), &temp_var.gradPhiLower()}};
140  }
141 
142  for (MooseIndex(3) type_index = 0; type_index < 3; ++type_index)
143  {
144  prepareMatrixTagLower(_assembly, ivar, jvar, jacobian_types[type_index]);
145 
147  if (jvariable.isVector())
148  {
149  _vector_phi = vector_phis[type_index];
150  _vector_grad_phi = vector_grad_phis[type_index];
151  }
152  else
153  {
154  _phi = phis[type_index];
155  _grad_phi = grad_phis[type_index];
156  }
157 
158  for (_i = 0; _i < test_space_size; _i++)
159  for (_j = 0; _j < shape_space_sizes[type_index]; _j++)
160  for (_qp = 0; _qp < _qrule_msm->n_points(); _qp++)
161  _local_ke(_i, _j) +=
162  _JxW_msm[_qp] * _coord[_qp] * computeQpJacobian(jacobian_types[type_index], jvar);
164  }
165  }
166 }
User for mortar methods.
const VectorVariablePhiValue * _vector_phi
The current shape functions for vector variables.
InputParameters validParams< MortarConstraintBase >()
void accumulateTaggedLocalResidual()
Local residual blocks will be appended by adding the current local kernel residual.
const VariablePhiValue * _phi
The current shape functions.
const MooseVariable & _slave_var
Reference to the slave variable.
unsigned int number() const
Get variable number coming from libMesh.
virtual Real computeQpJacobian(Moose::ConstraintJacobianType jacobian_type, unsigned int jvar)=0
compute the jacobian at the quadrature points
std::vector< std::pair< MooseVariableFEBase *, MooseVariableFEBase * > > & couplingEntries()
Definition: Assembly.h:762
MortarType
Definition: MooseTypes.h:536
const FieldVariablePhiValue & phiFace() const
const VariablePhiGradient * _grad_phi
The current shape function gradients.
const std::vector< dof_id_type > & dofIndicesLower() const final
Get dof indices for the current lower dimensional element (this is meaningful when performing mortar ...
const VectorVariablePhiGradient * _vector_grad_phi
The current shape function gradients for vector variables.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int _i
Definition: Constraint.h:75
MortarConstraint(const InputParameters &parameters)
DenseMatrix< Number > _local_ke
Holds residual entries as they are accumulated by this Kernel.
const VariableTestValue & _test
The shape functions corresponding to the lagrange multiplier variable.
const VariableTestValue & _test_slave
The shape functions corresponding to the slave interior primal variable.
const MooseVariable * _var
Pointer to the lagrange multipler variable. nullptr if none.
InputParameters validParams< MortarConstraint >()
void computeJacobian(Moose::MortarType mortar_type) override
compute the residual for the specified element type
unsigned int _j
Definition: Constraint.h:75
const std::vector< dof_id_type > & dofIndices() const final
Get local DoF indices.
void computeResidual(Moose::MortarType mortar_type) override
compute the residual for the specified element type
virtual const std::vector< dof_id_type > & dofIndices() const
Get local DoF indices.
const std::vector< Real > & _JxW_msm
The element Jacobian times weights.
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
const MooseVariable & _master_var
Reference to the master variable.
const VariableTestValue & _test_master
The shape functions corresponding to the master interior primal variable.
virtual const std::vector< dof_id_type > & dofIndicesLower() const =0
Get dof indices for the current lower dimensional element (this is meaningful when performing mortar ...
Assembly & _assembly
Definition: Constraint.h:72
const QBase *const & _qrule_msm
The quadrature rule.
void prepareVectorTagNeighbor(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual the according to active tags for DG and interface kernels...
virtual bool isVector() const =0
const std::vector< dof_id_type > & dofIndicesNeighbor() const final
Get neighbor DOF indices for currently selected element.
ConstraintJacobianType
Definition: MooseTypes.h:543
void prepareMatrixTagLower(Assembly &assembly, unsigned int ivar, unsigned int jvar, Moose::ConstraintJacobianType type)
Prepare data for computing the jacobian according to the ative tags for mortar.
DenseVector< Number > _local_re
Holds residual entries as they are accumulated by this Kernel.
virtual const std::vector< dof_id_type > & dofIndicesNeighbor() const =0
Get neighbor DOF indices for currently selected element.
void prepareVectorTagLower(Assembly &assembly, unsigned int ivar)
Prepare data for computing the residual according to active tags for mortar constraints.
void prepareVectorTag(Assembly &assembly, unsigned int ivar)
Prepare data for computing element residual the according to active tags.
const MooseArray< Real > & _coord
Member for handling change of coordinate systems (xyz, rz, spherical)
unsigned int _qp
Definition: Constraint.h:76
virtual Real computeQpResidual(Moose::MortarType mortar_type)=0
compute the residual at the quadrature points