www.mooseframework.org
SparsityBasedContactConstraint.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 
9 
10 // MOOSE includes
11 #include "MooseVariable.h"
12 
13 #include "libmesh/petsc_macro.h"
14 #include "libmesh/petsc_matrix.h"
15 
16 template <>
17 InputParameters
19 {
20  InputParameters params = validParams<NodeFaceConstraint>();
21  return params;
22 }
23 
24 void
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 }
virtual void getConnectedDofIndices()
Gets the indices for all dofs conected to the constraint Get indices for all the slave Jacobian colum...
InputParameters validParams< SparsityBasedContactConstraint >()