Adding test case for 1D mortar method

Closes idaholab#4211

framework/include/bcs/OneDEqualValueConstraintBC.h

@@ -0,0 +1,47 @@
+/****************************************************************/
+/*               DO NOT MODIFY THIS HEADER                      */
+/* MOOSE - Multiphysics Object Oriented Simulation Environment  */
+/*                                                              */
+/*           (c) 2010 Battelle Energy Alliance, LLC             */
+/*                   ALL RIGHTS RESERVED                        */
+/*                                                              */
+/*          Prepared by Battelle Energy Alliance, LLC           */
+/*            Under Contract No. DE-AC07-05ID14517              */
+/*            With the U. S. Department of Energy               */
+/*                                                              */
+/*            See COPYRIGHT for full restrictions               */
+/****************************************************************/
+
+#ifndef ONEDEQUALVALUECONSTRAINTBC_H
+#define ONEDEQUALVALUECONSTRAINTBC_H
+
+#include "IntegratedBC.h"
+
+class OneDEqualValueConstraintBC;
+
+template<>
+InputParameters validParams<OneDEqualValueConstraintBC>();
+
+
+/**
+ * This is the \int \lambda \d g term from the mortar method. This can connect two 1D domains only.
+ *
+ * For higher dimensions, you should use face-face constraints.
+ */
+class OneDEqualValueConstraintBC : public IntegratedBC
+{
+public:
+  OneDEqualValueConstraintBC(const std::string & name, InputParameters parameters);
+  virtual ~OneDEqualValueConstraintBC();
+
+protected:
+  virtual Real computeQpResidual();
+  virtual Real computeQpJacobian();
+  virtual Real computeQpOffDiagJacobian(unsigned jvar);
+
+  VariableValue & _lambda;
+  unsigned int _lambda_var_number;
+  unsigned int _component;
+};
+
+#endif // ONEDEQUALVALUECONSTRAINTBC_H

framework/include/kernels/NodalEqualValueConstraint.h

@@ -0,0 +1,45 @@
+/****************************************************************/
+/*               DO NOT MODIFY THIS HEADER                      */
+/* MOOSE - Multiphysics Object Oriented Simulation Environment  */
+/*                                                              */
+/*           (c) 2010 Battelle Energy Alliance, LLC             */
+/*                   ALL RIGHTS RESERVED                        */
+/*                                                              */
+/*          Prepared by Battelle Energy Alliance, LLC           */
+/*            Under Contract No. DE-AC07-05ID14517              */
+/*            With the U. S. Department of Energy               */
+/*                                                              */
+/*            See COPYRIGHT for full restrictions               */
+/****************************************************************/
+
+#ifndef NODALEQUALVALUECONSTRAINT_H
+#define NODALEQUALVALUECONSTRAINT_H
+
+#include "NodalScalarKernel.h"
+
+class NodalEqualValueConstraint;
+
+template<>
+InputParameters validParams<NodalEqualValueConstraint>();
+
+/**
+ * Constraint to enforce equal values (in 1D)
+ */
+class NodalEqualValueConstraint : public NodalScalarKernel
+{
+public:
+  NodalEqualValueConstraint(const std::string & name, InputParameters parameters);
+  virtual ~NodalEqualValueConstraint();
+
+  virtual void computeResidual();
+  virtual void computeJacobian();
+
+protected:
+  std::vector<Real> _normals;
+
+  std::vector<unsigned int> _val_number;
+  std::vector<VariableValue *> _value;
+};
+
+
+#endif /* NODALEQUALVALUECONSTRAINT_H */

framework/src/base/Moose.C

@@ -70,6 +70,7 @@
 #include "WeakGradientBC.h"
 #include "DiffusionFluxBC.h"
 #include "PostprocessorDirichletBC.h"
+#include "OneDEqualValueConstraintBC.h"
 
 // auxkernels
 #include "ConstantAux.h"
@@ -233,6 +234,7 @@
 // ScalarKernels
 #include "ODETimeDerivative.h"
 #include "FunctionScalarAux.h"
+#include "NodalEqualValueConstraint.h"
 
 // indicators
 #include "AnalyticalIndicator.h"
@@ -420,6 +422,7 @@ registerObjects(Factory & factory)
   registerBoundaryCondition(WeakGradientBC);
   registerBoundaryCondition(DiffusionFluxBC);
   registerBoundaryCondition(PostprocessorDirichletBC);
+  registerBoundaryCondition(OneDEqualValueConstraintBC);
 
   // dirac kernels
   registerDiracKernel(ConstantPointSource);
@@ -579,6 +582,7 @@ registerObjects(Factory & factory)
 
   // Scalar kernels
   registerScalarKernel(ODETimeDerivative);
+  registerScalarKernel(NodalEqualValueConstraint);
 
   // indicators
   registerIndicator(AnalyticalIndicator);

framework/src/bcs/OneDEqualValueConstraintBC.C

@@ -0,0 +1,63 @@
+/****************************************************************/
+/*               DO NOT MODIFY THIS HEADER                      */
+/* MOOSE - Multiphysics Object Oriented Simulation Environment  */
+/*                                                              */
+/*           (c) 2010 Battelle Energy Alliance, LLC             */
+/*                   ALL RIGHTS RESERVED                        */
+/*                                                              */
+/*          Prepared by Battelle Energy Alliance, LLC           */
+/*            Under Contract No. DE-AC07-05ID14517              */
+/*            With the U. S. Department of Energy               */
+/*                                                              */
+/*            See COPYRIGHT for full restrictions               */
+/****************************************************************/
+
+#include "OneDEqualValueConstraintBC.h"
+
+template<>
+InputParameters validParams<OneDEqualValueConstraintBC>()
+{
+  InputParameters params = validParams<IntegratedBC>();
+  params.addRequiredCoupledVar("lambda", "Lagrange multiplier");
+  params.addRequiredParam<unsigned int>("component", "Component of the Lagrange multiplier");
+  return params;
+}
+
+
+OneDEqualValueConstraintBC::OneDEqualValueConstraintBC(const std::string & name, InputParameters parameters) :
+    IntegratedBC(name, parameters),
+    _lambda(coupledScalarValue("lambda")),
+    _lambda_var_number(coupledScalar("lambda")),
+    _component(getParam<unsigned int>("component"))
+{
+}
+
+OneDEqualValueConstraintBC::~OneDEqualValueConstraintBC()
+{
+}
+
+Real
+OneDEqualValueConstraintBC::computeQpResidual()
+{
+  return _lambda[_component] * _normals[_qp](0) * _test[_i][_qp];
+}
+
+Real
+OneDEqualValueConstraintBC::computeQpJacobian()
+{
+  return 0.;
+}
+
+Real
+OneDEqualValueConstraintBC::computeQpOffDiagJacobian(unsigned jvar)
+{
+  if (jvar == _lambda_var_number)
+  {
+    if (_j == _component)
+      return _normals[_qp](0) * _test[_i][_qp];
+    else
+      return 0.;
+  }
+  else
+    return 0.;
+}

framework/src/kernels/NodalEqualValueConstraint.C

@@ -0,0 +1,74 @@
+/****************************************************************/
+/*               DO NOT MODIFY THIS HEADER                      */
+/* MOOSE - Multiphysics Object Oriented Simulation Environment  */
+/*                                                              */
+/*           (c) 2010 Battelle Energy Alliance, LLC             */
+/*                   ALL RIGHTS RESERVED                        */
+/*                                                              */
+/*          Prepared by Battelle Energy Alliance, LLC           */
+/*            Under Contract No. DE-AC07-05ID14517              */
+/*            With the U. S. Department of Energy               */
+/*                                                              */
+/*            See COPYRIGHT for full restrictions               */
+/****************************************************************/
+
+#include "NodalEqualValueConstraint.h"
+
+template<>
+InputParameters validParams<NodalEqualValueConstraint>()
+{
+  InputParameters params = validParams<NodalScalarKernel>();
+  params.addRequiredParam<std::vector<Real> >("normals", "node normals");
+  params.addRequiredCoupledVar("var", "Variables to put the constraint on");
+  return params;
+}
+
+NodalEqualValueConstraint::NodalEqualValueConstraint(const std::string & name, InputParameters parameters) :
+    NodalScalarKernel(name, parameters),
+    _normals(getParam<std::vector<Real> >("normals"))
+{
+  if (_normals.size() != 2)
+    mooseError("Component '" << this->name() << "' can only have 2 connections.");
+
+  unsigned int n = coupledComponents("var");
+  _value.resize(n);
+  _val_number.resize(n);
+  for (unsigned int k = 0; k < n; k++)
+  {
+    _value[k] = &coupledValue("var", k);
+    _val_number[k] = coupled("var", k);
+  }
+}
+
+NodalEqualValueConstraint::~NodalEqualValueConstraint()
+{
+}
+
+void
+NodalEqualValueConstraint::computeResidual()
+{
+  // LM residuals
+  DenseVector<Number> & lm_re = _assembly.residualBlock(_var.number());
+
+  for (unsigned int k = 0; k < _value.size(); k++)
+    lm_re(k) = (*_value[k])[0] * _normals[0] + (*_value[k])[1] * _normals[1];
+}
+
+void
+NodalEqualValueConstraint::computeJacobian()
+{
+  // do the diagonal block
+  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
+  // put zeroes on the diagonal (we have to do it, otherwise PETSc will complain!)
+  for (unsigned int i = 0; i < ke.m(); i++)
+    for (unsigned int j = 0; j < ke.n(); j++)
+      ke(i, j) = 0.;
+
+  for (unsigned int k = 0; k < _value.size(); k++)
+  {
+    DenseMatrix<Number> & ken = _assembly.jacobianBlock(_var.number(), _val_number[k]);
+
+    ken(k, 0) = _normals[0];
+    ken(k, 1) = _normals[1];
+  }
+}

test/tests/mortar/1d/1d.i

@@ -0,0 +1,79 @@
+[Mesh]
+  file = 2-lines.e
+  construct_side_list_from_node_list = true
+[]
+
+[Variables]
+  [./u]
+    order = FIRST
+    family = LAGRANGE
+    block = '1 2'
+  [../]
+
+  [./lm]
+    order = FIRST
+    family = SCALAR
+  [../]
+[]
+
+[Kernels]
+  [./diff]
+    type = Diffusion
+    variable = u
+  [../]
+[]
+
+[ScalarKernels]
+  [./ced]
+    type = NodalEqualValueConstraint
+    variable = lm
+    var = u
+    nodes = '2 3'
+    normals = '1 -1'
+  [../]
+[]
+
+[BCs]
+  [./left]
+    type = DirichletBC
+    variable = u
+    boundary = '1'
+    value = 1
+  [../]
+
+  [./right]
+    type = DirichletBC
+    variable = u
+    boundary = '2'
+    value = 3
+  [../]
+
+  [./evc]
+    type = OneDEqualValueConstraintBC
+    variable = u
+    boundary = '100 101'
+    lambda = lm
+    component = 0
+  [../]
+[]
+
+[Preconditioning]
+  [./fmp]
+    type = SMP
+    full = true
+    solve_type = 'NEWTON'
+  [../]
+[]
+
+[Executioner]
+  type = Steady
+[]
+
+[Outputs]
+  output_initial = true
+  exodus = true
+  [./console]
+    type = Console
+    perf_log = true
+  [../]
+[]

test/tests/mortar/1d/2-lines.jou

@@ -0,0 +1,29 @@
+create vertex 0    0 0
+create vertex 1    0 0
+create curve vertex 1 2
+
+create vertex 1.1 0 0
+create vertex 2.1 0 0
+create curve vertex 3 4
+
+curve 1 interval 2
+curve 1 scheme bias factor 1.0
+mesh curve 1
+
+curve 2 interval 2
+curve 2 scheme bias factor 1.0
+mesh curve 2
+
+set duplicate block elements off
+block 1 curve 1
+block 1 element type BEAM2
+
+set duplicate block elements off
+block 2 curve 2
+block 2 element type BEAM2
+
+nodeset 1 vertex 1
+nodeset 2 vertex 4
+
+nodeset 100 vertex 2
+nodeset 101 vertex 3

test/tests/mortar/1d/tests

@@ -0,0 +1,10 @@
+[Tests]
+  [./test]
+    type = 'Exodiff'
+    input = '1d.i'
+    exodiff = '1d_out.e'
+    max_parallel = 1
+    max_threads = 1
+    compiler = 'GCC CLANG'
+  [../]
+[]