libMesh
Public Member Functions | Public Attributes | Protected Member Functions | List of all members
libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectVertices Struct Reference

#include <generic_projector.h>

Inheritance diagram for libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectVertices:
[legend]

Public Member Functions

 ProjectVertices (GenericProjector &p)
 
 ProjectVertices (ProjectVertices &p_v, Threads::split)
 
void operator() (const node_range &range)
 
void insert_id (dof_id_type id, const FValue &val, processor_id_type pid)
 
void insert_ids (const std::vector< dof_id_type > &ids, const std::vector< FValue > &vals, processor_id_type pid)
 
void construct_projection (const std::vector< dof_id_type > &dof_indices_var, const std::vector< unsigned int > &involved_dofs, unsigned int var_component, const Node *node, const FEBase &fe)
 
void insert_id (dof_id_type id, const FValue &val, processor_id_type pid)
 
void insert_ids (const std::vector< dof_id_type > &ids, const std::vector< FValue > &vals, processor_id_type pid)
 
void find_dofs_to_send (const Node &node, const Elem &elem, unsigned short node_num, const var_set &vars)
 
void join (const SubFunctor &other)
 

Public Attributes

ProjectionAction action
 
FFunctor f
 
std::unique_ptr< GFunctor > g
 
const Systemsystem
 
FEMContext context
 
std::vector< FEContinuityconts
 
GenericProjectorprojector
 
std::unordered_map< dof_id_type, std::pair< FValue, processor_id_type > > new_ids_to_push
 
std::unordered_map< dof_id_type, FValue > new_ids_to_save
 

Protected Member Functions

void construct_projection (const std::vector< dof_id_type > &dof_indices_var, const std::vector< unsigned int > &involved_dofs, unsigned int var_component, const Node *node, const FEBase &fe)
 

Detailed Description

template<typename FFunctor, typename GFunctor, typename FValue, typename ProjectionAction>
struct libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectVertices

Definition at line 262 of file generic_projector.h.

Constructor & Destructor Documentation

◆ ProjectVertices() [1/2]

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectVertices::ProjectVertices ( GenericProjector p)

Definition at line 263 of file generic_projector.h.

263 : SubProjector(p) {}

◆ ProjectVertices() [2/2]

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectVertices::ProjectVertices ( ProjectVertices p_v,
Threads::split   
)

Definition at line 265 of file generic_projector.h.

265 : SubProjector(p_v.projector) {}

Member Function Documentation

◆ construct_projection() [1/2]

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
void libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubProjector::construct_projection ( const std::vector< dof_id_type > &  dof_indices_var,
const std::vector< unsigned int > &  involved_dofs,
unsigned int  var_component,
const Node node,
const FEBase fe 
)
protectedinherited

Definition at line 2478 of file generic_projector.h.

References libMesh::C_ONE, libMesh::DenseMatrix< T >::cholesky_solve(), libMesh::FEAbstract::get_continuity(), libMesh::FEGenericBase< OutputType >::get_dphi(), libMesh::FEAbstract::get_JxW(), libMesh::FEGenericBase< OutputType >::get_phi(), libMesh::DenseVector< T >::get_values(), libMesh::FEAbstract::get_xyz(), libMesh::DofObject::invalid_processor_id, libMesh::DofObject::processor_id(), and libMesh::DenseVector< T >::size().

2483 {
2484  const std::vector<Real> & JxW = fe.get_JxW();
2485  const std::vector<std::vector<Real>> & phi = fe.get_phi();
2486  const std::vector<std::vector<RealGradient>> * dphi = nullptr;
2487  const std::vector<Point> & xyz_values = fe.get_xyz();
2488  const FEContinuity cont = fe.get_continuity();
2489  const std::unordered_map<dof_id_type, FValue> & ids_to_save =
2490  this->projector.ids_to_save;
2491 
2492  if (cont == C_ONE)
2493  dphi = &(fe.get_dphi());
2494 
2495  const unsigned int n_involved_dofs =
2496  cast_int<unsigned int>(involved_dofs.size());
2497 
2498  std::vector<dof_id_type> free_dof_ids;
2499  DenseVector<FValue> Uinvolved(n_involved_dofs);
2500  std::vector<char> dof_is_fixed(n_involved_dofs, false); // bools
2501 
2502  for (auto i : IntRange<unsigned int>(0, n_involved_dofs))
2503  {
2504  const dof_id_type id = dof_indices_var[involved_dofs[i]];
2505  auto iter = ids_to_save.find(id);
2506  if (iter == ids_to_save.end())
2507  free_dof_ids.push_back(id);
2508  else
2509  {
2510  dof_is_fixed[i] = true;
2511  Uinvolved(i) = iter->second;
2512  }
2513  }
2514 
2515  const unsigned int free_dofs = free_dof_ids.size();
2516 
2517  // There may be nothing to project
2518  if (!free_dofs)
2519  return;
2520 
2521  // The element matrix and RHS for projections.
2522  // Note that Ke is always real-valued, whereas
2523  // Fe may be complex valued if complex number
2524  // support is enabled
2525  DenseMatrix<Real> Ke(free_dofs, free_dofs);
2526  DenseVector<FValue> Fe(free_dofs);
2527  // The new degree of freedom coefficients to solve for
2528  DenseVector<FValue> Ufree(free_dofs);
2529 
2530  const unsigned int n_qp =
2531  cast_int<unsigned int>(xyz_values.size());
2532 
2533  // Loop over the quadrature points
2534  for (unsigned int qp=0; qp<n_qp; qp++)
2535  {
2536  // solution at the quadrature point
2537  FValue fineval = f.eval_at_point(context,
2538  var_component,
2539  xyz_values[qp],
2540  system.time);
2541  // solution grad at the quadrature point
2542  VectorValue<FValue> finegrad;
2543  if (cont == C_ONE)
2544  finegrad = g->eval_at_point(context,
2545  var_component,
2546  xyz_values[qp],
2547  system.time);
2548 
2549  // Form edge projection matrix
2550  for (unsigned int sidei=0, freei=0;
2551  sidei != n_involved_dofs; ++sidei)
2552  {
2553  unsigned int i = involved_dofs[sidei];
2554  // fixed DoFs aren't test functions
2555  if (dof_is_fixed[sidei])
2556  continue;
2557  for (unsigned int sidej=0, freej=0;
2558  sidej != n_involved_dofs; ++sidej)
2559  {
2560  unsigned int j = involved_dofs[sidej];
2561  if (dof_is_fixed[sidej])
2562  Fe(freei) -= phi[i][qp] * phi[j][qp] *
2563  JxW[qp] * Uinvolved(sidej);
2564  else
2565  Ke(freei,freej) += phi[i][qp] *
2566  phi[j][qp] * JxW[qp];
2567  if (cont == C_ONE)
2568  {
2569  if (dof_is_fixed[sidej])
2570  Fe(freei) -= ( (*dphi)[i][qp] *
2571  (*dphi)[j][qp] ) *
2572  JxW[qp] * Uinvolved(sidej);
2573  else
2574  Ke(freei,freej) += ( (*dphi)[i][qp] *
2575  (*dphi)[j][qp] )
2576  * JxW[qp];
2577  }
2578  if (!dof_is_fixed[sidej])
2579  freej++;
2580  }
2581  Fe(freei) += phi[i][qp] * fineval * JxW[qp];
2582  if (cont == C_ONE)
2583  Fe(freei) += (finegrad * (*dphi)[i][qp] ) *
2584  JxW[qp];
2585  freei++;
2586  }
2587  }
2588 
2589  Ke.cholesky_solve(Fe, Ufree);
2590 
2591  // Transfer new edge solutions to element
2592  const processor_id_type pid = node ?
2593  node->processor_id() : DofObject::invalid_processor_id;
2594  insert_ids(free_dof_ids, Ufree.get_values(), pid);
2595 }
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
Definition: system.h:1545
uint8_t processor_id_type
Definition: id_types.h:104
static const processor_id_type invalid_processor_id
An invalid processor_id to distinguish DoFs that have not been assigned to a processor.
Definition: dof_object.h:404
std::unordered_map< dof_id_type, FValue > ids_to_save
FEContinuity
defines an enum for finite element types to libmesh_assert a certain level (or type? Hcurl?) of continuity.
void insert_ids(const std::vector< dof_id_type > &ids, const std::vector< FValue > &vals, processor_id_type pid)
uint8_t dof_id_type
Definition: id_types.h:67

◆ construct_projection() [2/2]

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
void libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubProjector::construct_projection

Definition at line 2478 of file generic_projector.h.

2483 {
2484  const std::vector<Real> & JxW = fe.get_JxW();
2485  const std::vector<std::vector<Real>> & phi = fe.get_phi();
2486  const std::vector<std::vector<RealGradient>> * dphi = nullptr;
2487  const std::vector<Point> & xyz_values = fe.get_xyz();
2488  const FEContinuity cont = fe.get_continuity();
2489  const std::unordered_map<dof_id_type, FValue> & ids_to_save =
2490  this->projector.ids_to_save;
2491 
2492  if (cont == C_ONE)
2493  dphi = &(fe.get_dphi());
2494 
2495  const unsigned int n_involved_dofs =
2496  cast_int<unsigned int>(involved_dofs.size());
2497 
2498  std::vector<dof_id_type> free_dof_ids;
2499  DenseVector<FValue> Uinvolved(n_involved_dofs);
2500  std::vector<char> dof_is_fixed(n_involved_dofs, false); // bools
2501 
2502  for (auto i : IntRange<unsigned int>(0, n_involved_dofs))
2503  {
2504  const dof_id_type id = dof_indices_var[involved_dofs[i]];
2505  auto iter = ids_to_save.find(id);
2506  if (iter == ids_to_save.end())
2507  free_dof_ids.push_back(id);
2508  else
2509  {
2510  dof_is_fixed[i] = true;
2511  Uinvolved(i) = iter->second;
2512  }
2513  }
2514 
2515  const unsigned int free_dofs = free_dof_ids.size();
2516 
2517  // There may be nothing to project
2518  if (!free_dofs)
2519  return;
2520 
2521  // The element matrix and RHS for projections.
2522  // Note that Ke is always real-valued, whereas
2523  // Fe may be complex valued if complex number
2524  // support is enabled
2525  DenseMatrix<Real> Ke(free_dofs, free_dofs);
2526  DenseVector<FValue> Fe(free_dofs);
2527  // The new degree of freedom coefficients to solve for
2528  DenseVector<FValue> Ufree(free_dofs);
2529 
2530  const unsigned int n_qp =
2531  cast_int<unsigned int>(xyz_values.size());
2532 
2533  // Loop over the quadrature points
2534  for (unsigned int qp=0; qp<n_qp; qp++)
2535  {
2536  // solution at the quadrature point
2537  FValue fineval = f.eval_at_point(context,
2538  var_component,
2539  xyz_values[qp],
2540  system.time);
2541  // solution grad at the quadrature point
2542  VectorValue<FValue> finegrad;
2543  if (cont == C_ONE)
2544  finegrad = g->eval_at_point(context,
2545  var_component,
2546  xyz_values[qp],
2547  system.time);
2548 
2549  // Form edge projection matrix
2550  for (unsigned int sidei=0, freei=0;
2551  sidei != n_involved_dofs; ++sidei)
2552  {
2553  unsigned int i = involved_dofs[sidei];
2554  // fixed DoFs aren't test functions
2555  if (dof_is_fixed[sidei])
2556  continue;
2557  for (unsigned int sidej=0, freej=0;
2558  sidej != n_involved_dofs; ++sidej)
2559  {
2560  unsigned int j = involved_dofs[sidej];
2561  if (dof_is_fixed[sidej])
2562  Fe(freei) -= phi[i][qp] * phi[j][qp] *
2563  JxW[qp] * Uinvolved(sidej);
2564  else
2565  Ke(freei,freej) += phi[i][qp] *
2566  phi[j][qp] * JxW[qp];
2567  if (cont == C_ONE)
2568  {
2569  if (dof_is_fixed[sidej])
2570  Fe(freei) -= ( (*dphi)[i][qp] *
2571  (*dphi)[j][qp] ) *
2572  JxW[qp] * Uinvolved(sidej);
2573  else
2574  Ke(freei,freej) += ( (*dphi)[i][qp] *
2575  (*dphi)[j][qp] )
2576  * JxW[qp];
2577  }
2578  if (!dof_is_fixed[sidej])
2579  freej++;
2580  }
2581  Fe(freei) += phi[i][qp] * fineval * JxW[qp];
2582  if (cont == C_ONE)
2583  Fe(freei) += (finegrad * (*dphi)[i][qp] ) *
2584  JxW[qp];
2585  freei++;
2586  }
2587  }
2588 
2589  Ke.cholesky_solve(Fe, Ufree);
2590 
2591  // Transfer new edge solutions to element
2592  const processor_id_type pid = node ?
2593  node->processor_id() : DofObject::invalid_processor_id;
2594  insert_ids(free_dof_ids, Ufree.get_values(), pid);
2595 }
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
Definition: system.h:1545
uint8_t processor_id_type
Definition: id_types.h:104
static const processor_id_type invalid_processor_id
An invalid processor_id to distinguish DoFs that have not been assigned to a processor.
Definition: dof_object.h:404
std::unordered_map< dof_id_type, FValue > ids_to_save
FEContinuity
defines an enum for finite element types to libmesh_assert a certain level (or type? Hcurl?) of continuity.
void insert_ids(const std::vector< dof_id_type > &ids, const std::vector< FValue > &vals, processor_id_type pid)
uint8_t dof_id_type
Definition: id_types.h:67

◆ find_dofs_to_send()

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
void libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::find_dofs_to_send ( const Node node,
const Elem elem,
unsigned short  node_num,
const var_set vars 
)
inherited

Definition at line 2322 of file generic_projector.h.

References libMesh::Elem::active(), libMesh::Variable::active_on_subdomain(), libMesh::DofMap::dof_indices(), libMesh::MeshBase::elem_ref(), libMesh::DofObject::id(), libMesh::DofObject::invalid_processor_id, libMesh::Parallel::Utils::is_sorted(), libMesh::Elem::node_ptr(), libMesh::DofObject::processor_id(), and libMesh::Elem::subdomain_id().

2323 {
2324  libmesh_assert (&node == elem.node_ptr(node_num));
2325 
2326  // Any ghosted node in our range might have an owner who needs our
2327  // data
2328  const processor_id_type owner = node.processor_id();
2329  if (owner != system.processor_id())
2330  {
2331  const MeshBase & mesh = system.get_mesh();
2332  const DofMap & dof_map = system.get_dof_map();
2333 
2334  // But let's check and see if we can be certain the owner can
2335  // compute any or all of its own dof coefficients on that node.
2336  std::vector<dof_id_type> node_dof_ids, patch_dof_ids;
2337  for (const auto & var : vars)
2338  {
2339  const Variable & variable = system.variable(var);
2340 
2341  if (!variable.active_on_subdomain(elem.subdomain_id()))
2342  continue;
2343 
2344  dof_map.dof_indices(elem, node_num, node_dof_ids, var);
2345  }
2346  libmesh_assert(std::is_sorted(node_dof_ids.begin(),
2347  node_dof_ids.end()));
2348  const std::vector<dof_id_type> & patch =
2349  (*this->projector.nodes_to_elem)[node.id()];
2350  for (const auto & elem_id : patch)
2351  {
2352  const Elem & patch_elem = mesh.elem_ref(elem_id);
2353  if (!patch_elem.active() || owner != patch_elem.processor_id())
2354  continue;
2355  dof_map.dof_indices(&patch_elem, patch_dof_ids);
2356  std::sort(patch_dof_ids.begin(), patch_dof_ids.end());
2357 
2358  // Remove any node_dof_ids that we see can be calculated on
2359  // this element
2360  std::vector<dof_id_type> diff_ids(node_dof_ids.size());
2361  auto it = std::set_difference(node_dof_ids.begin(), node_dof_ids.end(),
2362  patch_dof_ids.begin(), patch_dof_ids.end(), diff_ids.begin());
2363  diff_ids.resize(it-diff_ids.begin());
2364  node_dof_ids.swap(diff_ids);
2365  if (node_dof_ids.empty())
2366  break;
2367  }
2368 
2369  // Give ids_to_push default invalid pid: not yet computed
2370  for (auto id : node_dof_ids)
2372  }
2373 }
const Variable & variable(unsigned int var) const
Return a constant reference to Variable var.
Definition: system.h:2167
const MeshBase & get_mesh() const
Definition: system.h:2067
uint8_t processor_id_type
Definition: id_types.h:104
bool is_sorted(const std::vector< KeyType > &v)
static const processor_id_type invalid_processor_id
An invalid processor_id to distinguish DoFs that have not been assigned to a processor.
Definition: dof_object.h:404
std::unordered_map< dof_id_type, std::vector< dof_id_type > > * nodes_to_elem
std::unordered_map< dof_id_type, std::pair< FValue, processor_id_type > > new_ids_to_push
processor_id_type processor_id() const
const DofMap & get_dof_map() const
Definition: system.h:2083

◆ insert_id() [1/2]

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
void libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::insert_id ( dof_id_type  id,
const FValue &  val,
processor_id_type  pid 
)
inherited

Definition at line 1112 of file generic_projector.h.

References libMesh::DofObject::invalid_processor_id.

1113 {
1114  auto iter = new_ids_to_push.find(id);
1115  if (iter == new_ids_to_push.end())
1116  action.insert(id, val);
1117  else
1118  {
1119  libmesh_assert(pid != DofObject::invalid_processor_id);
1120  iter->second = std::make_pair(val, pid);
1121  }
1122  if (!this->projector.done_saving_ids)
1123  {
1124  libmesh_assert(!new_ids_to_save.count(id));
1125  new_ids_to_save[id] = val;
1126  }
1127 }
std::unordered_map< dof_id_type, FValue > new_ids_to_save
static const processor_id_type invalid_processor_id
An invalid processor_id to distinguish DoFs that have not been assigned to a processor.
Definition: dof_object.h:404
std::unordered_map< dof_id_type, std::pair< FValue, processor_id_type > > new_ids_to_push

◆ insert_id() [2/2]

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
void libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::insert_id

Definition at line 1112 of file generic_projector.h.

1113 {
1114  auto iter = new_ids_to_push.find(id);
1115  if (iter == new_ids_to_push.end())
1116  action.insert(id, val);
1117  else
1118  {
1119  libmesh_assert(pid != DofObject::invalid_processor_id);
1120  iter->second = std::make_pair(val, pid);
1121  }
1122  if (!this->projector.done_saving_ids)
1123  {
1124  libmesh_assert(!new_ids_to_save.count(id));
1125  new_ids_to_save[id] = val;
1126  }
1127 }
std::unordered_map< dof_id_type, FValue > new_ids_to_save
static const processor_id_type invalid_processor_id
An invalid processor_id to distinguish DoFs that have not been assigned to a processor.
Definition: dof_object.h:404
std::unordered_map< dof_id_type, std::pair< FValue, processor_id_type > > new_ids_to_push

◆ insert_ids() [1/2]

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
void libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::insert_ids ( const std::vector< dof_id_type > &  ids,
const std::vector< FValue > &  vals,
processor_id_type  pid 
)
inherited

Definition at line 1133 of file generic_projector.h.

References libMesh::index_range(), and libMesh::DofObject::invalid_processor_id.

1134 {
1135  libmesh_assert_equal_to(ids.size(), vals.size());
1136  for (auto i : index_range(ids))
1137  {
1138  const dof_id_type id = ids[i];
1139  const FValue & val = vals[i];
1140 
1141  auto iter = new_ids_to_push.find(id);
1142  if (iter == new_ids_to_push.end())
1143  action.insert(id, val);
1144  else
1145  {
1146  libmesh_assert(pid != DofObject::invalid_processor_id);
1147  iter->second = std::make_pair(val, pid);
1148  }
1149  if (!this->projector.done_saving_ids)
1150  {
1151  libmesh_assert(!new_ids_to_save.count(id));
1152  new_ids_to_save[id] = val;
1153  }
1154  }
1155 }
IntRange< std::size_t > index_range(const std::vector< T > &vec)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:104
std::unordered_map< dof_id_type, FValue > new_ids_to_save
static const processor_id_type invalid_processor_id
An invalid processor_id to distinguish DoFs that have not been assigned to a processor.
Definition: dof_object.h:404
std::unordered_map< dof_id_type, std::pair< FValue, processor_id_type > > new_ids_to_push
uint8_t dof_id_type
Definition: id_types.h:67

◆ insert_ids() [2/2]

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
void libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::insert_ids

Definition at line 1133 of file generic_projector.h.

1134 {
1135  libmesh_assert_equal_to(ids.size(), vals.size());
1136  for (auto i : index_range(ids))
1137  {
1138  const dof_id_type id = ids[i];
1139  const FValue & val = vals[i];
1140 
1141  auto iter = new_ids_to_push.find(id);
1142  if (iter == new_ids_to_push.end())
1143  action.insert(id, val);
1144  else
1145  {
1146  libmesh_assert(pid != DofObject::invalid_processor_id);
1147  iter->second = std::make_pair(val, pid);
1148  }
1149  if (!this->projector.done_saving_ids)
1150  {
1151  libmesh_assert(!new_ids_to_save.count(id));
1152  new_ids_to_save[id] = val;
1153  }
1154  }
1155 }
IntRange< std::size_t > index_range(const std::vector< T > &vec)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:104
std::unordered_map< dof_id_type, FValue > new_ids_to_save
static const processor_id_type invalid_processor_id
An invalid processor_id to distinguish DoFs that have not been assigned to a processor.
Definition: dof_object.h:404
std::unordered_map< dof_id_type, std::pair< FValue, processor_id_type > > new_ids_to_push
uint8_t dof_id_type
Definition: id_types.h:67

◆ join()

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
void libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::join ( const SubFunctor other)
inherited

Definition at line 1161 of file generic_projector.h.

References libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::new_ids_to_push, and libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::new_ids_to_save.

1162 {
1163  new_ids_to_push.insert(other.new_ids_to_push.begin(), other.new_ids_to_push.end());
1164  new_ids_to_save.insert(other.new_ids_to_save.begin(), other.new_ids_to_save.end());
1165 }
std::unordered_map< dof_id_type, FValue > new_ids_to_save
std::unordered_map< dof_id_type, std::pair< FValue, processor_id_type > > new_ids_to_push

◆ operator()()

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
void libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectVertices::operator() ( const node_range range)

Definition at line 1637 of file generic_projector.h.

References libMesh::C_ONE, libMesh::C_ZERO, libMesh::Elem::dim(), libMesh::DISCONTINUOUS, libMesh::DofObject::dof_number(), libMesh::FEInterface::extra_hanging_dofs(), libMesh::FEType::family, libMesh::Elem::get_node_index(), libMesh::HERMITE, libMesh::Elem::hmin(), libMesh::DofObject::id(), libMesh::Elem::is_vertex_on_parent(), libMesh::DofObject::n_comp(), libMesh::FEInterface::n_dofs_at_node(), libMesh::Elem::parent(), libMesh::Elem::point(), libMesh::DofObject::processor_id(), libMesh::Real, libMesh::SCALAR, libMesh::TOLERANCE, libMesh::Variable::type(), libMesh::Elem::type(), and libMesh::Elem::which_child_am_i().

1638 {
1639  LOG_SCOPE ("project_vertices","GenericProjector");
1640 
1641  const unsigned int sys_num = system.number();
1642 
1643  // Variables with extra hanging dofs can't safely use eval_at_node
1644  // in as many places as variables without can.
1645  std::vector<unsigned short> extra_hanging_dofs;
1646  for (auto v_num : this->projector.variables)
1647  {
1648  if (extra_hanging_dofs.size() <= v_num)
1649  extra_hanging_dofs.resize(v_num+1, false);
1650  extra_hanging_dofs[v_num] =
1652  }
1653 
1654  for (const auto & v_pair : range)
1655  {
1656  const Node & vertex = *v_pair.first;
1657  const Elem & elem = *std::get<0>(v_pair.second);
1658  const unsigned int n = std::get<1>(v_pair.second);
1659  const var_set & vertex_vars = std::get<2>(v_pair.second);
1660 
1661  context.pre_fe_reinit(system, &elem);
1662 
1663  this->find_dofs_to_send(vertex, elem, n, vertex_vars);
1664 
1665  // Look at all the variables we're supposed to interpolate from
1666  // this element on this vertex
1667  for (const auto & var : vertex_vars)
1668  {
1669  const Variable & variable = system.variable(var);
1670  const FEType & base_fe_type = variable.type();
1671  const unsigned int var_component =
1673 
1674  if (base_fe_type.family == SCALAR)
1675  continue;
1676 
1677  const FEContinuity & cont = this->conts[var];
1678  if (cont == DISCONTINUOUS)
1679  {
1680  libmesh_assert_equal_to(vertex.n_comp(sys_num, var), 0);
1681  }
1682  else if (cont == C_ZERO)
1683  {
1684  libmesh_assert(vertex.n_comp(sys_num, var));
1685  const dof_id_type id = vertex.dof_number(sys_num, var, 0);
1686  // C_ZERO elements have a single nodal value DoF at vertices
1687  const FValue val = f.eval_at_node
1688  (context, var_component, /*dim=*/ 0, // Don't care w/C0
1689  vertex, extra_hanging_dofs[var], system.time);
1690  insert_id(id, val, vertex.processor_id());
1691  }
1692  else if (cont == C_ONE)
1693  {
1694  libmesh_assert(vertex.n_comp(sys_num, var));
1695  const dof_id_type first_id = vertex.dof_number(sys_num, var, 0);
1696 
1697  // C_ONE elements have a single nodal value and dim
1698  // gradient values at vertices, as well as cross
1699  // gradients for HERMITE. We need to have an element in
1700  // hand to figure out dim and to have in case this
1701  // vertex is a new vertex.
1702  const int dim = elem.dim();
1703 #ifndef NDEBUG
1704  // For now all C1 elements at a vertex had better have
1705  // the same dimension. If anyone hits these asserts let
1706  // me know; we could probably support a mixed-dimension
1707  // mesh IFF the 2D elements were all parallel to xy and
1708  // the 1D elements all parallel to x.
1709  for (const auto e_id : (*this->projector.nodes_to_elem)[vertex.id()])
1710  {
1711  const Elem & e = system.get_mesh().elem_ref(e_id);
1712  libmesh_assert_equal_to(dim, e.dim());
1713  }
1714 #endif
1715 #ifdef LIBMESH_ENABLE_AMR
1716  bool is_parent_vertex = false;
1717  if (elem.parent())
1718  {
1719  const int i_am_child =
1720  elem.parent()->which_child_am_i(&elem);
1721  is_parent_vertex =
1722  elem.parent()->is_vertex_on_parent(i_am_child, n);
1723  }
1724 #else
1725  const bool is_parent_vertex = false;
1726 #endif
1727 
1728  // The hermite element vertex shape functions are weird
1729  if (base_fe_type.family == HERMITE)
1730  {
1731  const FValue val =
1732  f.eval_at_node(context,
1733  var_component,
1734  dim,
1735  vertex,
1736  extra_hanging_dofs[var],
1737  system.time);
1738  insert_id(first_id, val, vertex.processor_id());
1739 
1740  VectorValue<FValue> grad =
1741  is_parent_vertex ?
1742  g->eval_at_node(context,
1743  var_component,
1744  dim,
1745  vertex,
1746  extra_hanging_dofs[var],
1747  system.time) :
1748  g->eval_at_point(context,
1749  var_component,
1750  vertex,
1751  system.time);
1752  // x derivative
1753  insert_id(first_id+1, grad(0),
1754  vertex.processor_id());
1755  if (dim > 1)
1756  {
1757  // We'll finite difference mixed derivatives
1758  Real delta_x = TOLERANCE * elem.hmin();
1759 
1760  Point nxminus = elem.point(n),
1761  nxplus = elem.point(n);
1762  nxminus(0) -= delta_x;
1763  nxplus(0) += delta_x;
1764  VectorValue<FValue> gxminus =
1765  g->eval_at_point(context,
1766  var_component,
1767  nxminus,
1768  system.time);
1769  VectorValue<FValue> gxplus =
1770  g->eval_at_point(context,
1771  var_component,
1772  nxplus,
1773  system.time);
1774  // y derivative
1775  insert_id(first_id+2, grad(1),
1776  vertex.processor_id());
1777  // xy derivative
1778  insert_id(first_id+3,
1779  (gxplus(1) - gxminus(1)) / 2. / delta_x,
1780  vertex.processor_id());
1781 
1782  if (dim > 2)
1783  {
1784  // z derivative
1785  insert_id(first_id+4, grad(2),
1786  vertex.processor_id());
1787  // xz derivative
1788  insert_id(first_id+5,
1789  (gxplus(2) - gxminus(2)) / 2. / delta_x,
1790  vertex.processor_id());
1791 
1792  // We need new points for yz
1793  Point nyminus = elem.point(n),
1794  nyplus = elem.point(n);
1795  nyminus(1) -= delta_x;
1796  nyplus(1) += delta_x;
1797  VectorValue<FValue> gyminus =
1798  g->eval_at_point(context,
1799  var_component,
1800  nyminus,
1801  system.time);
1802  VectorValue<FValue> gyplus =
1803  g->eval_at_point(context,
1804  var_component,
1805  nyplus,
1806  system.time);
1807  // yz derivative
1808  insert_id(first_id+6,
1809  (gyplus(2) - gyminus(2)) / 2. / delta_x,
1810  vertex.processor_id());
1811  // Getting a 2nd order xyz is more tedious
1812  Point nxmym = elem.point(n),
1813  nxmyp = elem.point(n),
1814  nxpym = elem.point(n),
1815  nxpyp = elem.point(n);
1816  nxmym(0) -= delta_x;
1817  nxmym(1) -= delta_x;
1818  nxmyp(0) -= delta_x;
1819  nxmyp(1) += delta_x;
1820  nxpym(0) += delta_x;
1821  nxpym(1) -= delta_x;
1822  nxpyp(0) += delta_x;
1823  nxpyp(1) += delta_x;
1824  VectorValue<FValue> gxmym =
1825  g->eval_at_point(context,
1826  var_component,
1827  nxmym,
1828  system.time);
1829  VectorValue<FValue> gxmyp =
1830  g->eval_at_point(context,
1831  var_component,
1832  nxmyp,
1833  system.time);
1834  VectorValue<FValue> gxpym =
1835  g->eval_at_point(context,
1836  var_component,
1837  nxpym,
1838  system.time);
1839  VectorValue<FValue> gxpyp =
1840  g->eval_at_point(context,
1841  var_component,
1842  nxpyp,
1843  system.time);
1844  FValue gxzplus = (gxpyp(2) - gxmyp(2))
1845  / 2. / delta_x;
1846  FValue gxzminus = (gxpym(2) - gxmym(2))
1847  / 2. / delta_x;
1848  // xyz derivative
1849  insert_id(first_id+7,
1850  (gxzplus - gxzminus) / 2. / delta_x,
1851  vertex.processor_id());
1852  }
1853  }
1854  }
1855  else
1856  {
1857  // Currently other C_ONE elements have a single nodal
1858  // value shape function and nodal gradient component
1859  // shape functions
1860  libmesh_assert_equal_to
1862  (dim, base_fe_type, elem.type(),
1863  elem.get_node_index(&vertex)),
1864  (unsigned int)(1 + dim));
1865  const FValue val =
1866  f.eval_at_node(context, var_component, dim,
1867  vertex, extra_hanging_dofs[var],
1868  system.time);
1869  insert_id(first_id, val, vertex.processor_id());
1870  VectorValue<FValue> grad =
1871  is_parent_vertex ?
1872  g->eval_at_node(context, var_component, dim,
1873  vertex, extra_hanging_dofs[var],
1874  system.time) :
1875  g->eval_at_point(context, var_component, vertex,
1876  system.time);
1877  for (int i=0; i!= dim; ++i)
1878  insert_id(first_id + i + 1, grad(i),
1879  vertex.processor_id());
1880  }
1881  }
1882  else
1883  libmesh_error_msg("Unknown continuity " << cont);
1884  }
1885  }
1886 }
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
Definition: system.h:1545
const Elem * parent() const
Definition: elem.h:2453
unsigned int variable_scalar_number(const std::string &var, unsigned int component) const
Definition: system.h:2198
const Variable & variable(unsigned int var) const
Return a constant reference to Variable var.
Definition: system.h:2167
void insert_id(dof_id_type id, const FValue &val, processor_id_type pid)
static const Real TOLERANCE
const MeshBase & get_mesh() const
Definition: system.h:2067
unsigned int number() const
Definition: system.h:2059
const std::vector< unsigned int > & variables
unsigned int which_child_am_i(const Elem *e) const
Definition: elem.h:2581
void find_dofs_to_send(const Node &node, const Elem &elem, unsigned short node_num, const var_set &vars)
std::unordered_map< dof_id_type, std::vector< dof_id_type > > * nodes_to_elem
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
FEContinuity
defines an enum for finite element types to libmesh_assert a certain level (or type? Hcurl?) of continuity.
virtual const Elem & elem_ref(const dof_id_type i) const
Definition: mesh_base.h:520
std::set< unsigned int > var_set
static unsigned int n_dofs_at_node(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
virtual void pre_fe_reinit(const System &, const Elem *e)
Reinitializes local data vectors/matrices on the current geometric element.
std::vector< FEContinuity > conts
uint8_t dof_id_type
Definition: id_types.h:67
static bool extra_hanging_dofs(const FEType &fe_t)
const FEType & type() const
Definition: variable.h:119

Member Data Documentation

◆ action

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
ProjectionAction libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::action

Definition at line 188 of file generic_projector.h.

◆ context

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
FEMContext libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::context

Definition at line 193 of file generic_projector.h.

◆ conts

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
std::vector<FEContinuity> libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::conts

Definition at line 196 of file generic_projector.h.

◆ f

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
FFunctor libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::f

Definition at line 189 of file generic_projector.h.

◆ g

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
std::unique_ptr<GFunctor> libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubProjector::g

Definition at line 222 of file generic_projector.h.

◆ new_ids_to_push

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
std::unordered_map<dof_id_type, std::pair<FValue, processor_id_type> > libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::new_ids_to_push
inherited

◆ new_ids_to_save

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
std::unordered_map<dof_id_type, FValue> libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::new_ids_to_save
inherited

◆ projector

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
GenericProjector& libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::projector
inherited

Definition at line 149 of file generic_projector.h.

◆ system

template<typename FFunctor , typename GFunctor , typename FValue , typename ProjectionAction >
const System& libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::system

Definition at line 198 of file generic_projector.h.


The documentation for this struct was generated from the following file: