libMesh::HPCoarsenTest Class Reference

This class uses the error estimate given by different types of derefinement (h coarsening or p reduction) to choose between h refining and p elevation. More...

#include <hp_coarsentest.h>

Inheritance diagram for libMesh::HPCoarsenTest:

Public Member Functions
	HPCoarsenTest ()
	Constructor. More...
	HPCoarsenTest (const HPCoarsenTest &)=delete
	This class cannot be (default) copy constructed/assigned because it has unique_ptr members. More...
HPCoarsenTest &	operator= (const HPCoarsenTest &)=delete
	HPCoarsenTest (HPCoarsenTest &&)=default
	Defaulted move ctor, move assignment operator, and destructor. More...
HPCoarsenTest &	operator= (HPCoarsenTest &&)=default
virtual	~HPCoarsenTest ()=default
virtual void	select_refinement (System &system) override
	This pure virtual function must be redefined in derived classes to take a mesh flagged for h refinement and potentially change the desired refinement type. More...
void	extra_quadrature_order (const int extraorder)
	Increases or decreases the order of the quadrature rule used for numerical integration. More...

Public Attributes
Real	p_weight
	Because the coarsening test seems to always choose p refinement, we're providing an option to make h refinement more likely. More...
std::vector< float >	component_scale
	This vector can be used to "scale" certain variables in a system. More...

Protected Member Functions
void	add_projection (const System &, const Elem *, unsigned int var)
	The helper function which adds individual fine element data to the coarse element projection. More...

Protected Attributes
Elem *	coarse
	The coarse element on which a solution projection is cached. More...
std::vector< dof_id_type >	dof_indices
	Global DOF indices for fine elements. More...
std::unique_ptr< FEBase >	fe
	The finite element objects for fine and coarse elements. More...
std::unique_ptr< FEBase >	fe_coarse
const std::vector< std::vector< Real > > *	phi
	The shape functions and their derivatives. More...
const std::vector< std::vector< Real > > *	phi_coarse
const std::vector< std::vector< RealGradient > > *	dphi
const std::vector< std::vector< RealGradient > > *	dphi_coarse
const std::vector< std::vector< RealTensor > > *	d2phi
const std::vector< std::vector< RealTensor > > *	d2phi_coarse
const std::vector< Real > *	JxW
	Mapping jacobians. More...
const std::vector< Point > *	xyz_values
	Quadrature locations. More...
std::vector< Point >	coarse_qpoints
std::unique_ptr< QBase >	qrule
	The quadrature rule for the fine element. More...
DenseMatrix< Number >	Ke
	Linear system for projections. More...
DenseVector< Number >	Fe
DenseVector< Number >	Uc
	Coefficients for projected coarse and projected p-derefined solutions. More...
DenseVector< Number >	Up
int	_extra_order
	Extra order to use for quadrature rule. More...

Detailed Description

This class uses the error estimate given by different types of derefinement (h coarsening or p reduction) to choose between h refining and p elevation.

Currently we assume that a set of elements has already been flagged for h refinement, and we may want to change some of those elements to be flagged for p refinement.

This code is currently experimental and will not produce optimal hp meshes without significant improvement.

Author

Roy H. Stogner

Date

2006

Definition at line 67 of file hp_coarsentest.h.

Constructor & Destructor Documentation

HPCoarsenTest() [1/3]

libMesh::HPCoarsenTest::HPCoarsenTest ( )

inline

Constructor.

Definition at line 74 of file hp_coarsentest.h.

                  : p_weight(1.0), _extra_order(1)
  {
    libmesh_experimental();
  }

HPCoarsenTest() [2/3]

libMesh::HPCoarsenTest::HPCoarsenTest ( const HPCoarsenTest &  )

delete

This class cannot be (default) copy constructed/assigned because it has unique_ptr members.

Explicitly deleting these functions is the best way to document this fact.

HPCoarsenTest() [3/3]

libMesh::HPCoarsenTest::HPCoarsenTest ( HPCoarsenTest &&  )

default

Defaulted move ctor, move assignment operator, and destructor.

~HPCoarsenTest()

virtual libMesh::HPCoarsenTest::~HPCoarsenTest ( )

virtualdefault

Member Function Documentation

add_projection()

void libMesh::HPCoarsenTest::add_projection ( const System &  system, const Elem *  elem, unsigned int  var  )

protected

The helper function which adds individual fine element data to the coarse element projection.

Definition at line 48 of file hp_coarsentest.C.

References libMesh::Elem::active(), libMesh::TypeVector< T >::add_scaled(), libMesh::TypeTensor< T >::add_scaled(), libMesh::C_ONE, libMesh::C_ZERO, libMesh::Elem::child_ref_range(), coarse, coarse_qpoints, libMesh::TypeTensor< T >::contract(), libMesh::System::current_solution(), d2phi, d2phi_coarse, dof_indices, libMesh::DofMap::dof_indices(), dphi, fe, Fe, fe_coarse, libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::index_range(), libMesh::FEMap::inverse_map(), Ke, libMesh::libmesh_assert(), libMesh::make_range(), libMesh::MeshBase::mesh_dimension(), phi_coarse, qrule, libMesh::DenseVector< T >::resize(), libMesh::DenseMatrix< T >::resize(), libMesh::DenseVector< T >::size(), libMesh::Elem::subactive(), Uc, xyz_values, libMesh::DenseMatrix< T >::zero(), libMesh::DenseVector< T >::zero(), and libMesh::zero.

Referenced by select_refinement().

{
  // If we have children, we need to add their projections instead
  if (!elem->active())
    {
      libmesh_assert(!elem->subactive());
      for (auto & child : elem->child_ref_range())
        this->add_projection(system, &child, var);
      return;
    }

  // The DofMap for this system
  const DofMap & dof_map = system.get_dof_map();

  const FEContinuity cont = fe->get_continuity();

  fe->reinit(elem);

  dof_map.dof_indices(elem, dof_indices, var);

  const unsigned int n_dofs =
    cast_int<unsigned int>(dof_indices.size());

  FEMap::inverse_map (system.get_mesh().mesh_dimension(), coarse,
                      *xyz_values, coarse_qpoints);

  fe_coarse->reinit(coarse, &coarse_qpoints);

  const unsigned int n_coarse_dofs =
    cast_int<unsigned int>(phi_coarse->size());

  if (Uc.size() == 0)
    {
      Ke.resize(n_coarse_dofs, n_coarse_dofs);
      Ke.zero();
      Fe.resize(n_coarse_dofs);
      Fe.zero();
      Uc.resize(n_coarse_dofs);
      Uc.zero();
    }
  libmesh_assert_equal_to (Uc.size(), phi_coarse->size());

  // Loop over the quadrature points
  for (auto qp : make_range(qrule->n_points()))
    {
      // The solution value at the quadrature point
      Number val = libMesh::zero;
      Gradient grad;
      Tensor hess;

      for (unsigned int i=0; i != n_dofs; i++)
        {
          dof_id_type dof_num = dof_indices[i];
          val += (*phi)[i][qp] *
            system.current_solution(dof_num);
          if (cont == C_ZERO || cont == C_ONE)
            grad.add_scaled((*dphi)[i][qp],system.current_solution(dof_num));
          if (cont == C_ONE)
            hess.add_scaled((*d2phi)[i][qp], system.current_solution(dof_num));
        }

      // The projection matrix and vector
      for (auto i : index_range(Fe))
        {
          Fe(i) += (*JxW)[qp] *
            (*phi_coarse)[i][qp]*val;
          if (cont == C_ZERO || cont == C_ONE)
            Fe(i) += (*JxW)[qp] *
              (grad*(*dphi_coarse)[i][qp]);
          if (cont == C_ONE)
            Fe(i) += (*JxW)[qp] *
              hess.contract((*d2phi_coarse)[i][qp]);

          for (auto j : index_range(Fe))
            {
              Ke(i,j) += (*JxW)[qp] *
                (*phi_coarse)[i][qp]*(*phi_coarse)[j][qp];
              if (cont == C_ZERO || cont == C_ONE)
                Ke(i,j) += (*JxW)[qp] *
                  (*dphi_coarse)[i][qp]*(*dphi_coarse)[j][qp];
              if (cont == C_ONE)
                Ke(i,j) += (*JxW)[qp] *
                  ((*d2phi_coarse)[i][qp].contract((*d2phi_coarse)[j][qp]));
            }
        }
    }
}

extra_quadrature_order()

void libMesh::HPCoarsenTest::extra_quadrature_order ( const int  extraorder )

inline

Increases or decreases the order of the quadrature rule used for numerical integration.

The default extraorder is 1, because properly integrating L2 error requires integrating the squares of terms with order p+1, and 2p+2 is 1 higher than what we default to using for reasonable mass matrix integration.

Definition at line 115 of file hp_coarsentest.h.

References _extra_order.

  { _extra_order = extraorder; }

operator=() [1/2]

HPCoarsenTest& libMesh::HPCoarsenTest::operator= ( const HPCoarsenTest &  )

delete

operator=() [2/2]

HPCoarsenTest& libMesh::HPCoarsenTest::operator= ( HPCoarsenTest &&  )

default

select_refinement()

void libMesh::HPCoarsenTest::select_refinement ( System &  system )

overridevirtual

This pure virtual function must be redefined in derived classes to take a mesh flagged for h refinement and potentially change the desired refinement type.

Implements libMesh::HPSelector.

Definition at line 138 of file hp_coarsentest.C.

References _extra_order, add_projection(), libMesh::TypeVector< T >::add_scaled(), libMesh::TypeTensor< T >::add_scaled(), libMesh::FEGenericBase< OutputType >::build(), libMesh::C_ONE, libMesh::C_ZERO, libMesh::DenseMatrix< T >::cholesky_solve(), coarse, coarse_qpoints, libMesh::HPSelector::component_scale, libMesh::TypeTensor< T >::contract(), libMesh::System::current_solution(), d2phi, d2phi_coarse, libMesh::FEType::default_quadrature_rule(), dim, libMesh::DISCONTINUOUS, libMesh::Elem::DO_NOTHING, dof_indices, libMesh::DofMap::dof_indices(), libMesh::DofObject::dof_number(), dphi, dphi_coarse, libMesh::ErrorVectorReal, libMesh::SmoothnessEstimator::estimate_smoothness(), fe, Fe, fe_coarse, libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::Elem::hack_p_level(), libMesh::index_range(), libMesh::FEMap::inverse_map(), JxW, Ke, libMesh::libmesh_assert(), libMesh::MeshRefinement::make_flags_parallel_consistent(), libMesh::make_range(), libMesh::ErrorVector::mean(), mesh, libMesh::FEInterface::n_dofs(), n_nodes, n_vars, libMesh::System::n_vars(), libMesh::TensorTools::norm_sq(), libMesh::TypeVector< T >::norm_sq(), libMesh::TypeTensor< T >::norm_sq(), libMesh::System::number(), libMesh::Elem::p_level(), p_weight, libMesh::Elem::parent(), phi, phi_coarse, qrule, libMesh::Real, libMesh::Elem::REFINE, libMesh::DenseVector< T >::resize(), libMesh::DenseMatrix< T >::resize(), libMesh::TypeVector< T >::subtract_scaled(), libMesh::TypeTensor< T >::subtract_scaled(), Uc, Up, libMesh::DofMap::variable_type(), xyz_values, libMesh::DenseMatrix< T >::zero(), libMesh::DenseVector< T >::zero(), and libMesh::zero.

Referenced by main().

{
  LOG_SCOPE("select_refinement()", "HPCoarsenTest");

  // The current mesh
  MeshBase & mesh = system.get_mesh();

  // The dimensionality of the mesh
  const unsigned int dim = mesh.mesh_dimension();

  // The number of variables in the system
  const unsigned int n_vars = system.n_vars();

  // The DofMap for this system
  const DofMap & dof_map = system.get_dof_map();

  // The system number (for doing bad hackery)
  const unsigned int sys_num = system.number();

  // Check for a valid component_scale
  if (!component_scale.empty())
    {
      libmesh_error_msg_if(component_scale.size() != n_vars,
                           "ERROR: component_scale is the wrong size:\n"
                           << " component_scale.size()="
                           << component_scale.size()
                           << "\n n_vars="
                           << n_vars);
    }
  else
    {
      // No specified scaling.  Scale all variables by one.
      component_scale.resize (n_vars, 1.0);
    }

  // Estimates smoothness of solution on each cell
  // via Legendre coefficient decay rate.
  ErrorVector smoothness;
  SmoothnessEstimator estimate_smoothness;
  estimate_smoothness.estimate_smoothness(system, smoothness);

  // Resize the error_per_cell vectors to handle
  // the number of elements, initialize them to 0.
  std::vector<ErrorVectorReal> h_error_per_cell(mesh.max_elem_id(), 0.);
  std::vector<ErrorVectorReal> p_error_per_cell(mesh.max_elem_id(), 0.);

  // Loop over all the variables in the system
  for (unsigned int var=0; var<n_vars; var++)
    {
      // Possibly skip this variable
      if (!component_scale.empty())
        if (component_scale[var] == 0.0) continue;

      // The type of finite element to use for this variable
      const FEType & fe_type = dof_map.variable_type (var);

      // Finite element objects for a fine (and probably a coarse)
      // element will be needed
      fe = FEBase::build (dim, fe_type);
      fe_coarse = FEBase::build (dim, fe_type);

      // Any cached coarse element results have expired
      coarse = nullptr;
      unsigned int cached_coarse_p_level = 0;

      const FEContinuity cont = fe->get_continuity();
      libmesh_assert (cont == DISCONTINUOUS || cont == C_ZERO ||
                      cont == C_ONE);

      // Build an appropriate quadrature rule
      qrule = fe_type.default_quadrature_rule(dim, _extra_order);

      // Tell the refined finite element about the quadrature
      // rule.  The coarse finite element need not know about it
      fe->attach_quadrature_rule (qrule.get());

      // We will always do the integration
      // on the fine elements.  Get their Jacobian values, etc..
      JxW = &(fe->get_JxW());
      xyz_values = &(fe->get_xyz());

      // The shape functions
      phi = &(fe->get_phi());
      phi_coarse = &(fe_coarse->get_phi());

      // The shape function derivatives
      if (cont == C_ZERO || cont == C_ONE)
        {
          dphi = &(fe->get_dphi());
          dphi_coarse = &(fe_coarse->get_dphi());
        }

      // The shape function second derivatives
      if (cont == C_ONE)
        {
#ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
          d2phi = &(fe->get_d2phi());
          d2phi_coarse = &(fe_coarse->get_d2phi());
#else
          libmesh_error_msg("Minimization of H2 error without second derivatives is not possible.");
#endif
        }

      // Iterate over all the active elements in the mesh
      // that live on this processor.
      for (auto & elem : mesh.active_local_element_ptr_range())
        {
          // We're only checking elements that are already flagged for h
          // refinement
          if (elem->refinement_flag() != Elem::REFINE)
            continue;

          const dof_id_type e_id = elem->id();

          // Find the projection onto the parent element,
          // if necessary
          if (elem->parent() &&
              (coarse != elem->parent() ||
               cached_coarse_p_level != elem->p_level()))
            {
              Uc.resize(0);

              coarse = elem->parent();
              cached_coarse_p_level = elem->p_level();

              unsigned int old_parent_level = coarse->p_level();
              coarse->hack_p_level(elem->p_level());

              this->add_projection(system, coarse, var);

              coarse->hack_p_level(old_parent_level);

              // Solve the h-coarsening projection problem
              Ke.cholesky_solve(Fe, Uc);
            }

          fe->reinit(elem);

          // Get the DOF indices for the fine element
          dof_map.dof_indices (elem, dof_indices, var);

          // The number of quadrature points
          const unsigned int n_qp = qrule->n_points();

          // The number of DOFS on the fine element
          const unsigned int n_dofs =
            cast_int<unsigned int>(dof_indices.size());

          // The number of nodes on the fine element
          const unsigned int n_nodes = elem->n_nodes();

          // The average element value (used as an ugly hack
          // when we have nothing p-coarsened to compare to)
          // Real average_val = 0.;
          Number average_val = 0.;

          // Calculate this variable's contribution to the p
          // refinement error

          if (elem->p_level() == 0)
            {
              unsigned int n_vertices = 0;
              for (unsigned int n = 0; n != n_nodes; ++n)
                if (elem->is_vertex(n))
                  {
                    n_vertices++;
                    const Node & node = elem->node_ref(n);
                    average_val += system.current_solution
                      (node.dof_number(sys_num,var,0));
                  }
              average_val /= n_vertices;
            }
          else
            {
              unsigned int old_elem_level = elem->p_level();
              elem->hack_p_level(old_elem_level - 1);

              fe_coarse->reinit(elem, &(qrule->get_points()));

              const unsigned int n_coarse_dofs =
                cast_int<unsigned int>(phi_coarse->size());

              elem->hack_p_level(old_elem_level);

              Ke.resize(n_coarse_dofs, n_coarse_dofs);
              Ke.zero();
              Fe.resize(n_coarse_dofs);
              Fe.zero();

              // Loop over the quadrature points
              for (auto qp : make_range(qrule->n_points()))
                {
                  // The solution value at the quadrature point
                  Number val = libMesh::zero;
                  Gradient grad;
                  Tensor hess;

                  for (unsigned int i=0; i != n_dofs; i++)
                    {
                      dof_id_type dof_num = dof_indices[i];
                      val += (*phi)[i][qp] *
                        system.current_solution(dof_num);
                      if (cont == C_ZERO || cont == C_ONE)
                        grad.add_scaled((*dphi)[i][qp], system.current_solution(dof_num));
                      if (cont == C_ONE)
                        hess.add_scaled((*d2phi)[i][qp], system.current_solution(dof_num));
                    }

                  // The projection matrix and vector
                  for (auto i : index_range(Fe))
                    {
                      Fe(i) += (*JxW)[qp] *
                        (*phi_coarse)[i][qp]*val;
                      if (cont == C_ZERO || cont == C_ONE)
                        Fe(i) += (*JxW)[qp] *
                          grad * (*dphi_coarse)[i][qp];
                      if (cont == C_ONE)
                        Fe(i) += (*JxW)[qp] *
                          hess.contract((*d2phi_coarse)[i][qp]);

                      for (auto j : index_range(Fe))
                        {
                          Ke(i,j) += (*JxW)[qp] *
                            (*phi_coarse)[i][qp]*(*phi_coarse)[j][qp];
                          if (cont == C_ZERO || cont == C_ONE)
                            Ke(i,j) += (*JxW)[qp] *
                              (*dphi_coarse)[i][qp]*(*dphi_coarse)[j][qp];
                          if (cont == C_ONE)
                            Ke(i,j) += (*JxW)[qp] *
                              ((*d2phi_coarse)[i][qp].contract((*d2phi_coarse)[j][qp]));
                        }
                    }
                }

              // Solve the p-coarsening projection problem
              Ke.cholesky_solve(Fe, Up);
            }

          // loop over the integration points on the fine element
          for (unsigned int qp=0; qp<n_qp; qp++)
            {
              Number value_error = 0.;
              Gradient grad_error;
              Tensor hessian_error;
              for (unsigned int i=0; i<n_dofs; i++)
                {
                  const dof_id_type dof_num = dof_indices[i];
                  value_error += (*phi)[i][qp] *
                    system.current_solution(dof_num);
                  if (cont == C_ZERO || cont == C_ONE)
                    grad_error.add_scaled((*dphi)[i][qp], system.current_solution(dof_num));
                  if (cont == C_ONE)
                    hessian_error.add_scaled((*d2phi)[i][qp], system.current_solution(dof_num));
                }
              if (elem->p_level() == 0)
                {
                  value_error -= average_val;
                }
              else
                {
                  for (auto i : index_range(Up))
                    {
                      value_error -= (*phi_coarse)[i][qp] * Up(i);
                      if (cont == C_ZERO || cont == C_ONE)
                        grad_error.subtract_scaled((*dphi_coarse)[i][qp], Up(i));
                      if (cont == C_ONE)
                        hessian_error.subtract_scaled((*d2phi_coarse)[i][qp], Up(i));
                    }
                }

              p_error_per_cell[e_id] += static_cast<ErrorVectorReal>
                (component_scale[var] *
                 (*JxW)[qp] * TensorTools::norm_sq(value_error));
              if (cont == C_ZERO || cont == C_ONE)
                p_error_per_cell[e_id] += static_cast<ErrorVectorReal>
                  (component_scale[var] *
                   (*JxW)[qp] * grad_error.norm_sq());
              if (cont == C_ONE)
                p_error_per_cell[e_id] += static_cast<ErrorVectorReal>
                  (component_scale[var] *
                   (*JxW)[qp] * hessian_error.norm_sq());
            }

          // Calculate this variable's contribution to the h
          // refinement error

          if (!elem->parent())
            {
              // For now, we'll always start with an h refinement
              h_error_per_cell[e_id] =
                std::numeric_limits<ErrorVectorReal>::max() / 2;
            }
          else
            {
              FEMap::inverse_map (dim, coarse, *xyz_values,
                                  coarse_qpoints);

              unsigned int old_parent_level = coarse->p_level();
              coarse->hack_p_level(elem->p_level());

              fe_coarse->reinit(coarse, &coarse_qpoints);

              coarse->hack_p_level(old_parent_level);

              // The number of DOFS on the coarse element
              unsigned int n_coarse_dofs =
                cast_int<unsigned int>(phi_coarse->size());

              // Loop over the quadrature points
              for (unsigned int qp=0; qp<n_qp; qp++)
                {
                  // The solution difference at the quadrature point
                  Number value_error = libMesh::zero;
                  Gradient grad_error;
                  Tensor hessian_error;

                  for (unsigned int i=0; i != n_dofs; ++i)
                    {
                      const dof_id_type dof_num = dof_indices[i];
                      value_error += (*phi)[i][qp] *
                        system.current_solution(dof_num);
                      if (cont == C_ZERO || cont == C_ONE)
                        grad_error.add_scaled((*dphi)[i][qp], system.current_solution(dof_num));
                      if (cont == C_ONE)
                        hessian_error.add_scaled((*d2phi)[i][qp], system.current_solution(dof_num));
                    }

                  for (unsigned int i=0; i != n_coarse_dofs; ++i)
                    {
                      value_error -= (*phi_coarse)[i][qp] * Uc(i);
                      if (cont == C_ZERO || cont == C_ONE)
                        grad_error.subtract_scaled((*dphi_coarse)[i][qp], Uc(i));
                      if (cont == C_ONE)
                        hessian_error.subtract_scaled((*d2phi_coarse)[i][qp], Uc(i));
                    }

                  h_error_per_cell[e_id] += static_cast<ErrorVectorReal>
                    (component_scale[var] *
                     (*JxW)[qp] * TensorTools::norm_sq(value_error));
                  if (cont == C_ZERO || cont == C_ONE)
                    h_error_per_cell[e_id] += static_cast<ErrorVectorReal>
                      (component_scale[var] *
                       (*JxW)[qp] * grad_error.norm_sq());
                  if (cont == C_ONE)
                    h_error_per_cell[e_id] += static_cast<ErrorVectorReal>
                      (component_scale[var] *
                       (*JxW)[qp] * hessian_error.norm_sq());
                }

            }
        }
    }

  // Now that we've got our approximations for p_error and h_error, let's see
  // if we want to switch any h refinement flags to p refinement

  // Iterate over all the active elements in the mesh
  // that live on this processor.
  for (auto & elem : mesh.active_local_element_ptr_range())
    {
      // We're only checking elements that are already flagged for h
      // refinement
      if (elem->refinement_flag() != Elem::REFINE)
        continue;

      const dof_id_type e_id = elem->id();

      unsigned int dofs_per_elem = 0, dofs_per_p_elem = 0;

      // Loop over all the variables in the system
      for (unsigned int var=0; var<n_vars; var++)
        {
          // The type of finite element to use for this variable
          const FEType & fe_type = dof_map.variable_type (var);

          // FIXME: we're overestimating the number of DOFs added by h
          // refinement

          // Compute number of DOFs for elem at current p_level()
          dofs_per_elem +=
            FEInterface::n_dofs(fe_type, elem);

          // Compute number of DOFs for elem at current p_level() + 1
          dofs_per_p_elem +=
            FEInterface::n_dofs(fe_type, elem->p_level() + 1, elem);
        }

      const unsigned int new_h_dofs = dofs_per_elem *
        (elem->n_children() - 1);

      const unsigned int new_p_dofs = dofs_per_p_elem -
        dofs_per_elem;

      /*
        libMesh::err << "Cell " << e_id << ": h = " << elem->hmax()
        << ", p = " << elem->p_level() + 1 << "," << std::endl
        << "     h_error = " << h_error_per_cell[e_id]
        << ", p_error = " << p_error_per_cell[e_id] << std::endl
        << "     new_h_dofs = " << new_h_dofs
        << ", new_p_dofs = " << new_p_dofs << std::endl;
      */
      const Real p_value =
        std::sqrt(p_error_per_cell[e_id]) * p_weight / new_p_dofs;
      const Real h_value =
        std::sqrt(h_error_per_cell[e_id]) /
        static_cast<Real>(new_h_dofs);
      if (smoothness[elem->id()] > smoothness.mean() && p_value > h_value)
        {
          elem->set_p_refinement_flag(Elem::REFINE);
          elem->set_refinement_flag(Elem::DO_NOTHING);
        }
    }

  // libMesh::MeshRefinement will now assume that users have set
  // refinement flags consistently on all processors, but all we've
  // done so far is set refinement flags on local elements.  Let's
  // make sure that flags on geometrically ghosted elements are all
  // set to whatever their owners decided.
  MeshRefinement(mesh).make_flags_parallel_consistent();
}

Member Data Documentation

_extra_order

int libMesh::HPCoarsenTest::_extra_order

protected

Extra order to use for quadrature rule.

Definition at line 178 of file hp_coarsentest.h.

Referenced by extra_quadrature_order(), and select_refinement().

coarse

Elem* libMesh::HPCoarsenTest::coarse

protected

The coarse element on which a solution projection is cached.

Definition at line 128 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

coarse_qpoints

std::vector<Point> libMesh::HPCoarsenTest::coarse_qpoints

protected

Definition at line 156 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

component_scale

std::vector<float> libMesh::HPSelector::component_scale

inherited

This vector can be used to "scale" certain variables in a system.

If the mask is not empty, the consideration given to each component's h and p error estimates will be scaled by component_scale[c].

Definition at line 82 of file hp_selector.h.

Referenced by select_refinement().

d2phi

const std::vector<std::vector<RealTensor> >* libMesh::HPCoarsenTest::d2phi

protected

Definition at line 145 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

d2phi_coarse

const std::vector<std::vector<RealTensor> > * libMesh::HPCoarsenTest::d2phi_coarse

protected

Definition at line 145 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

dof_indices

std::vector<dof_id_type> libMesh::HPCoarsenTest::dof_indices

protected

Global DOF indices for fine elements.

Definition at line 133 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

dphi

const std::vector<std::vector<RealGradient> >* libMesh::HPCoarsenTest::dphi

protected

Definition at line 144 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

dphi_coarse

const std::vector<std::vector<RealGradient> > * libMesh::HPCoarsenTest::dphi_coarse

protected

Definition at line 144 of file hp_coarsentest.h.

Referenced by select_refinement().

fe

std::unique_ptr<FEBase> libMesh::HPCoarsenTest::fe

protected

The finite element objects for fine and coarse elements.

Definition at line 138 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

Fe

DenseVector<Number> libMesh::HPCoarsenTest::Fe

protected

Definition at line 167 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

fe_coarse

std::unique_ptr<FEBase> libMesh::HPCoarsenTest::fe_coarse

protected

Definition at line 138 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

JxW

const std::vector<Real>* libMesh::HPCoarsenTest::JxW

protected

Mapping jacobians.

Definition at line 150 of file hp_coarsentest.h.

Referenced by select_refinement().

Ke

DenseMatrix<Number> libMesh::HPCoarsenTest::Ke

protected

Linear system for projections.

Definition at line 166 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

p_weight

Real libMesh::HPCoarsenTest::p_weight

Because the coarsening test seems to always choose p refinement, we're providing an option to make h refinement more likely.

Definition at line 106 of file hp_coarsentest.h.

Referenced by select_refinement().

phi

const std::vector<std::vector<Real> >* libMesh::HPCoarsenTest::phi

protected

The shape functions and their derivatives.

Definition at line 143 of file hp_coarsentest.h.

Referenced by select_refinement().

phi_coarse

const std::vector<std::vector<Real> > * libMesh::HPCoarsenTest::phi_coarse

protected

Definition at line 143 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

qrule

std::unique_ptr<QBase> libMesh::HPCoarsenTest::qrule

protected

The quadrature rule for the fine element.

Definition at line 161 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

Uc

DenseVector<Number> libMesh::HPCoarsenTest::Uc

protected

Coefficients for projected coarse and projected p-derefined solutions.

Definition at line 172 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

Up

DenseVector<Number> libMesh::HPCoarsenTest::Up

protected

Definition at line 173 of file hp_coarsentest.h.

Referenced by select_refinement().

xyz_values

const std::vector<Point>* libMesh::HPCoarsenTest::xyz_values

protected

Quadrature locations.

Definition at line 155 of file hp_coarsentest.h.

Referenced by add_projection(), and select_refinement().

---

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

• include/error_estimation/hp_coarsentest.h
• src/error_estimation/hp_coarsentest.C