libMesh
Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | List of all members
OverlappingCouplingGhostingTest Class Reference
Inheritance diagram for OverlappingCouplingGhostingTest:
[legend]

Public Member Functions

 CPPUNIT_TEST_SUITE (OverlappingCouplingGhostingTest)
 
 CPPUNIT_TEST (testSparsityCouplingMatrix)
 
 CPPUNIT_TEST (testSparsityNullCouplingMatrix)
 
 CPPUNIT_TEST (testSparsityNullCouplingMatrixUnifRef)
 
 CPPUNIT_TEST_SUITE_END ()
 
void setUp ()
 
void tearDown ()
 
void testSparsityCouplingMatrix ()
 
void testSparsityNullCouplingMatrix ()
 
void testSparsityNullCouplingMatrixUnifRef ()
 

Protected Member Functions

void build_quad_mesh (unsigned int n_refinements=0)
 
void init (MeshBase &mesh)
 
void clear ()
 
void setup_coupling_matrix (std::unique_ptr< CouplingMatrix > &coupling)
 

Protected Attributes

std::unique_ptr< MeshBase_mesh
 
std::unique_ptr< EquationSystems_es
 

Private Member Functions

void run_sparsity_pattern_test (const unsigned int n_refinements, bool build_coupling_matrix)
 

Detailed Description

Definition at line 621 of file overlapping_coupling_test.C.

Member Function Documentation

◆ build_quad_mesh()

void OverlappingTestBase::build_quad_mesh ( unsigned int  n_refinements = 0)
inlineprotectedinherited

Definition at line 238 of file overlapping_coupling_test.C.

239  {
240  // We are making assumptions in various places about the presence
241  // of the elements on the current processor so we're restricting to
242  // ReplicatedMesh for now.
243  _mesh = libmesh_make_unique<ReplicatedMesh>(*TestCommWorld);
244 
245  _mesh->set_mesh_dimension(2);
246 
247  _mesh->add_point( Point(0.0,0.0),0 );
248  _mesh->add_point( Point(1.0,0.0),1 );
249  _mesh->add_point( Point(1.0,1.0),2 );
250  _mesh->add_point( Point(0.0,1.0),3 );
251 
252  {
253  Elem* elem = _mesh->add_elem( new Quad4 );
254  elem->set_id(0);
255  elem->subdomain_id() = 1;
256 
257  for (unsigned int n=0; n<4; n++)
258  elem->set_node(n) = _mesh->node_ptr(n);
259  }
260 
261  _mesh->add_point( Point(1.0,2.0),4 );
262  _mesh->add_point( Point(0.0,2.0),5 );
263 
264  {
265  Elem* elem = _mesh->add_elem( new Quad4 );
266  elem->set_id(1);
267  elem->subdomain_id() = 1;
268 
269  elem->set_node(0) = _mesh->node_ptr(3);
270  elem->set_node(1) = _mesh->node_ptr(2);
271  elem->set_node(2) = _mesh->node_ptr(4);
272  elem->set_node(3) = _mesh->node_ptr(5);
273  }
274 
275  _mesh->add_point( Point(0.0,0.0),6 );
276  _mesh->add_point( Point(1.0,0.0),7 );
277  _mesh->add_point( Point(1.0,2.0),8 );
278  _mesh->add_point( Point(0.0,2.0),9 );
279 
280  {
281  Elem* elem = _mesh->add_elem( new Quad4 );
282  elem->set_id(2);
283  elem->subdomain_id() = 2;
284 
285  elem->set_node(0) = _mesh->node_ptr(6);
286  elem->set_node(1) = _mesh->node_ptr(7);
287  elem->set_node(2) = _mesh->node_ptr(8);
288  elem->set_node(3) = _mesh->node_ptr(9);
289  }
290 
291  _mesh->partitioner() = libmesh_make_unique<OverlappingTestPartitioner>();
292 
293  _mesh->prepare_for_use();
294 
295 #ifdef LIBMESH_ENABLE_AMR
296  if (n_refinements > 0)
297  {
298  MeshRefinement refine(*_mesh);
299  refine.uniformly_refine(n_refinements);
300  }
301 #endif // LIBMESH_ENABLE_AMR
302  }

References libMesh::DofObject::set_id(), libMesh::Elem::set_node(), libMesh::Elem::subdomain_id(), TestCommWorld, and libMesh::MeshRefinement::uniformly_refine().

◆ clear()

void OverlappingTestBase::clear ( )
inlineprotectedinherited

Definition at line 324 of file overlapping_coupling_test.C.

325  {
326  _es.reset();
327  _mesh.reset();
328  }

◆ CPPUNIT_TEST() [1/3]

OverlappingCouplingGhostingTest::CPPUNIT_TEST ( testSparsityCouplingMatrix  )

◆ CPPUNIT_TEST() [2/3]

OverlappingCouplingGhostingTest::CPPUNIT_TEST ( testSparsityNullCouplingMatrix  )

◆ CPPUNIT_TEST() [3/3]

OverlappingCouplingGhostingTest::CPPUNIT_TEST ( testSparsityNullCouplingMatrixUnifRef  )

◆ CPPUNIT_TEST_SUITE()

OverlappingCouplingGhostingTest::CPPUNIT_TEST_SUITE ( OverlappingCouplingGhostingTest  )

◆ CPPUNIT_TEST_SUITE_END()

OverlappingCouplingGhostingTest::CPPUNIT_TEST_SUITE_END ( )

◆ init()

void OverlappingTestBase::init ( MeshBase mesh)
inlineprotectedinherited

Definition at line 304 of file overlapping_coupling_test.C.

305  {
306  _es = libmesh_make_unique<EquationSystems>(*_mesh);
307  LinearImplicitSystem & sys = _es->add_system<LinearImplicitSystem> ("SimpleSystem");
308 
309  std::set<subdomain_id_type> sub_one;
310  sub_one.insert(1);
311 
312  std::set<subdomain_id_type> sub_two;
313  sub_two.insert(2);
314 
315  sys.add_variable("U", FIRST, LAGRANGE, &sub_two);
316  sys.add_variable("L", FIRST, LAGRANGE, &sub_two);
317 
318  sys.add_variable("V", FIRST, LAGRANGE, &sub_one);
319  sys.add_variable("p", FIRST, LAGRANGE, &sub_one);
320 
321  _es->init();
322  }

References libMesh::System::add_variable(), libMesh::FIRST, and libMesh::LAGRANGE.

◆ run_sparsity_pattern_test()

void OverlappingCouplingGhostingTest::run_sparsity_pattern_test ( const unsigned int  n_refinements,
bool  build_coupling_matrix 
)
inlineprivate

Definition at line 658 of file overlapping_coupling_test.C.

659  {
660  this->build_quad_mesh(n_refinements);
661  this->init(*_mesh);
662 
663  std::unique_ptr<CouplingMatrix> coupling_matrix;
664  if (build_coupling_matrix)
665  this->setup_coupling_matrix(coupling_matrix);
666 
667  LinearImplicitSystem & system = _es->get_system<LinearImplicitSystem>("SimpleSystem");
668 
669  // If we don't add this coupling functor and properly recompute the
670  // sparsity pattern, then PETSc will throw a malloc error when we
671  // try to assemble into the global matrix
672  OverlappingCouplingFunctor coupling_functor(system);
673  coupling_functor.set_coupling_matrix(coupling_matrix);
674 
675  DofMap & dof_map = system.get_dof_map();
676  dof_map.add_coupling_functor(coupling_functor);
677  dof_map.reinit_send_list(system.get_mesh());
678 
679  // Update current local solution
681 
682  system.current_local_solution->init(system.n_dofs(), system.n_local_dofs(),
683  dof_map.get_send_list(), false,
685 
686  system.solution->localize(*(system.current_local_solution),dof_map.get_send_list());
687 
688  // Now that we've added the coupling functor, we need
689  // to recompute the sparsity
690  dof_map.clear_sparsity();
691  dof_map.compute_sparsity(system.get_mesh());
692 
693  // Now that we've recomputed the sparsity pattern, we need
694  // to reinitialize the system matrix.
695  libMesh::SparseMatrix<libMesh::Number> & matrix = system.get_matrix("System Matrix");
696  libmesh_assert(dof_map.is_attached(matrix));
697  matrix.init();
698 
699  std::unique_ptr<PointLocatorBase> point_locator = _mesh->sub_point_locator();
700 
701  const unsigned int u_var = system.variable_number("U");
702  const unsigned int v_var = system.variable_number("V");
703 
704  DenseMatrix<Number> K12, K21;
705 
706  FEMContext subdomain_one_context(system);
707  FEMContext subdomain_two_context(system);
708 
709  // Add normally coupled parts of the matrix
710  for (const auto & elem : _mesh->active_local_subdomain_elements_ptr_range(1))
711  {
712  subdomain_one_context.pre_fe_reinit(system,elem);
713  subdomain_one_context.elem_fe_reinit();
714 
715  std::vector<dof_id_type> & rows = subdomain_one_context.get_dof_indices();
716 
717  // Fill with ones in case PETSc ignores the zeros at some point
718  std::fill( subdomain_one_context.get_elem_jacobian().get_values().begin(),
719  subdomain_one_context.get_elem_jacobian().get_values().end(),
720  1);
721 
722  // Insert the Jacobian for the dofs for this element
723  system.matrix->add_matrix( subdomain_one_context.get_elem_jacobian(), rows );
724  }
725 
726  for (const auto & elem : _mesh->active_local_subdomain_elements_ptr_range(2))
727  {
728  // A little extra unit testing on the range iterator
729  CPPUNIT_ASSERT_EQUAL(2, (int)elem->subdomain_id());
730 
731  const std::vector<libMesh::Point> & qpoints = subdomain_two_context.get_element_fe(u_var)->get_xyz();
732 
733  // Setup the context for the current element
734  subdomain_two_context.pre_fe_reinit(system,elem);
735  subdomain_two_context.elem_fe_reinit();
736 
737  // We're only assembling rows for the dofs on subdomain 2 (U,L), so
738  // the current element will have all those dof_indices.
739  std::vector<dof_id_type> & rows = subdomain_two_context.get_dof_indices();
740 
741  std::fill( subdomain_two_context.get_elem_jacobian().get_values().begin(),
742  subdomain_two_context.get_elem_jacobian().get_values().end(),
743  1);
744 
745  // Insert the Jacobian for the normally coupled dofs for this element
746  system.matrix->add_matrix( subdomain_two_context.get_elem_jacobian(), rows );
747 
748  std::set<subdomain_id_type> allowed_subdomains;
749  allowed_subdomains.insert(1);
750 
751  // Now loop over the quadrature points and find the subdomain-one element that overlaps
752  // with the current subdomain-two element and then add a local element matrix with
753  // the coupling to the global matrix to try and trip any issues with sparsity pattern
754  // construction
755  for ( const auto & qp : qpoints )
756  {
757  const Elem * overlapping_elem = (*point_locator)( qp, &allowed_subdomains );
758  CPPUNIT_ASSERT(overlapping_elem);
759 
760  // Setup the context for the overlapping element
761  subdomain_one_context.pre_fe_reinit(system,overlapping_elem);
762  subdomain_one_context.elem_fe_reinit();
763 
764  // We're only coupling to the "V" variable so only need those dof indices
765  std::vector<dof_id_type> & v_indices = subdomain_one_context.get_dof_indices(v_var);
766  std::vector<dof_id_type> columns(rows);
767  columns.insert( columns.end(), v_indices.begin(), v_indices.end() );
768 
769  // This will also zero the matrix so we can just insert zeros for this test
770  K21.resize( rows.size(), columns.size() );
771 
772  std::fill(K21.get_values().begin(), K21.get_values().end(), 1);
773 
774  // Now adding this local matrix to the global would trip a PETSc
775  // malloc error if the sparsity pattern hasn't been correctly
776  // built to include the overlapping coupling.
777  system.matrix->add_matrix (K21, rows, columns);
778 
779  // Now add the other part of the overlapping coupling
780  K12.resize(v_indices.size(), rows.size());
781  std::fill(K12.get_values().begin(), K12.get_values().end(), 1);
782  system.matrix->add_matrix(K12,v_indices,rows);
783  }
784  } // end element loop
785 
786  // We need to make sure to close the matrix for this test. There could still
787  // be PETSc malloc errors tripped here if we didn't allocate the off-processor
788  // part of the sparsity pattern correctly.
789  system.matrix->close();
790  }

References libMesh::DofMap::add_coupling_functor(), libMesh::SparseMatrix< T >::add_matrix(), libMesh::NumericVector< T >::build(), libMesh::DofMap::clear_sparsity(), libMesh::SparseMatrix< T >::close(), libMesh::ParallelObject::comm(), libMesh::DofMap::compute_sparsity(), libMesh::System::current_local_solution, libMesh::FEMContext::elem_fe_reinit(), libMesh::DiffContext::get_dof_indices(), libMesh::System::get_dof_map(), libMesh::DiffContext::get_elem_jacobian(), libMesh::FEMContext::get_element_fe(), libMesh::ImplicitSystem::get_matrix(), libMesh::System::get_mesh(), libMesh::DofMap::get_send_list(), libMesh::DenseMatrix< T >::get_values(), libMesh::GHOSTED, libMesh::TriangleWrapper::init(), libMesh::SparseMatrix< T >::init(), libMesh::DofMap::is_attached(), libMesh::libmesh_assert(), libMesh::ImplicitSystem::matrix, libMesh::System::n_dofs(), libMesh::System::n_local_dofs(), libMesh::FEMContext::pre_fe_reinit(), libMesh::DofMap::reinit_send_list(), libMesh::DenseMatrix< T >::resize(), libMesh::System::solution, and libMesh::System::variable_number().

◆ setUp()

void OverlappingCouplingGhostingTest::setUp ( )
inline

Definition at line 635 of file overlapping_coupling_test.C.

636  {}

◆ setup_coupling_matrix()

void OverlappingTestBase::setup_coupling_matrix ( std::unique_ptr< CouplingMatrix > &  coupling)
inlineprotectedinherited

Definition at line 330 of file overlapping_coupling_test.C.

331  {
332  System & system = _es->get_system("SimpleSystem");
333 
334  coupling = libmesh_make_unique<CouplingMatrix>(system.n_vars());
335 
336  const unsigned int u_var = system.variable_number("U");
337  const unsigned int l_var = system.variable_number("L");
338  const unsigned int v_var = system.variable_number("V");
339  const unsigned int p_var = system.variable_number("p");
340 
341  // Only adding the overlapping couplings since the primary ones should
342  // be there by default.
343  (*coupling)(u_var,v_var) = true;
344  (*coupling)(l_var,v_var) = true;
345  (*coupling)(l_var,p_var) = true;
346  (*coupling)(v_var,u_var) = true;
347  (*coupling)(v_var,l_var) = true;
348  }

References libMesh::System::n_vars(), and libMesh::System::variable_number().

◆ tearDown()

void OverlappingCouplingGhostingTest::tearDown ( )
inline

Definition at line 638 of file overlapping_coupling_test.C.

639  { this->clear(); }

◆ testSparsityCouplingMatrix()

void OverlappingCouplingGhostingTest::testSparsityCouplingMatrix ( )
inline

Definition at line 641 of file overlapping_coupling_test.C.

642  {
643  this->run_sparsity_pattern_test(0, true);
644  }

◆ testSparsityNullCouplingMatrix()

void OverlappingCouplingGhostingTest::testSparsityNullCouplingMatrix ( )
inline

Definition at line 646 of file overlapping_coupling_test.C.

647  {
648  this->run_sparsity_pattern_test(0, false);
649  }

◆ testSparsityNullCouplingMatrixUnifRef()

void OverlappingCouplingGhostingTest::testSparsityNullCouplingMatrixUnifRef ( )
inline

Definition at line 651 of file overlapping_coupling_test.C.

652  {
653  this->run_sparsity_pattern_test(1, false);
654  }

Member Data Documentation

◆ _es

std::unique_ptr<EquationSystems> OverlappingTestBase::_es
protectedinherited

Definition at line 236 of file overlapping_coupling_test.C.

◆ _mesh

std::unique_ptr<MeshBase> OverlappingTestBase::_mesh
protectedinherited

Definition at line 234 of file overlapping_coupling_test.C.


The documentation for this class was generated from the following file:
libMesh::SparseMatrix::close
virtual void close()=0
Calls the SparseMatrix's internal assembly routines, ensuring that the values are consistent across p...
libMesh::System
Manages consistently variables, degrees of freedom, and coefficient vectors.
Definition: system.h:100
libMesh::System::n_vars
unsigned int n_vars() const
Definition: system.h:2155
libMesh::DofMap::compute_sparsity
void compute_sparsity(const MeshBase &)
Computes the sparsity pattern for the matrices corresponding to proc_id and sends that data to Linear...
Definition: dof_map.C:1768
libMesh::DofMap::is_attached
bool is_attached(SparseMatrix< Number > &matrix)
Matrices should not be attached more than once.
Definition: dof_map.C:310
OverlappingCouplingGhostingTest::run_sparsity_pattern_test
void run_sparsity_pattern_test(const unsigned int n_refinements, bool build_coupling_matrix)
Definition: overlapping_coupling_test.C:658
libMesh::SparseMatrix::init
virtual void init(const numeric_index_type m, const numeric_index_type n, const numeric_index_type m_l, const numeric_index_type n_l, const numeric_index_type nnz=30, const numeric_index_type noz=10, const numeric_index_type blocksize=1)=0
Initialize SparseMatrix with the specified sizes.
libMesh::DofObject::set_id
dof_id_type & set_id()
Definition: dof_object.h:776
libMesh::ParallelObject::comm
const Parallel::Communicator & comm() const
Definition: parallel_object.h:94
libMesh::DenseMatrix< Number >
libMesh::GHOSTED
Definition: enum_parallel_type.h:37
libMesh::DenseMatrix::get_values
std::vector< T > & get_values()
Definition: dense_matrix.h:371
libMesh::SparseMatrix
Generic sparse matrix.
Definition: vector_fe_ex5.C:45
libMesh::NumericVector::build
static std::unique_ptr< NumericVector< T > > build(const Parallel::Communicator &comm, const SolverPackage solver_package=libMesh::default_solver_package())
Builds a NumericVector on the processors in communicator comm using the linear solver package specifi...
Definition: numeric_vector.C:49
libMesh::DenseMatrix::resize
void resize(const unsigned int new_m, const unsigned int new_n)
Resize the matrix.
Definition: dense_matrix.h:822
libMesh::MeshRefinement
Implements (adaptive) mesh refinement algorithms for a MeshBase.
Definition: mesh_refinement.h:61
OverlappingCouplingFunctor
Definition: overlapping_coupling_test.C:37
libMesh::libmesh_assert
libmesh_assert(ctx)
libMesh::System::n_local_dofs
dof_id_type n_local_dofs() const
Definition: system.C:187
OverlappingTestBase::_mesh
std::unique_ptr< MeshBase > _mesh
Definition: overlapping_coupling_test.C:234
libMesh::System::get_mesh
const MeshBase & get_mesh() const
Definition: system.h:2083
libMesh::System::add_variable
unsigned int add_variable(const std::string &var, const FEType &type, const std::set< subdomain_id_type > *const active_subdomains=nullptr)
Adds the variable var to the list of variables for this system.
Definition: system.C:1069
libMesh::System::current_local_solution
std::unique_ptr< NumericVector< Number > > current_local_solution
All the values I need to compute my contribution to the simulation at hand.
Definition: system.h:1551
OverlappingTestBase::_es
std::unique_ptr< EquationSystems > _es
Definition: overlapping_coupling_test.C:236
TestCommWorld
libMesh::Parallel::Communicator * TestCommWorld
Definition: driver.C:111
libMesh::Point
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
libMesh::ImplicitSystem::matrix
SparseMatrix< Number > * matrix
The system matrix.
Definition: implicit_system.h:393
libMesh::Quad4
The QUAD4 is an element in 2D composed of 4 nodes.
Definition: face_quad4.h:51
libMesh::DofMap::get_send_list
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:496
libMesh::Elem::set_node
virtual Node *& set_node(const unsigned int i)
Definition: elem.h:2059
libMesh::DofMap::reinit_send_list
void reinit_send_list(MeshBase &mesh)
Clears the _send_list vector and then rebuilds it.
Definition: dof_map.C:1690
OverlappingTestBase::setup_coupling_matrix
void setup_coupling_matrix(std::unique_ptr< CouplingMatrix > &coupling)
Definition: overlapping_coupling_test.C:330
OverlappingTestBase::build_quad_mesh
void build_quad_mesh(unsigned int n_refinements=0)
Definition: overlapping_coupling_test.C:238
libMesh::System::solution
std::unique_ptr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1539
libMesh::SparseMatrix::add_matrix
virtual void add_matrix(const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols)=0
Add the full matrix dm to the SparseMatrix.
libMesh::DofMap
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:176
OverlappingTestBase::init
void init(MeshBase &mesh)
Definition: overlapping_coupling_test.C:304
libMesh::DofMap::clear_sparsity
void clear_sparsity()
Clears the sparsity pattern.
Definition: dof_map.C:1800
libMesh::Elem::subdomain_id
subdomain_id_type subdomain_id() const
Definition: elem.h:2069
libMesh::Elem
This is the base class from which all geometric element types are derived.
Definition: elem.h:100
libMesh::System::variable_number
unsigned short int variable_number(const std::string &var) const
Definition: system.C:1232
libMesh::System::get_dof_map
const DofMap & get_dof_map() const
Definition: system.h:2099
libMesh::System::n_dofs
dof_id_type n_dofs() const
Definition: system.C:150
libMesh::LAGRANGE
Definition: enum_fe_family.h:36
libMesh::LinearImplicitSystem
Manages consistently variables, degrees of freedom, coefficient vectors, matrices and linear solvers ...
Definition: linear_implicit_system.h:55
libMesh::ImplicitSystem::get_matrix
const SparseMatrix< Number > & get_matrix(const std::string &mat_name) const
Definition: implicit_system.C:262
OverlappingTestBase::clear
void clear()
Definition: overlapping_coupling_test.C:324
libMesh::FIRST
Definition: enum_order.h:42
libMesh::DofMap::add_coupling_functor
void add_coupling_functor(GhostingFunctor &coupling_functor, bool to_mesh=true)
Adds a functor which can specify coupling requirements for creation of sparse matrices.
Definition: dof_map.C:1838
libMesh::FEMContext
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:62