Classes
struct	PeriodicQuadFunction

class	PeriodicBCTest

Functions
Number	quadratic_solution (const Point &p, const Parameters &, const std::string &, const std::string &)

void	periodic_bc_test_poisson (EquationSystems &es, const std::string &)

	CPPUNIT_TEST_SUITE_REGISTRATION (PeriodicBCTest)

Variables
const boundary_id_type	top_id = 50

const boundary_id_type	bottom_id = 60

const boundary_id_type	side_id = 70

Function Documentation

◆ CPPUNIT_TEST_SUITE_REGISTRATION()

CPPUNIT_TEST_SUITE_REGISTRATION ( PeriodicBCTest )

◆ periodic_bc_test_poisson()

void periodic_bc_test_poisson	(	EquationSystems &	es,
		const std::string &
	)

Definition at line 74 of file periodic_bc_test.C.

{
  const MeshBase& mesh = es.get_mesh();
  LinearImplicitSystem& system = es.get_system<LinearImplicitSystem>("PBCSys");
  const DofMap& dof_map = system.get_dof_map();
 
  // Interior FE to integrate smooth part of poisson problem
  FEType fe_type = dof_map.variable_type(0);
  std::unique_ptr<FEBase> fe (FEBase::build(2, fe_type));
  QGauss qrule (2, fe_type.default_quadrature_order());
  fe->attach_quadrature_rule(&qrule);
 
  const std::vector<Real> & JxW = fe->get_JxW();
  const std::vector<std::vector<Real>> & phi = fe->get_phi();
  const std::vector<std::vector<RealGradient>> & dphi = fe->get_dphi();
 
  // Face FE to integrate delta function part of forcing
  std::unique_ptr<FEBase> fe_face (FEBase::build(2, fe_type));
  QGauss qface (1, fe_type.default_quadrature_order());
  fe_face->attach_quadrature_rule(&qface);
 
  const std::vector<Real> & JxW_face = fe_face->get_JxW();
  const std::vector<std::vector<Real>> & phi_face = fe_face->get_phi();
  std::unique_ptr<const Elem> elem_side;
 
  DenseMatrix<Number> Ke;
  DenseVector<Number> Fe;
 
  std::vector<dof_id_type> dof_indices;
 
  for (const Elem * elem : mesh.active_local_element_ptr_range())
    {
      dof_map.dof_indices (elem, dof_indices);
      const unsigned int n_dofs = dof_indices.size();
 
      Ke.resize (n_dofs, n_dofs);
      Fe.resize (n_dofs);
 
      fe->reinit (elem);
 
      // Integrate the poisson problem over the element interior, where we have
      // a nice f=2 forcing function aside from the discontinuity
      for (unsigned int qp=0; qp<qrule.n_points(); qp++)
        {
          for (unsigned int i=0; i != n_dofs; i++)
            {
              for (unsigned int j=0; j != n_dofs; j++)
                Ke(i,j) += -JxW[qp]*(dphi[i][qp]*dphi[j][qp]);
 
              Fe(i) += JxW[qp]*phi[i][qp]*2;
            }
        }
 
      for(unsigned int s=0; s != elem->n_sides(); ++s)
      {
        elem->build_side_ptr(elem_side, s);
        Point centroid = elem_side->vertex_average();
        if (std::abs(centroid(1) - 1) < TOLERANCE)
          {
            fe_face->reinit(elem, s);
 
            for (unsigned int qp=0; qp<qface.n_points(); qp++)
              {
                // delta function divided by 2 since we hit it from
                // both sides
                for (unsigned int i=0; i != n_dofs; i++)
                  Fe(i) += JxW_face[qp]*phi_face[i][qp]*(-4);
              }
          }
      }
 
      dof_map.heterogenously_constrain_element_matrix_and_vector(Ke, Fe, dof_indices);
      system.matrix->add_matrix (Ke, dof_indices);
      system.rhs->add_vector    (Fe, dof_indices);
    }
 
  system.rhs->close();
  system.matrix->close();
}

References libMesh::SparseMatrix< T >::add_matrix(), libMesh::NumericVector< T >::add_vector(), libMesh::FEGenericBase< OutputType >::build(), libMesh::NumericVector< T >::close(), libMesh::SparseMatrix< T >::close(), libMesh::FEType::default_quadrature_order(), libMesh::System::get_dof_map(), libMesh::EquationSystems::get_mesh(), libMesh::EquationSystems::get_system(), libMesh::ImplicitSystem::matrix, mesh, libMesh::QBase::n_points(), libMesh::DenseVector< T >::resize(), libMesh::DenseMatrix< T >::resize(), libMesh::ExplicitSystem::rhs, libMesh::DenseVector< T >::size(), and libMesh::TOLERANCE.

Referenced by PeriodicBCTest::testPeriodicBC().

◆ quadratic_solution()

Number quadratic_solution	(	const Point &	p,
		const Parameters &	,
		const std::string &	,
		const std::string &
	)

Definition at line 23 of file periodic_bc_test.C.

{
  const Real & y = LIBMESH_DIM > 1 ? p(1) : 0;
 
  // Discontinuity at y=1 from the forcing function there.
  // Periodic on -3 < y < 2
  return (y > 1) ? (y-3)*(y-3) : (y+1)*(y+1);
}