www.mooseframework.org
Public Member Functions | List of all members
SparsityBasedContactConstraint Class Reference

#include <SparsityBasedContactConstraint.h>

Inheritance diagram for SparsityBasedContactConstraint:
[legend]

Public Member Functions

 SparsityBasedContactConstraint (const InputParameters &parameters)
 
virtual ~SparsityBasedContactConstraint ()
 
virtual Real computeQpSlaveValue ()
 
virtual Real computeQpResidual (Moose::ConstraintType)
 
virtual Real computeQpJacobian (Moose::ConstraintJacobianType)
 
virtual void getConnectedDofIndices ()
 Gets the indices for all dofs conected to the constraint Get indices for all the slave Jacobian columns, not only those based on the mesh connectivity. More...
 

Detailed Description

Definition at line 20 of file SparsityBasedContactConstraint.h.

Constructor & Destructor Documentation

SparsityBasedContactConstraint::SparsityBasedContactConstraint ( const InputParameters &  parameters)
inline

Definition at line 23 of file SparsityBasedContactConstraint.h.

24  : NodeFaceConstraint(parameters){};
virtual SparsityBasedContactConstraint::~SparsityBasedContactConstraint ( )
inlinevirtual

Definition at line 25 of file SparsityBasedContactConstraint.h.

25 {}

Member Function Documentation

virtual Real SparsityBasedContactConstraint::computeQpJacobian ( Moose::ConstraintJacobianType  )
inlinevirtual

Reimplemented in GluedContactConstraint.

Definition at line 39 of file SparsityBasedContactConstraint.h.

40  {
41  mooseError("Unimplemented pure virtual method. SparsityBasedContactConstraint should only be "
42  "used as a base class");
43  return 0;
44  }
virtual Real SparsityBasedContactConstraint::computeQpResidual ( Moose::ConstraintType  )
inlinevirtual

Reimplemented in GluedContactConstraint.

Definition at line 33 of file SparsityBasedContactConstraint.h.

34  {
35  mooseError("Unimplemented pure virtual method. SparsityBasedContactConstraint should only be "
36  "used as a base class");
37  return 0;
38  }
virtual Real SparsityBasedContactConstraint::computeQpSlaveValue ( )
inlinevirtual

Reimplemented in GluedContactConstraint.

Definition at line 27 of file SparsityBasedContactConstraint.h.

28  {
29  mooseError("Unimplemented pure virtual method. SparsityBasedContactConstraint should only be "
30  "used as a base class");
31  return 0;
32  }
void SparsityBasedContactConstraint::getConnectedDofIndices ( )
virtual

Gets the indices for all dofs conected to the constraint Get indices for all the slave Jacobian columns, not only those based on the mesh connectivity.

Definition at line 25 of file SparsityBasedContactConstraint.C.

Referenced by computeQpJacobian().

26 {
27 #if defined(LIBMESH_HAVE_PETSC) && !PETSC_VERSION_LESS_THAN(3, 3, 0)
28  _connected_dof_indices.clear();
29 
30  // An ugly hack: have to extract the row and look at it's sparsity structure, since otherwise I
31  // won't get the dofs connected to this row by virtue of intervariable coupling
32  // This is easier than sifting through all coupled variables, selecting those active on the
33  // current subdomain, dealing with the scalar variables, etc.
34  // Also, importantly, this will miss coupling to variables that might have introduced by prior
35  // constraints similar to this one!
36  PetscMatrix<Number> * petsc_jacobian = dynamic_cast<PetscMatrix<Number> *>(_jacobian);
37  mooseAssert(petsc_jacobian, "Expected a PETSc matrix");
38  Mat jac = petsc_jacobian->mat();
39  PetscErrorCode ierr;
40  PetscInt ncols;
41  const PetscInt * cols;
42  ierr = MatGetRow(jac, static_cast<PetscInt>(_var.nodalDofIndex()), &ncols, &cols, PETSC_NULL);
43  CHKERRABORT(_communicator.get(), ierr);
44  bool debug = false;
45  if (debug)
46  {
47  libMesh::out << "_connected_dof_indices: adding " << ncols << " dofs from Jacobian row["
48  << _var.nodalDofIndex() << "] = [";
49  }
50  for (PetscInt i = 0; i < ncols; ++i)
51  {
52  if (debug)
53  {
54  libMesh::out << cols[i] << " ";
55  }
56  _connected_dof_indices.push_back(cols[i]);
57  }
58  if (debug)
59  {
60  libMesh::out << "]\n";
61  }
62  ierr = MatRestoreRow(jac, static_cast<PetscInt>(_var.nodalDofIndex()), &ncols, &cols, PETSC_NULL);
63  CHKERRABORT(_communicator.get(), ierr);
64 #else
65  NodeFaceConstraint::getConnectedDofIndices();
66 #endif
67 }

The documentation for this class was generated from the following files: