DofObject::set_extra_integer(const unsigned int index,
                             const dof_id_type value)
{
  libmesh_assert_less(index, this->n_extra_integers());
  libmesh_assert_less(this->n_pseudo_systems(), _idx_buf.size());

  const unsigned int start_idx_i = this->start_idx_ints();

  libmesh_assert_less(start_idx_i+index, _idx_buf.size());
  _idx_buf[start_idx_i+index] = value;
}

   * A clone() is needed because GhostingFunctor can not be shared between
   * different meshes. The operations in GhostingFunctor are mesh dependent.
   */
  virtual std::unique_ptr<GhostingFunctor> clone () const override
  { return std::make_unique<OverlapCoupling>(*this); }

  // Change coupling matrix after construction
  void set_dof_coupling(const CouplingMatrix * dof_coupling)

  /**
   * Constructor.
   */
  PointNeighborCoupling(const PointNeighborCoupling & other) :
    GhostingFunctor(other),
    _dof_coupling(other._dof_coupling),
#ifdef LIBMESH_ENABLE_PERIODIC
    _periodic_bcs(other._periodic_bcs),
#endif
    _n_levels(other._n_levels)
  {}

  /**
   * A clone() is needed because GhostingFunctor can not be shared between
   * different meshes. The operations in  GhostingFunctor are mesh dependent.
   */
  virtual std::unique_ptr<GhostingFunctor> clone () const override
  { return std::make_unique<PointNeighborCoupling>(*this); }

  // Change coupling matrix after construction
  void set_dof_coupling(const CouplingMatrix * dof_coupling)

// ------------------------------------------------------------
// Elem class member functions
inline
const std::set<unsigned char> & MeshBase::elem_dimensions() const
{
  libmesh_assert(_preparation.has_cached_elem_data);
  return _elem_dims;
}


inline
const std::set<Order> & MeshBase::elem_default_orders() const
{
  libmesh_assert(_preparation.has_cached_elem_data);
  return _elem_default_orders;
}


inline
Order MeshBase::supported_nodal_order() const
{
  libmesh_assert(_preparation.has_cached_elem_data);
  return _supported_nodal_order;
}

inline
unsigned int MeshBase::spatial_dimension () const
{
  libmesh_assert(_preparation.has_cached_elem_data);

  return cast_int<unsigned int>(_spatial_dimension);
}

template <typename T>

// Declare inline no-ops for assertions with NDEBUG defined
#ifdef NDEBUG

inline void libmesh_assert_valid_is_prepared (const MeshBase &) {}

inline void libmesh_assert_equal_n_systems (const MeshBase &) {}

   * Copy/move ctor, copy/move assignment operator, and destructor are
   * all explicitly defaulted for this simple class.
   */
  QBase (const QBase &) = default;
  QBase (QBase &&) = default;
  QBase & operator= (const QBase &) = default;
  QBase & operator= (QBase &&) = default;

   * Copy/move ctor, copy/move assignment operator, and destructor are
   * all explicitly defaulted for this simple class.
   */
  QSimpson (const QSimpson &) = default;
  QSimpson (QSimpson &&) = default;
  QSimpson & operator= (const QSimpson &) = default;
  QSimpson & operator= (QSimpson &&) = default;

  // Swap out the new mid elem node with the previous one
  if (mid_elem_node)
  {
    libmesh_assert(returnval->n_nodes() == this->n_nodes());
    returnval->set_node(returnval->n_nodes() - 1, this->_nodelinks_data.back());
  }

  returnval->set_id() = this->id();

  if (this->interior_parent())
    {
      if (!rhs.interior_parent())
        return false;
      if (this->interior_parent()->id() !=
          rhs.interior_parent()->id())
        return false;
    }
  else if (rhs.interior_parent())

OverlapCoupling::OverlapCoupling(const OverlapCoupling & other) :
  GhostingFunctor(other),
  _dof_coupling(other._dof_coupling)
{
  for (auto i : make_range(2))
    if (other._q_rules[i].get())
      _q_rules[i] = other._q_rules[i]->clone();
}

  return *this;
}

std::string_view DistributedMesh::subclass_first_difference_from (const MeshBase & other_mesh_base) const
{
  const DistributedMesh * dist_mesh_ptr =
    dynamic_cast<const DistributedMesh *>(&other_mesh_base);
  if (!dist_mesh_ptr)
    return "DistributedMesh class";
  const DistributedMesh & other_mesh = *dist_mesh_ptr;

#define CHECK_MEMBER(member_name) \
  if (member_name != other_mesh.member_name) \
    return #member_name;

  CHECK_MEMBER(_is_serial);
  CHECK_MEMBER(_is_serial_on_proc_0);
  CHECK_MEMBER(_deleted_coarse_elements);
  CHECK_MEMBER(_n_nodes);
  CHECK_MEMBER(_n_elem);
  CHECK_MEMBER(_max_node_id);
  CHECK_MEMBER(_max_elem_id);
  // We expect these things to change in a prepare_for_use();
  // they're conceptually "mutable"...
/*

  CHECK_MEMBER(_next_unpartitioned_unique_id);
#endif
*/
  if (!this->nodes_and_elements_equal(other_mesh))
    return "nodes and/or elements";

  if (_extra_ghost_elems.size() !=
      other_mesh._extra_ghost_elems.size())
    return "_extra_ghost_elems size";
  for (auto & elem : _extra_ghost_elems)
    {
      libmesh_assert(this->query_elem_ptr(elem->id()) == elem);
      const Elem * other_elem = other_mesh.query_elem_ptr(elem->id());
      if (!other_elem ||
          !other_mesh._extra_ghost_elems.count(const_cast<Elem *>(other_elem)))
        return "_extra_ghost_elems entry";
    }

  return "";
}

DistributedMesh::~DistributedMesh ()

      // We communicate valid neighbor links for newly-redistributed
      // elements, but we may still have remote_elem links that should
      // be corrected but whose correction wasn't communicated.
      this->find_neighbors(/*reset_remote_elements*/ false,
                           /*reset_current_list*/ false /*non-default*/,
                           /*assert_valid*/ false,
                           /*check_non_remote*/ false /*non-default!*/);

      // Is this necessary?  If we are called from prepare_for_use(), this will be called
      // anyway... but users can always call partition directly, in which case we do need

                      sender_could_become_owner)
                    {
                      if (it != repartitioned_node_pids.end() &&
                          pid < it->second)
                        it->second = pid;
                      else
                        repartitioned_node_pids[n] = pid;
                    }
                  else
                    if (it == repartitioned_node_pids.end())
                      repartitioned_node_pids[n] =
                        DofObject::invalid_processor_id;

                  repartitioned_node_sets_to_push[pid].insert(n);

  // _elemset_codes stores pointers to entries in _elemset_codes_inverse_map,
  // so it is not possible to simply copy it directly from other_mesh
  for (const auto & [set, code] : _elemset_codes_inverse_map)
    _elemset_codes.emplace(code, &set);
}

MeshBase& MeshBase::operator= (MeshBase && other_mesh)

}


void MeshBase::assert_equal_to (const MeshBase & other_mesh,
                                std::string_view failure_context) const
{
#ifndef NDEBUG
  LOG_SCOPE("assert_equal_to()", "MeshBase");

  std::string_view local_diff = first_difference_from(other_mesh);

  bool diff_found = !local_diff.empty();
  this->comm().max(diff_found);

  if (diff_found)
    {
      // Construct a user-friendly message to throw on pid 0
      std::set<std::string> unique_diffs;
      if (!local_diff.empty())
        unique_diffs.insert(std::string(local_diff));
      this->comm().set_union(unique_diffs);

      if (!this->processor_id())
        {
          std::string error_msg {failure_context};
          error_msg += "\nMeshes failed asserted equality in at least these aspects:\n";
          for (auto & diff : unique_diffs)
            {
              error_msg += diff;
              error_msg += '\n';
            }
          libmesh_assert_msg(!diff_found, error_msg);
        }

      // We're not going to throw on other processors because we don't

#else
  libmesh_ignore(other_mesh, failure_context);
#endif // NDEBUG
}


bool MeshBase::locally_equals (const MeshBase & other_mesh) const
{
  const std::string_view diff = first_difference_from(other_mesh);
  return diff.empty();
}


std::string_view MeshBase::first_difference_from(const MeshBase & other_mesh) const
{
  // Check whether (almost) everything in the base is equal
  //

  if (member_name != other_mesh.member_name) \
    return #member_name;

  CHECK_MEMBER(_n_parts);
  CHECK_MEMBER(_default_mapping_type);
  CHECK_MEMBER(_default_mapping_data);
  CHECK_MEMBER(_preparation);
  CHECK_MEMBER(_count_lower_dim_elems_in_point_locator);

  // We should either both have our own interior parents or both not;
  // but if we both don't then we can't really assert anything else

  // pointing at two different copies of "the same" interior mesh.
  if ((_interior_mesh == this) !=
      (other_mesh._interior_mesh == &other_mesh))
    return "_interior_mesh";

  CHECK_MEMBER(_skip_noncritical_partitioning);
  CHECK_MEMBER(_skip_all_partitioning);
  CHECK_MEMBER(_skip_renumber_nodes_and_elements);
  CHECK_MEMBER(_skip_find_neighbors);
  CHECK_MEMBER(_skip_detect_interior_parents);
  CHECK_MEMBER(_allow_remote_element_removal);
  CHECK_MEMBER(_allow_node_and_elem_unique_id_overlap);
  CHECK_MEMBER(_spatial_dimension);
  CHECK_MEMBER(_point_locator_close_to_point_tol);
  CHECK_MEMBER(_block_id_to_name);
  CHECK_MEMBER(_elem_dims);
  CHECK_MEMBER(_elem_default_orders);
  CHECK_MEMBER(_supported_nodal_order);
  CHECK_MEMBER(_mesh_subdomains);
  CHECK_MEMBER(_all_elemset_ids);
  CHECK_MEMBER(_elem_integer_names);
  CHECK_MEMBER(_elem_integer_default_values);
  CHECK_MEMBER(_node_integer_names);
  CHECK_MEMBER(_node_integer_default_values);

  if (static_cast<bool>(_default_ghosting) != static_cast<bool>(other_mesh._default_ghosting))
    return "_default_ghosting";
  if (static_cast<bool>(_partitioner) != static_cast<bool>(other_mesh._partitioner))
    return "_partitioner";
  if (*boundary_info != *other_mesh.boundary_info)
    return "boundary_info";

  // First check whether the "existence" of the two pointers differs (one present, one absent)
  if (static_cast<bool>(_disjoint_neighbor_boundary_pairs) !=
      static_cast<bool>(other_mesh._disjoint_neighbor_boundary_pairs))
    return "_disjoint_neighbor_boundary_pairs existence";
  // If both exist, compare the contents (Weak Test: just compare sizes like `_ghosting_functors`)
  if (_disjoint_neighbor_boundary_pairs &&
      (_disjoint_neighbor_boundary_pairs->size() != other_mesh._disjoint_neighbor_boundary_pairs->size()))
    return "_disjoint_neighbor_boundary_pairs size";

  const constraint_rows_type & other_rows =
    other_mesh.get_constraint_rows();

      const dof_id_type node_id = node->id();
      const Node * other_node = other_mesh.query_node_ptr(node_id);
      if (!other_node)
        return "_constraint_rows node presence";

      auto it = other_rows.find(other_node);
      if (it == other_rows.end())
        return "_constraint_rows row presence";

      const auto & other_row = it->second;
      if (row.size() != other_row.size())
        return "_constraint_rows row size";

      for (auto i : index_range(row))
        {

              elem_pair.second !=
              other_elem_pair.second ||
              coef != other_coef)
            return "_constraint_rows row entry";
        }
    }

  for (const auto & [elemset_code, elemset_ptr] : this->_elemset_codes)
    if (const auto it = other_mesh._elemset_codes.find(elemset_code);
        it == other_mesh._elemset_codes.end() || *elemset_ptr != *it->second)
      return "_elemset_codes";

  // FIXME: we have no good way to compare ghosting functors, since
  // they're in a vector of pointers, and we have no way *at all*

  // we have the same number, is all.
  if (_ghosting_functors.size() !=
      other_mesh._ghosting_functors.size())
    return "_ghosting_functors size";

  // Same deal for partitioners.  We tested that we both have one or
  // both don't, but are they equivalent?  Let's guess "yes".

  // Now let the subclasses decide whether everything else is equal
  return this->subclass_first_difference_from(other_mesh);
}

  libmesh_assert(this->comm().verify(this->is_serial()));

  _preparation.libmesh_assert_consistent(this->comm());

#ifdef DEBUG
  // If we don't go into this method with valid constraint rows, we're

  // If this mesh thinks it's already  partially prepared, then in
  // optimized builds we'll trust it, but in debug builds we'll check.
  const bool was_partly_prepared = (_preparation == Preparation());
#endif

  // A distributed mesh may have processors with no elements (or

#ifdef DEBUG
  // The if() here avoids both unnecessary work *and* stack overflow
  if (was_partly_prepared)
    MeshTools::libmesh_assert_valid_is_prepared(*this);

  MeshTools::libmesh_assert_valid_boundary_ids(*this);
#ifdef LIBMESH_ENABLE_UNIQUE_ID

                --_elem_dims.end(), // --end() is valid if the set is non-empty
                std::ostream_iterator<unsigned int>(oss, ", "));
      oss << cast_int<unsigned int>(*_elem_dims.rbegin());
      oss << "}";
      if (!this->preparation().has_cached_elem_data)
        oss << " (may be out of date)";
      oss << '\n';
    }

  if (!_elem_default_orders.empty())

                     [](Order o)
                       { return Utility::enum_to_string<Order>(o); });
      oss << Utility::enum_to_string<Order>(*_elem_default_orders.rbegin());
      oss << "}";
      if (!this->preparation().has_cached_elem_data)
        oss << " (may be out of date)";
      oss << '\n';
    }

  oss << "  supported_nodal_order()=" << this->_supported_nodal_order;
  if (!this->preparation().has_cached_elem_data)
    oss << " (may be out of date)";
  oss << '\n';

  oss << "  spatial_dimension()=" << int(this->_spatial_dimension);
  if (!this->preparation().has_cached_elem_data)
    oss << " (may be out of date)";
  oss << '\n';

  oss << "  n_nodes()="               << this->n_nodes()                                      << '\n'
      << "    n_local_nodes()="       << this->n_local_nodes()                                << '\n'
      << "  n_elem()="                << this->n_elem()                                       << '\n'
      << "    n_local_elem()="        << this->n_local_elem()                                 << '\n';

MeshBase::Preparation::operator== (const Preparation & other) const
{
  if (is_partitioned != other.is_partitioned)
    return false;
  if (has_synched_id_counts != other.has_synched_id_counts)
    return false;
  if (has_neighbor_ptrs != other.has_neighbor_ptrs)
    return false;
  if (has_cached_elem_data != other.has_cached_elem_data)
    return false;
  if (has_interior_parent_ptrs != other.has_interior_parent_ptrs)
    return false;
  if (has_removed_remote_elements != other.has_removed_remote_elements)
    return false;
  if (has_removed_orphaned_nodes != other.has_removed_orphaned_nodes)
    return false;
  if (has_reinit_ghosting_functors != other.has_reinit_ghosting_functors)

void
MeshBase::Preparation::libmesh_assert_consistent (const Parallel::Communicator & libmesh_dbg_var(comm))
{
  libmesh_assert(comm.verify(is_partitioned));
  libmesh_assert(comm.verify(has_synched_id_counts));
  libmesh_assert(comm.verify(has_neighbor_ptrs));
  libmesh_assert(comm.verify(has_cached_elem_data));
  libmesh_assert(comm.verify(has_interior_parent_ptrs));
  libmesh_assert(comm.verify(has_removed_remote_elements));
  libmesh_assert(comm.verify(has_removed_orphaned_nodes));
  libmesh_assert(comm.verify(has_reinit_ghosting_functors));
  libmesh_assert(comm.verify(has_boundary_id_sets));
}


// Explicit instantiations for our template function

void VariationalMeshSmoother::setup()
{
  // Check for multiple dimensions
  if (!_mesh.preparation().has_cached_elem_data)
    _mesh.cache_elem_data();

  if (_mesh.elem_dimensions().size() > 1)
    libmesh_not_implemented_msg("Meshes containing elements of differing dimension are not yet supported.");

{
  // If we've been handed an unprepared mesh then we need to be made
  // aware of that and fix that; we're relying on neighbor pointers.
  MeshTools::libmesh_assert_valid_is_prepared(mesh);

  if (!mesh.is_prepared())
    mesh.prepare_for_use();

std::unique_ptr<MeshBase> reprepared_mesh_clone (const MeshBase & mesh)
{
  const MeshBase::Preparation prep = mesh.preparation();

  // If the mesh thinks it's prepared in some way, *re*-preparing in
  // that way shouldn't change a clone of it, as long as we disallow
  // repartitioning or renumbering or remote element removal.
  std::unique_ptr<MeshBase> mesh_clone = mesh.clone();

  const bool old_allow_renumbering = mesh_clone->allow_renumbering();
  const bool old_allow_remote_element_removal =
    mesh_clone->allow_remote_element_removal();
  const bool old_skip_partitioning = mesh_clone->skip_partitioning();
  mesh_clone->allow_renumbering(false);
  mesh_clone->allow_remote_element_removal(false);
  mesh_clone->skip_partitioning(true);

  // If the mesh thinks it's already completely prepared, test that
  if (prep)
    mesh_clone->prepare_for_use();
  // If the mesh thinks it's somewhat prepared, test each way it
  // thinks so.
  else
    {
      if (prep.has_synched_id_counts)
        mesh_clone->update_parallel_id_counts();

      if (prep.has_neighbor_ptrs)
        mesh_clone->find_neighbors();

      if (prep.has_cached_elem_data)
        mesh_clone->cache_elem_data();

      if (mesh.allow_detect_interior_parents() &&
          prep.has_interior_parent_ptrs)
        mesh_clone->detect_interior_parents();

      if (old_allow_remote_element_removal &&
          prep.has_removed_remote_elements)
        mesh_clone->delete_remote_elements();

      if (prep.has_removed_orphaned_nodes)
        mesh_clone->remove_orphaned_nodes();

      if (prep.has_boundary_id_sets)
        mesh_clone->get_boundary_info().regenerate_id_sets();

      // I don't know how we'll tell if this changes anything, but
      // we'll do it for completeness
      if (prep.has_reinit_ghosting_functors)
        mesh_clone->reinit_ghosting_functors();
    }

  // Restore original flag values
  mesh_clone->allow_renumbering(old_allow_renumbering);
  mesh_clone->allow_remote_element_removal(old_allow_remote_element_removal);
  mesh_clone->skip_partitioning(old_skip_partitioning);

  return mesh_clone;
}


}

{
  LOG_SCOPE("valid_is_prepared()", "MeshTools");

  const MeshBase::Preparation prep = mesh.preparation();

  // If the mesh doesn't think *anything* has been prepared, we have
  // nothing to check.
  if (prep == MeshBase::Preparation())
    return true;

  // If the mesh thinks it's partitioned, check.  These are counts,
  // not caches, so it's a real check.
  if (prep.is_partitioned)
    if (mesh.n_unpartitioned_elem() ||
        mesh.n_unpartitioned_nodes())
      return false;

  // If the mesh thinks it's prepared in some way, *re*-preparing in
  // that way shouldn't change a clone of it, as long as we disallow
  // repartitioning or renumbering or remote element removal.
  std::unique_ptr<MeshBase> mesh_clone = reprepared_mesh_clone(mesh);

  // Check whether the original and clone compare equal
  return (mesh == *mesh_clone);

#ifndef NDEBUG


void libmesh_assert_valid_is_prepared (const MeshBase & mesh)
{
  LOG_SCOPE("libmesh_assert_valid_is_prepared()", "MeshTools");

  const MeshBase::Preparation prep = mesh.preparation();

  // If the mesh doesn't think *anything* has been prepared, we have
  // nothing to check.
  if (prep == MeshBase::Preparation())
    return;

  // If the mesh thinks it's partitioned, check.  These are counts,
  // not caches, so it's a real check.
  if (prep.is_partitioned)
    libmesh_assert_msg
      (!mesh.n_unpartitioned_elem() && !mesh.n_unpartitioned_nodes(),
       "Mesh preparation().is_partitioned does not reflect mesh data.");

  // If the mesh thinks it's prepared in some way, *re*-preparing in
  // that way shouldn't change a clone of it, as long as we disallow
  // repartitioning or renumbering or remote element removal.
  std::unique_ptr<MeshBase> mesh_clone = reprepared_mesh_clone(mesh);

  mesh.assert_equal_to
    (*mesh_clone,
     "Mesh data does not match mesh preparation().\n"
     "Efficiently-prepared mesh != exhaustively-prepared clone.");
}

  LOG_SCOPE("libmesh_assert_parallel_consistent_procids()", "MeshTools");

  if (mesh.n_processors() == 1)
    return;

  libmesh_parallel_only(mesh.comm());

}


std::string_view ReplicatedMesh::subclass_first_difference_from (const MeshBase & other_mesh_base) const
{
  const ReplicatedMesh * rep_mesh_ptr =
    dynamic_cast<const ReplicatedMesh *>(&other_mesh_base);
  if (!rep_mesh_ptr)
    return "ReplicatedMesh class";
  const ReplicatedMesh & other_mesh = *rep_mesh_ptr;

#define CHECK_MEMBER(member_name) \
  if (member_name != other_mesh.member_name) \
    return #member_name;

  CHECK_MEMBER(_n_nodes);
  CHECK_MEMBER(_n_elem);
#ifdef LIBMESH_ENABLE_UNIQUE_ID
  CHECK_MEMBER(_next_unique_id);
#endif
  if (!this->nodes_and_elements_equal(other_mesh))
    return "nodes and/or elements";

  return "";
}

           });

        MeshBase * other_interior_mesh =
          const_cast<MeshBase *>(&other_mesh.interior_mesh());

        // If we don't already have interior parents, then we can just
        // use whatever interior_mesh we need for the incoming

        if (!existing_interior_parents)
          {
            if (other_interior_mesh == &other_mesh)
              this->set_interior_mesh(*this);
            else
              this->set_interior_mesh(*other_interior_mesh);
          }

        if (other_interior_mesh == &other_mesh &&

          for (auto & elem_pair : ip_map)
            elem_pair.second->set_interior_parent(
              this->elem_ptr(elem_pair.first->interior_parent()->id() + element_id_offset));
        else if (other_interior_mesh == _interior_mesh)
          for (auto & elem_pair : ip_map)
            {
              Elem * ip = const_cast<Elem *>(elem_pair.first->interior_parent());
              libmesh_assert(ip == remote_elem ||
                             ip == other_interior_mesh->elem_ptr(ip->id()));
              elem_pair.second->set_interior_parent(ip);
            }
        else
          libmesh_error_msg("Cannot copy boundary elements between meshes with different interior meshes");

      const bool skipped_partitioning = this->skip_partitioning();
      this->skip_partitioning(true);

      const Preparation old_preparation = this->preparation();
      this->prepare_for_use();

      //But in the long term, don't change our policies.

      // That prepare_for_use() call marked us as prepared, but we
      // specifically avoided some important preparation, so we might not
      // actually be prepared now.
      if (skip_find_neighbors)
        this->unset_has_neighbor_ptrs();

      const Preparation other_preparation = other_mesh.preparation();
      if (!old_preparation.is_partitioned ||
          !other_preparation.is_partitioned)
        this->unset_is_partitioned();
      if (!old_preparation.has_removed_orphaned_nodes ||
          !other_preparation.has_removed_orphaned_nodes)
        this->unset_has_removed_orphaned_nodes();
      if (!old_preparation.has_removed_remote_elements ||
          !other_preparation.has_removed_remote_elements)
        this->unset_has_removed_remote_elements();
    }

  // In general we've just invalidated just about everything, and we'd

            //
            // If we're only checking remote neighbors, after a
            // redistribution, then we'll skip the non-remote ones
            if ((element->neighbor_ptr(ms) == nullptr && check_non_remote) ||
                element->neighbor_ptr(ms) == remote_elem)
              {
                // Get the key for the side of this element.  Use the

  return QSIMPSON;
}

std::unique_ptr<QBase> QSimpson::clone() const
{
  return std::make_unique<QSimpson>(*this);
}

// See the files:

  // unprepared without re-preparing it.  At *this* point hopefully
  // everybody's done changing the mesh and we can make sure it's
  // prepared.
  if (!_mesh.is_prepared())
    _mesh.complete_preparation();

  // Tell all the \p DofObject entities how many systems
  // there are.

  const MeshBase & mesh = this->get_mesh();

  libmesh_assert(mesh.is_prepared());
#if defined(DEBUG) && !defined(LIBMESH_ENABLE_DEPRECATED)
  MeshTools::libmesh_assert_valid_is_prepared(mesh);
#endif

  const MeshBase & mesh = this->get_mesh();

  libmesh_assert(mesh.is_prepared());
#if defined(DEBUG) && !defined(LIBMESH_ENABLE_DEPRECATED)
  MeshTools::libmesh_assert_valid_is_prepared(mesh);
#endif

  const MeshBase & mesh = this->get_mesh();

  libmesh_assert(mesh.is_prepared());
#if defined(DEBUG) && !defined(LIBMESH_ENABLE_DEPRECATED)
  MeshTools::libmesh_assert_valid_is_prepared(mesh);
#endif

    es.parameters.get<unsigned int>("nonlinear solver maximum function evaluations");

  const double abs_resid_tol = parameters.have_parameter<Real>("nonlinear solver absolute residual tolerance") ?
    double(parameters.get<Real>("nonlinear solver absolute residual tolerance")) :
    double(es.parameters.get<Real>("nonlinear solver absolute residual tolerance"));

  const double rel_resid_tol = parameters.have_parameter<Real>("nonlinear solver relative residual tolerance") ?

                                       bool /*apply_heterogeneous_constraints*/,
                                       bool /*apply_no_constraints*/)
{
  libmesh_assert(this->get_mesh().is_prepared());
#if defined(DEBUG) && !defined(LIBMESH_ENABLE_DEPRECATED)
  MeshTools::libmesh_assert_valid_is_prepared(this->get_mesh());
#endif

  // Hopefully the user or the EquationSystems prepared the mesh, but
  // if not then we'd better do it now.
  if (!mesh.is_prepared())
    mesh.complete_preparation();

  // Distribute the degrees of freedom on the mesh
  auto total_dofs = _dof_map->distribute_dofs (mesh);

  // Log how long the user's assembly code takes
  LOG_SCOPE("assemble()", "System");

  libmesh_assert(this->get_mesh().is_prepared());
#if defined(DEBUG) && !defined(LIBMESH_ENABLE_DEPRECATED)
  MeshTools::libmesh_assert_valid_is_prepared(this->get_mesh());
#endif

  // Log how long the user's assembly code takes
  LOG_SCOPE("assemble_qoi()", "System");

  libmesh_assert(this->get_mesh().is_prepared());
#if defined(DEBUG) && !defined(LIBMESH_ENABLE_DEPRECATED)
  MeshTools::libmesh_assert_valid_is_prepared(this->get_mesh());
#endif

  // Log how long the user's assembly code takes
  LOG_SCOPE("assemble_qoi_derivative()", "System");

  libmesh_assert(this->get_mesh().is_prepared());
#if defined(DEBUG) && !defined(LIBMESH_ENABLE_DEPRECATED)
  MeshTools::libmesh_assert_valid_is_prepared(this->get_mesh());
#endif