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             : #pragma once
      11             : 
      12             : // MOOSE includes
      13             : #include "Constraint.h"
      14             : #include "NeighborCoupleableMooseVariableDependencyIntermediateInterface.h"
      15             : 
      16             : // libMesh forward declarations
      17             : namespace libMesh
      18             : {
      19             : template <typename T>
      20             : class SparseMatrix;
      21             : }
      22             : 
      23             : /**
      24             :  * A NodeFaceConstraint is used when you need to create constraints between
      25             :  * two surfaces in a mesh.  It works by allowing you to modify the residual
      26             :  * and jacobian entries on "this" side (the node side, also referred to as
      27             :  * the secondary side) and the "other" side (the face side, also referred to as
      28             :  * the primary side)
      29             :  *
      30             :  * This is common for contact algorithms and other constraints.
      31             :  */
      32             : class NodeFaceConstraint : public Constraint,
      33             :                            public NeighborCoupleableMooseVariableDependencyIntermediateInterface,
      34             :                            public NeighborMooseVariableInterface<Real>
      35             : {
      36             : public:
      37             :   static InputParameters validParams();
      38             : 
      39             :   NodeFaceConstraint(const InputParameters & parameters);
      40             :   virtual ~NodeFaceConstraint();
      41             : 
      42             :   /**
      43             :    * Compute the value the secondary node should have at the beginning of a timestep.
      44             :    */
      45             :   virtual void computeSecondaryValue(NumericVector<Number> & current_solution);
      46             : 
      47             :   /**
      48             :    * Computes the residual Nodal residual.
      49             :    */
      50             :   virtual void computeResidual() override;
      51             : 
      52             :   /**
      53             :    * Computes the jacobian for the current element.
      54             :    */
      55             :   virtual void computeJacobian() override;
      56             : 
      57             :   /**
      58             :    * Computes d-residual / d-jvar...
      59             :    */
      60             :   virtual void computeOffDiagJacobian(unsigned int jvar) override;
      61             : 
      62             :   /**
      63             :    * Gets the indices for all dofs connected to the constraint
      64             :    */
      65             :   virtual void getConnectedDofIndices(unsigned int var_num);
      66             : 
      67             :   /**
      68             :    * Whether or not this constraint should be applied.
      69             :    *
      70             :    * Get's called once per secondary node.
      71             :    */
      72       44852 :   virtual bool shouldApply() { return true; }
      73             : 
      74             :   /**
      75             :    * Whether or not the secondary's residual should be overwritten.
      76             :    *
      77             :    * When this returns true the secondary's residual as computed by the constraint will _replace_
      78             :    * the residual previously at that node for that variable.
      79             :    */
      80             :   virtual bool overwriteSecondaryResidual();
      81             : 
      82             :   /**
      83             :    * Whether or not the secondary's Jacobian row should be overwritten.
      84             :    *
      85             :    * When this returns true the secondary's Jacobian row as computed by the constraint will
      86             :    * _replace_ the residual previously at that node for that variable.
      87             :    */
      88        7392 :   virtual bool overwriteSecondaryJacobian() { return overwriteSecondaryResidual(); };
      89             : 
      90             :   /**
      91             :    * The variable on the Primary side of the domain.
      92             :    */
      93       11088 :   virtual MooseVariable & primaryVariable() { return _primary_var; }
      94             : 
      95             :   /**
      96             :    * The primary boundary ID for this constraint.
      97             :    */
      98       96052 :   BoundaryID primaryBoundary() const { return _primary; }
      99             : 
     100             :   /**
     101             :    * The secondary boundary ID for this constraint.
     102             :    */
     103       96052 :   BoundaryID secondaryBoundary() const { return _secondary; }
     104             : 
     105             :   /**
     106             :    * The variable number that this object operates on.
     107             :    */
     108       58012 :   const MooseVariable & variable() const override { return _var; }
     109             : 
     110             :   Real secondaryResidual() const;
     111             : 
     112             :   void residualSetup() override;
     113             : 
     114             :   /**
     115             :    * @returns material property dependencies
     116             :    */
     117             :   virtual const std::unordered_set<unsigned int> & getMatPropDependencies() const;
     118             : 
     119             :   /**
     120             :    * Whether (contact) constraint is of 'explicit dynamics' type.
     121             :    */
     122        5576 :   virtual bool isExplicitConstraint() const { return false; }
     123             : 
     124             :   /**
     125             :    * Allows for overwriting boundary variables (explicit dynamics contact).
     126             :    */
     127           0 :   virtual void overwriteBoundaryVariables(NumericVector<Number> & /*soln*/,
     128             :                                           const Node & /*secondary_node*/) const
     129             :   {
     130           0 :   }
     131             : 
     132             :   /**
     133             :    * @returns IDs of the subdomains that connect to the secondary boundary
     134             :    */
     135             :   std::set<SubdomainID> getSecondaryConnectedBlocks() const;
     136             : 
     137             : protected:
     138             :   // TODO: Make this protected or add an accessor
     139             :   // Do the same for all the other public members
     140             :   const SparseMatrix<Number> * _jacobian;
     141             : 
     142             :   /**
     143             :    * Compute the value the secondary node should have at the beginning of a timestep.
     144             :    */
     145             :   virtual Real computeQpSecondaryValue() = 0;
     146             : 
     147             :   /**
     148             :    * This is the virtual that derived classes should override for computing the residual on
     149             :    * neighboring element.
     150             :    */
     151             :   virtual Real computeQpResidual(Moose::ConstraintType type) = 0;
     152             : 
     153             :   /**
     154             :    * This is the virtual that derived classes should override for computing the Jacobian on
     155             :    * neighboring element.
     156             :    */
     157             :   virtual Real computeQpJacobian(Moose::ConstraintJacobianType type) = 0;
     158             : 
     159             :   /**
     160             :    * This is the virtual that derived classes should override for computing the off-diag Jacobian.
     161             :    */
     162           0 :   virtual Real computeQpOffDiagJacobian(Moose::ConstraintJacobianType /*type*/,
     163             :                                         unsigned int /*jvar*/)
     164             :   {
     165           0 :     return 0;
     166             :   }
     167             : 
     168             :   /// coupling interface:
     169           0 :   virtual const VariableValue & coupledSecondaryValue(const std::string & var_name,
     170             :                                                       unsigned int comp = 0)
     171             :   {
     172           0 :     return coupledValue(var_name, comp);
     173             :   }
     174           0 :   virtual const VariableValue & coupledSecondaryValueOld(const std::string & var_name,
     175             :                                                          unsigned int comp = 0)
     176             :   {
     177           0 :     return coupledValueOld(var_name, comp);
     178             :   }
     179           0 :   virtual const VariableValue & coupledSecondaryValueOlder(const std::string & var_name,
     180             :                                                            unsigned int comp = 0)
     181             :   {
     182           0 :     return coupledValueOlder(var_name, comp);
     183             :   }
     184             : 
     185           0 :   virtual const VariableGradient & coupledSecondaryGradient(const std::string & var_name,
     186             :                                                             unsigned int comp = 0)
     187             :   {
     188           0 :     return coupledGradient(var_name, comp);
     189             :   }
     190           0 :   virtual const VariableGradient & coupledSecondaryGradientOld(const std::string & var_name,
     191             :                                                                unsigned int comp = 0)
     192             :   {
     193           0 :     return coupledGradientOld(var_name, comp);
     194             :   }
     195           0 :   virtual const VariableGradient & coupledSecondaryGradientOlder(const std::string & var_name,
     196             :                                                                  unsigned int comp = 0)
     197             :   {
     198           0 :     return coupledGradientOlder(var_name, comp);
     199             :   }
     200             : 
     201           0 :   virtual const VariableSecond & coupledSecondarySecond(const std::string & var_name,
     202             :                                                         unsigned int comp = 0)
     203             :   {
     204           0 :     return coupledSecond(var_name, comp);
     205             :   }
     206             : 
     207           0 :   virtual const VariableValue & coupledPrimaryValue(const std::string & var_name,
     208             :                                                     unsigned int comp = 0)
     209             :   {
     210           0 :     return coupledNeighborValue(var_name, comp);
     211             :   }
     212           0 :   virtual const VariableValue & coupledPrimaryValueOld(const std::string & var_name,
     213             :                                                        unsigned int comp = 0)
     214             :   {
     215           0 :     return coupledNeighborValueOld(var_name, comp);
     216             :   }
     217           0 :   virtual const VariableValue & coupledPrimaryValueOlder(const std::string & var_name,
     218             :                                                          unsigned int comp = 0)
     219             :   {
     220           0 :     return coupledNeighborValueOlder(var_name, comp);
     221             :   }
     222             : 
     223           0 :   virtual const VariableGradient & coupledPrimaryGradient(const std::string & var_name,
     224             :                                                           unsigned int comp = 0)
     225             :   {
     226           0 :     return coupledNeighborGradient(var_name, comp);
     227             :   }
     228           0 :   virtual const VariableGradient & coupledPrimaryGradientOld(const std::string & var_name,
     229             :                                                              unsigned int comp = 0)
     230             :   {
     231           0 :     return coupledNeighborGradientOld(var_name, comp);
     232             :   }
     233           0 :   virtual const VariableGradient & coupledPrimaryGradientOlder(const std::string & var_name,
     234             :                                                                unsigned int comp = 0)
     235             :   {
     236           0 :     return coupledNeighborGradientOlder(var_name, comp);
     237             :   }
     238             : 
     239           0 :   virtual const VariableSecond & coupledPrimarySecond(const std::string & var_name,
     240             :                                                       unsigned int comp = 0)
     241             :   {
     242           0 :     return coupledNeighborSecond(var_name, comp);
     243             :   }
     244             : 
     245             :   /**
     246             :    * Builds the \p _boundary_ids data member and returns it
     247             :    * @returns a set containing the secondary and primary boundary IDs
     248             :    */
     249             :   const std::set<BoundaryID> & buildBoundaryIDs();
     250             : 
     251             :   /// Boundary ID for the secondary surface
     252             :   BoundaryID _secondary;
     253             :   /// Boundary ID for the primary surface
     254             :   BoundaryID _primary;
     255             : 
     256             :   MooseVariable & _var;
     257             : 
     258             :   const MooseArray<Point> & _primary_q_point;
     259             :   const QBase * const & _primary_qrule;
     260             : 
     261             :   /// the union of the secondary and primary boundary ids
     262             :   std::set<BoundaryID> _boundary_ids;
     263             : 
     264             :   PenetrationLocator & _penetration_locator;
     265             : 
     266             :   /// current node being processed
     267             :   const Node * const & _current_node;
     268             :   const Elem * const & _current_primary;
     269             : 
     270             :   /// Value of the unknown variable this BC is action on
     271             :   const VariableValue & _u_secondary;
     272             :   /// Shape function on the secondary side.  This will always
     273             :   VariablePhiValue _phi_secondary;
     274             :   /// Shape function on the secondary side.  This will always only have one entry and that entry will always be "1"
     275             :   VariableTestValue _test_secondary;
     276             : 
     277             :   /// Primary side variable
     278             :   MooseVariable & _primary_var;
     279             : 
     280             :   /// Number for the primary variable
     281             :   unsigned int _primary_var_num;
     282             : 
     283             :   /// Side shape function.
     284             :   const VariablePhiValue & _phi_primary;
     285             :   /// Gradient of side shape function
     286             :   const VariablePhiGradient & _grad_phi_primary;
     287             : 
     288             :   /// Side test function
     289             :   const VariableTestValue & _test_primary;
     290             :   /// Gradient of side shape function
     291             :   const VariableTestGradient & _grad_test_primary;
     292             : 
     293             :   /// Holds the current solution at the current quadrature point
     294             :   const VariableValue & _u_primary;
     295             :   /// Holds the current solution gradient at the current quadrature point
     296             :   const VariableGradient & _grad_u_primary;
     297             : 
     298             :   /// DOF map
     299             :   const DofMap & _dof_map;
     300             : 
     301             :   const std::map<dof_id_type, std::vector<dof_id_type>> & _node_to_elem_map;
     302             : 
     303             :   /**
     304             :    * Whether or not the secondary's residual should be overwritten.
     305             :    *
     306             :    * When this is true the secondary's residual as computed by the constraint will _replace_
     307             :    * the residual previously at that node for that variable.
     308             :    */
     309             :   bool _overwrite_secondary_residual;
     310             : 
     311             :   /// JxW on the primary face
     312             :   const MooseArray<Real> & _primary_JxW;
     313             : 
     314             :   /// Whether the secondary residual has been computed
     315             :   bool _secondary_residual_computed;
     316             : 
     317             :   /// The value of the secondary residual
     318             :   Real _secondary_residual;
     319             : 
     320             :   /// An empty material property dependency set for use with \p getMatPropDependencies
     321             :   const std::unordered_set<unsigned int> _empty_mat_prop_deps;
     322             : 
     323             :   std::vector<dof_id_type> _connected_dof_indices;
     324             : 
     325             :   /// The Jacobian corresponding to the derivatives of the neighbor/primary residual with respect to
     326             :   /// the elemental/secondary degrees of freedom.  We want to manually manipulate Kne because of the
     327             :   /// dependence of the primary residuals on dofs from all elements connected to the secondary node
     328             :   /// (e.g. those held by _connected_dof_indices)
     329             :   DenseMatrix<Number> _Kne;
     330             : 
     331             :   /// The Jacobian corresponding to the derivatives of the elemental/secondary residual with respect to
     332             :   /// the elemental/secondary degrees of freedom.  We want to manually manipulate Kee because of the
     333             :   /// dependence of the secondary/primary residuals on // dofs from all elements connected to the secondary
     334             :   /// node (e.g. those held by _connected_dof_indices) // and because when we're overwriting the
     335             :   /// secondary residual we traditionally want to use a different // scaling factor from the one
     336             :   /// associated with interior physics
     337             :   DenseMatrix<Number> _Kee;
     338             : 
     339             :   /// The Jacobian corresponding to the derivatives of the elemental/secondary residual with respect to
     340             :   /// the neighbor/primary degrees of freedom.  We want to manually manipulate Ken because when we're
     341             :   /// overwriting the secondary residual we traditionally want to use a different scaling factor from the
     342             :   /// one associated with interior physics
     343             :   DenseMatrix<Number> _Ken;
     344             : 
     345             :   friend class NonlinearSystemBase;
     346             : };
     347             : 
     348             : inline const std::unordered_set<unsigned int> &
     349        2725 : NodeFaceConstraint::getMatPropDependencies() const
     350             : {
     351        2725 :   return _empty_mat_prop_deps;
     352             : }