libMesh::BoundaryProjectSolution Class Reference

Public Member Functions

 BoundaryProjectSolution (const std::set< boundary_id_type > &b_in, const std::vector< unsigned int > &variables_in, const System &system_in, FunctionBase< Number > *f_in, FunctionBase< Gradient > *g_in, const Parameters &parameters_in, NumericVector< Number > &new_v_in)
 
 BoundaryProjectSolution (const BoundaryProjectSolution &in)
 
void operator() (const ConstElemRange &range) const
 

Private Attributes

const std::set< boundary_id_type > & b
 
const std::vector< unsigned int > & variables
 
const Systemsystem
 
UniquePtr< FunctionBase< Number > > f
 
UniquePtr< FunctionBase< Gradient > > g
 
const Parametersparameters
 
NumericVector< Number > & new_vector
 

Detailed Description

This class implements projecting an arbitrary boundary function to the current mesh. This may be executed in parallel on multiple threads.

Definition at line 521 of file system_projection.C.

Constructor & Destructor Documentation

libMesh::BoundaryProjectSolution::BoundaryProjectSolution ( const std::set< boundary_id_type > &  b_in,
const std::vector< unsigned int > &  variables_in,
const System system_in,
FunctionBase< Number > *  f_in,
FunctionBase< Gradient > *  g_in,
const Parameters parameters_in,
NumericVector< Number > &  new_v_in 
)
inline

Definition at line 533 of file system_projection.C.

References libMesh::libmesh_assert().

539  :
540  b(b_in),
541  variables(variables_in),
542  system(system_in),
543  f(f_in ? f_in->clone() : UniquePtr<FunctionBase<Number> >()),
544  g(g_in ? g_in->clone() : UniquePtr<FunctionBase<Gradient> >()),
545  parameters(parameters_in),
546  new_vector(new_v_in)
547  {
548  libmesh_assert(f.get());
549  f->init();
550  if (g.get())
551  g->init();
552  }
UniquePtr< FunctionBase< Number > > f
UniquePtr< FunctionBase< Gradient > > g
const std::vector< unsigned int > & variables
libmesh_assert(j)
NumericVector< Number > & new_vector
virtual UniquePtr< FunctionBase< Output > > clone() const =0
const std::set< boundary_id_type > & b
libMesh::BoundaryProjectSolution::BoundaryProjectSolution ( const BoundaryProjectSolution in)
inline

Definition at line 554 of file system_projection.C.

References libMesh::libmesh_assert(), and libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::operator()().

554  :
555  b(in.b),
556  variables(in.variables),
557  system(in.system),
558  f(in.f.get() ? in.f->clone() : UniquePtr<FunctionBase<Number> >()),
559  g(in.g.get() ? in.g->clone() : UniquePtr<FunctionBase<Gradient> >()),
560  parameters(in.parameters),
561  new_vector(in.new_vector)
562  {
563  libmesh_assert(f.get());
564  f->init();
565  if (g.get())
566  g->init();
567  }
UniquePtr< FunctionBase< Number > > f
UniquePtr< FunctionBase< Gradient > > g
const std::vector< unsigned int > & variables
libmesh_assert(j)
NumericVector< Number > & new_vector
const std::set< boundary_id_type > & b

Member Function Documentation

void libMesh::BoundaryProjectSolution::operator() ( const ConstElemRange range) const

This method projects an arbitrary boundary solution to the current mesh. The input function f gives the arbitrary solution, while the new_vector (which should already be correctly sized) gives the solution (to be computed) on the current mesh.

Definition at line 2128 of file system_projection.C.

References std::abs(), libMesh::Variable::active_on_subdomain(), libMesh::StoredRange< iterator_type, object_type >::begin(), libMesh::BoundaryInfo::boundary_ids(), libMesh::FEGenericBase< OutputType >::build(), libMesh::C_ONE, libMesh::C_ZERO, libMesh::DenseMatrix< T >::cholesky_solve(), libMesh::FEType::default_quadrature_rule(), libMesh::DISCONTINUOUS, libMesh::DofMap::dof_indices(), libMesh::FEInterface::dofs_on_edge(), libMesh::FEInterface::dofs_on_side(), libMesh::StoredRange< iterator_type, object_type >::end(), libMesh::FEType::family, libMesh::MeshBase::get_boundary_info(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::HERMITE, libMesh::Elem::is_edge_on_side(), libMesh::Elem::is_node_on_side(), libMesh::Elem::is_vertex(), libMesh::libmesh_assert(), libmesh_nullptr, libMesh::MeshBase::mesh_dimension(), libMesh::FEInterface::n_dofs_at_node(), libMesh::Elem::n_edges(), n_nodes, libMesh::Elem::n_nodes(), libMesh::Elem::n_sides(), libMesh::Elem::point(), libMesh::DenseVector< T >::resize(), libMesh::DenseMatrix< T >::resize(), libMesh::SCALAR, libMesh::Threads::spin_mtx, libMesh::Elem::subdomain_id(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::system, libMesh::System::time, libMesh::TOLERANCE, libMesh::Variable::type(), libMesh::Elem::type(), libMesh::DofMap::variable(), libMesh::System::variable_scalar_number(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::variables, libMesh::DenseMatrix< T >::zero(), and libMesh::DenseVector< T >::zero().

2129 {
2130  // We need data to project
2131  libmesh_assert(f.get());
2132 
2140  // The dimensionality of the current mesh
2141  const unsigned int dim = system.get_mesh().mesh_dimension();
2142 
2143  // The DofMap for this system
2144  const DofMap & dof_map = system.get_dof_map();
2145 
2146  // Boundary info for the current mesh
2147  const BoundaryInfo & boundary_info =
2149 
2150  // The element matrix and RHS for projections.
2151  // Note that Ke is always real-valued, whereas
2152  // Fe may be complex valued if complex number
2153  // support is enabled
2154  DenseMatrix<Real> Ke;
2155  DenseVector<Number> Fe;
2156  // The new element coefficients
2157  DenseVector<Number> Ue;
2158 
2159 
2160  // Loop over all the variables we've been requested to project
2161  for (std::size_t v=0; v!=variables.size(); v++)
2162  {
2163  const unsigned int var = variables[v];
2164 
2165  const Variable & variable = dof_map.variable(var);
2166 
2167  const FEType & fe_type = variable.type();
2168 
2169  if (fe_type.family == SCALAR)
2170  continue;
2171 
2172  const unsigned int var_component =
2174 
2175  // Get FE objects of the appropriate type
2176  UniquePtr<FEBase> fe (FEBase::build(dim, fe_type));
2177 
2178  // Prepare variables for projection
2179  UniquePtr<QBase> qedgerule (fe_type.default_quadrature_rule(1));
2180  UniquePtr<QBase> qsiderule (fe_type.default_quadrature_rule(dim-1));
2181 
2182  // The values of the shape functions at the quadrature
2183  // points
2184  const std::vector<std::vector<Real> > & phi = fe->get_phi();
2185 
2186  // The gradients of the shape functions at the quadrature
2187  // points on the child element.
2188  const std::vector<std::vector<RealGradient> > * dphi = libmesh_nullptr;
2189 
2190  const FEContinuity cont = fe->get_continuity();
2191 
2192  if (cont == C_ONE)
2193  {
2194  // We'll need gradient data for a C1 projection
2195  libmesh_assert(g.get());
2196 
2197  const std::vector<std::vector<RealGradient> > &
2198  ref_dphi = fe->get_dphi();
2199  dphi = &ref_dphi;
2200  }
2201 
2202  // The Jacobian * quadrature weight at the quadrature points
2203  const std::vector<Real> & JxW =
2204  fe->get_JxW();
2205 
2206  // The XYZ locations of the quadrature points
2207  const std::vector<Point> & xyz_values =
2208  fe->get_xyz();
2209 
2210  // The global DOF indices
2211  std::vector<dof_id_type> dof_indices;
2212  // Side/edge DOF indices
2213  std::vector<unsigned int> side_dofs;
2214 
2215  // Container to catch IDs passed back from BoundaryInfo.
2216  std::vector<boundary_id_type> bc_ids;
2217 
2218  // Iterate over all the elements in the range
2219  for (ConstElemRange::const_iterator elem_it=range.begin(); elem_it != range.end(); ++elem_it)
2220  {
2221  const Elem * elem = *elem_it;
2222 
2223  // Per-subdomain variables don't need to be projected on
2224  // elements where they're not active
2225  if (!variable.active_on_subdomain(elem->subdomain_id()))
2226  continue;
2227 
2228  // Find out which nodes, edges and sides are on a requested
2229  // boundary:
2230  std::vector<bool> is_boundary_node(elem->n_nodes(), false),
2231  is_boundary_edge(elem->n_edges(), false),
2232  is_boundary_side(elem->n_sides(), false);
2233  for (unsigned char s=0; s != elem->n_sides(); ++s)
2234  {
2235  // First see if this side has been requested
2236  boundary_info.boundary_ids (elem, s, bc_ids);
2237  bool do_this_side = false;
2238  for (std::vector<boundary_id_type>::iterator i=bc_ids.begin();
2239  i!=bc_ids.end(); ++i)
2240  if (b.count(*i))
2241  {
2242  do_this_side = true;
2243  break;
2244  }
2245  if (!do_this_side)
2246  continue;
2247 
2248  is_boundary_side[s] = true;
2249 
2250  // Then see what nodes and what edges are on it
2251  for (unsigned int n=0; n != elem->n_nodes(); ++n)
2252  if (elem->is_node_on_side(n,s))
2253  is_boundary_node[n] = true;
2254  for (unsigned int e=0; e != elem->n_edges(); ++e)
2255  if (elem->is_edge_on_side(e,s))
2256  is_boundary_edge[e] = true;
2257  }
2258 
2259  // Update the DOF indices for this element based on
2260  // the current mesh
2261  dof_map.dof_indices (elem, dof_indices, var);
2262 
2263  // The number of DOFs on the element
2264  const unsigned int n_dofs =
2265  cast_int<unsigned int>(dof_indices.size());
2266 
2267  // Fixed vs. free DoFs on edge/face projections
2268  std::vector<char> dof_is_fixed(n_dofs, false); // bools
2269  std::vector<int> free_dof(n_dofs, 0);
2270 
2271  // The element type
2272  const ElemType elem_type = elem->type();
2273 
2274  // The number of nodes on the new element
2275  const unsigned int n_nodes = elem->n_nodes();
2276 
2277  // Zero the interpolated values
2278  Ue.resize (n_dofs); Ue.zero();
2279 
2280  // In general, we need a series of
2281  // projections to ensure a unique and continuous
2282  // solution. We start by interpolating boundary nodes, then
2283  // hold those fixed and project boundary edges, then hold
2284  // those fixed and project boundary faces,
2285 
2286  // Interpolate node values first
2287  unsigned int current_dof = 0;
2288  for (unsigned int n=0; n!= n_nodes; ++n)
2289  {
2290  // FIXME: this should go through the DofMap,
2291  // not duplicate dof_indices code badly!
2292  const unsigned int nc =
2293  FEInterface::n_dofs_at_node (dim, fe_type, elem_type,
2294  n);
2295  if (!elem->is_vertex(n) || !is_boundary_node[n])
2296  {
2297  current_dof += nc;
2298  continue;
2299  }
2300  if (cont == DISCONTINUOUS)
2301  {
2302  libmesh_assert_equal_to (nc, 0);
2303  }
2304  // Assume that C_ZERO elements have a single nodal
2305  // value shape function
2306  else if (cont == C_ZERO)
2307  {
2308  libmesh_assert_equal_to (nc, 1);
2309  Ue(current_dof) = f->component(var_component,
2310  elem->point(n),
2311  system.time);
2312  dof_is_fixed[current_dof] = true;
2313  current_dof++;
2314  }
2315  // The hermite element vertex shape functions are weird
2316  else if (fe_type.family == HERMITE)
2317  {
2318  Ue(current_dof) = f->component(var_component,
2319  elem->point(n),
2320  system.time);
2321  dof_is_fixed[current_dof] = true;
2322  current_dof++;
2323  Gradient grad = g->component(var_component,
2324  elem->point(n),
2325  system.time);
2326  // x derivative
2327  Ue(current_dof) = grad(0);
2328  dof_is_fixed[current_dof] = true;
2329  current_dof++;
2330  if (dim > 1)
2331  {
2332  // We'll finite difference mixed derivatives
2333  Point nxminus = elem->point(n),
2334  nxplus = elem->point(n);
2335  nxminus(0) -= TOLERANCE;
2336  nxplus(0) += TOLERANCE;
2337  Gradient gxminus = g->component(var_component,
2338  nxminus,
2339  system.time);
2340  Gradient gxplus = g->component(var_component,
2341  nxplus,
2342  system.time);
2343  // y derivative
2344  Ue(current_dof) = grad(1);
2345  dof_is_fixed[current_dof] = true;
2346  current_dof++;
2347  // xy derivative
2348  Ue(current_dof) = (gxplus(1) - gxminus(1))
2349  / 2. / TOLERANCE;
2350  dof_is_fixed[current_dof] = true;
2351  current_dof++;
2352 
2353  if (dim > 2)
2354  {
2355  // z derivative
2356  Ue(current_dof) = grad(2);
2357  dof_is_fixed[current_dof] = true;
2358  current_dof++;
2359  // xz derivative
2360  Ue(current_dof) = (gxplus(2) - gxminus(2))
2361  / 2. / TOLERANCE;
2362  dof_is_fixed[current_dof] = true;
2363  current_dof++;
2364  // We need new points for yz
2365  Point nyminus = elem->point(n),
2366  nyplus = elem->point(n);
2367  nyminus(1) -= TOLERANCE;
2368  nyplus(1) += TOLERANCE;
2369  Gradient gyminus = g->component(var_component,
2370  nyminus,
2371  system.time);
2372  Gradient gyplus = g->component(var_component,
2373  nyplus,
2374  system.time);
2375  // xz derivative
2376  Ue(current_dof) = (gyplus(2) - gyminus(2))
2377  / 2. / TOLERANCE;
2378  dof_is_fixed[current_dof] = true;
2379  current_dof++;
2380  // Getting a 2nd order xyz is more tedious
2381  Point nxmym = elem->point(n),
2382  nxmyp = elem->point(n),
2383  nxpym = elem->point(n),
2384  nxpyp = elem->point(n);
2385  nxmym(0) -= TOLERANCE;
2386  nxmym(1) -= TOLERANCE;
2387  nxmyp(0) -= TOLERANCE;
2388  nxmyp(1) += TOLERANCE;
2389  nxpym(0) += TOLERANCE;
2390  nxpym(1) -= TOLERANCE;
2391  nxpyp(0) += TOLERANCE;
2392  nxpyp(1) += TOLERANCE;
2393  Gradient gxmym = g->component(var_component,
2394  nxmym,
2395  system.time);
2396  Gradient gxmyp = g->component(var_component,
2397  nxmyp,
2398  system.time);
2399  Gradient gxpym = g->component(var_component,
2400  nxpym,
2401  system.time);
2402  Gradient gxpyp = g->component(var_component,
2403  nxpyp,
2404  system.time);
2405  Number gxzplus = (gxpyp(2) - gxmyp(2))
2406  / 2. / TOLERANCE;
2407  Number gxzminus = (gxpym(2) - gxmym(2))
2408  / 2. / TOLERANCE;
2409  // xyz derivative
2410  Ue(current_dof) = (gxzplus - gxzminus)
2411  / 2. / TOLERANCE;
2412  dof_is_fixed[current_dof] = true;
2413  current_dof++;
2414  }
2415  }
2416  }
2417  // Assume that other C_ONE elements have a single nodal
2418  // value shape function and nodal gradient component
2419  // shape functions
2420  else if (cont == C_ONE)
2421  {
2422  libmesh_assert_equal_to (nc, 1 + dim);
2423  Ue(current_dof) = f->component(var_component,
2424  elem->point(n),
2425  system.time);
2426  dof_is_fixed[current_dof] = true;
2427  current_dof++;
2428  Gradient grad = g->component(var_component,
2429  elem->point(n),
2430  system.time);
2431  for (unsigned int i=0; i!= dim; ++i)
2432  {
2433  Ue(current_dof) = grad(i);
2434  dof_is_fixed[current_dof] = true;
2435  current_dof++;
2436  }
2437  }
2438  else
2439  libmesh_error_msg("Unknown continuity " << cont);
2440  }
2441 
2442  // In 3D, project any edge values next
2443  if (dim > 2 && cont != DISCONTINUOUS)
2444  for (unsigned int e=0; e != elem->n_edges(); ++e)
2445  {
2446  if (!is_boundary_edge[e])
2447  continue;
2448 
2449  FEInterface::dofs_on_edge(elem, dim, fe_type, e,
2450  side_dofs);
2451 
2452  // Some edge dofs are on nodes and already
2453  // fixed, others are free to calculate
2454  unsigned int free_dofs = 0;
2455  for (std::size_t i=0; i != side_dofs.size(); ++i)
2456  if (!dof_is_fixed[side_dofs[i]])
2457  free_dof[free_dofs++] = i;
2458 
2459  // There may be nothing to project
2460  if (!free_dofs)
2461  continue;
2462 
2463  Ke.resize (free_dofs, free_dofs); Ke.zero();
2464  Fe.resize (free_dofs); Fe.zero();
2465  // The new edge coefficients
2466  DenseVector<Number> Uedge(free_dofs);
2467 
2468  // Initialize FE data on the edge
2469  fe->attach_quadrature_rule (qedgerule.get());
2470  fe->edge_reinit (elem, e);
2471  const unsigned int n_qp = qedgerule->n_points();
2472 
2473  // Loop over the quadrature points
2474  for (unsigned int qp=0; qp<n_qp; qp++)
2475  {
2476  // solution at the quadrature point
2477  Number fineval = f->component(var_component,
2478  xyz_values[qp],
2479  system.time);
2480  // solution grad at the quadrature point
2481  Gradient finegrad;
2482  if (cont == C_ONE)
2483  finegrad = g->component(var_component,
2484  xyz_values[qp],
2485  system.time);
2486 
2487  // Form edge projection matrix
2488  for (std::size_t sidei=0, freei=0;
2489  sidei != side_dofs.size(); ++sidei)
2490  {
2491  unsigned int i = side_dofs[sidei];
2492  // fixed DoFs aren't test functions
2493  if (dof_is_fixed[i])
2494  continue;
2495  for (std::size_t sidej=0, freej=0;
2496  sidej != side_dofs.size(); ++sidej)
2497  {
2498  unsigned int j = side_dofs[sidej];
2499  if (dof_is_fixed[j])
2500  Fe(freei) -= phi[i][qp] * phi[j][qp] *
2501  JxW[qp] * Ue(j);
2502  else
2503  Ke(freei,freej) += phi[i][qp] *
2504  phi[j][qp] * JxW[qp];
2505  if (cont == C_ONE)
2506  {
2507  if (dof_is_fixed[j])
2508  Fe(freei) -= ((*dphi)[i][qp] *
2509  (*dphi)[j][qp]) *
2510  JxW[qp] * Ue(j);
2511  else
2512  Ke(freei,freej) += ((*dphi)[i][qp] *
2513  (*dphi)[j][qp])
2514  * JxW[qp];
2515  }
2516  if (!dof_is_fixed[j])
2517  freej++;
2518  }
2519  Fe(freei) += phi[i][qp] * fineval * JxW[qp];
2520  if (cont == C_ONE)
2521  Fe(freei) += (finegrad * (*dphi)[i][qp]) *
2522  JxW[qp];
2523  freei++;
2524  }
2525  }
2526 
2527  Ke.cholesky_solve(Fe, Uedge);
2528 
2529  // Transfer new edge solutions to element
2530  for (unsigned int i=0; i != free_dofs; ++i)
2531  {
2532  Number & ui = Ue(side_dofs[free_dof[i]]);
2534  std::abs(ui - Uedge(i)) < TOLERANCE);
2535  ui = Uedge(i);
2536  dof_is_fixed[side_dofs[free_dof[i]]] = true;
2537  }
2538  }
2539 
2540  // Project any side values (edges in 2D, faces in 3D)
2541  if (dim > 1 && cont != DISCONTINUOUS)
2542  for (unsigned int s=0; s != elem->n_sides(); ++s)
2543  {
2544  if (!is_boundary_side[s])
2545  continue;
2546 
2547  FEInterface::dofs_on_side(elem, dim, fe_type, s,
2548  side_dofs);
2549 
2550  // Some side dofs are on nodes/edges and already
2551  // fixed, others are free to calculate
2552  unsigned int free_dofs = 0;
2553  for (std::size_t i=0; i != side_dofs.size(); ++i)
2554  if (!dof_is_fixed[side_dofs[i]])
2555  free_dof[free_dofs++] = i;
2556 
2557  // There may be nothing to project
2558  if (!free_dofs)
2559  continue;
2560 
2561  Ke.resize (free_dofs, free_dofs); Ke.zero();
2562  Fe.resize (free_dofs); Fe.zero();
2563  // The new side coefficients
2564  DenseVector<Number> Uside(free_dofs);
2565 
2566  // Initialize FE data on the side
2567  fe->attach_quadrature_rule (qsiderule.get());
2568  fe->reinit (elem, s);
2569  const unsigned int n_qp = qsiderule->n_points();
2570 
2571  // Loop over the quadrature points
2572  for (unsigned int qp=0; qp<n_qp; qp++)
2573  {
2574  // solution at the quadrature point
2575  Number fineval = f->component(var_component,
2576  xyz_values[qp],
2577  system.time);
2578  // solution grad at the quadrature point
2579  Gradient finegrad;
2580  if (cont == C_ONE)
2581  finegrad = g->component(var_component,
2582  xyz_values[qp],
2583  system.time);
2584 
2585  // Form side projection matrix
2586  for (std::size_t sidei=0, freei=0;
2587  sidei != side_dofs.size(); ++sidei)
2588  {
2589  unsigned int i = side_dofs[sidei];
2590  // fixed DoFs aren't test functions
2591  if (dof_is_fixed[i])
2592  continue;
2593  for (std::size_t sidej=0, freej=0;
2594  sidej != side_dofs.size(); ++sidej)
2595  {
2596  unsigned int j = side_dofs[sidej];
2597  if (dof_is_fixed[j])
2598  Fe(freei) -= phi[i][qp] * phi[j][qp] *
2599  JxW[qp] * Ue(j);
2600  else
2601  Ke(freei,freej) += phi[i][qp] *
2602  phi[j][qp] * JxW[qp];
2603  if (cont == C_ONE)
2604  {
2605  if (dof_is_fixed[j])
2606  Fe(freei) -= ((*dphi)[i][qp] *
2607  (*dphi)[j][qp]) *
2608  JxW[qp] * Ue(j);
2609  else
2610  Ke(freei,freej) += ((*dphi)[i][qp] *
2611  (*dphi)[j][qp])
2612  * JxW[qp];
2613  }
2614  if (!dof_is_fixed[j])
2615  freej++;
2616  }
2617  Fe(freei) += (fineval * phi[i][qp]) * JxW[qp];
2618  if (cont == C_ONE)
2619  Fe(freei) += (finegrad * (*dphi)[i][qp]) *
2620  JxW[qp];
2621  freei++;
2622  }
2623  }
2624 
2625  Ke.cholesky_solve(Fe, Uside);
2626 
2627  // Transfer new side solutions to element
2628  for (unsigned int i=0; i != free_dofs; ++i)
2629  {
2630  Number & ui = Ue(side_dofs[free_dof[i]]);
2632  std::abs(ui - Uside(i)) < TOLERANCE);
2633  ui = Uside(i);
2634  dof_is_fixed[side_dofs[free_dof[i]]] = true;
2635  }
2636  }
2637 
2638  const dof_id_type
2639  first = new_vector.first_local_index(),
2640  last = new_vector.last_local_index();
2641 
2642  // Lock the new_vector since it is shared among threads.
2643  {
2644  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
2645 
2646  for (unsigned int i = 0; i < n_dofs; i++)
2647  if (dof_is_fixed[i] &&
2648  (dof_indices[i] >= first) &&
2649  (dof_indices[i] < last))
2650  {
2651  new_vector.set(dof_indices[i], Ue(i));
2652  }
2653  }
2654  } // end elem loop
2655  } // end variables loop
2656 }
static void dofs_on_edge(const Elem *const elem, const unsigned int dim, const FEType &fe_t, unsigned int e, std::vector< unsigned int > &di)
Definition: fe_interface.C:497
const BoundaryInfo & get_boundary_info() const
Definition: mesh_base.h:114
double abs(double a)
virtual numeric_index_type last_local_index() const =0
static void dofs_on_side(const Elem *const elem, const unsigned int dim, const FEType &fe_t, unsigned int s, std::vector< unsigned int > &di)
Definition: fe_interface.C:482
UniquePtr< FunctionBase< Number > > f
const class libmesh_nullptr_t libmesh_nullptr
static const Real TOLERANCE
UniquePtr< FunctionBase< Gradient > > g
const std::vector< unsigned int > & variables
libmesh_assert(j)
const dof_id_type n_nodes
Definition: tecplot_io.C:67
const MeshBase & get_mesh() const
Definition: system.h:2012
spin_mutex spin_mtx
Definition: threads.C:29
const DofMap & get_dof_map() const
Definition: system.h:2028
NumericVector< Number > & new_vector
NumberVectorValue Gradient
std::vector< object_type >::const_iterator const_iterator
Definition: stored_range.h:58
static unsigned int n_dofs_at_node(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
Definition: fe_interface.C:436
unsigned int mesh_dimension() const
Definition: mesh_base.C:147
virtual numeric_index_type first_local_index() const =0
unsigned int variable_scalar_number(const std::string &var, unsigned int component) const
Definition: system.h:2143
virtual void set(const numeric_index_type i, const T value)=0
const std::set< boundary_id_type > & b
uint8_t dof_id_type
Definition: id_types.h:64
static UniquePtr< FEGenericBase > build(const unsigned int dim, const FEType &type)

Member Data Documentation

const std::set<boundary_id_type>& libMesh::BoundaryProjectSolution::b
private

Definition at line 524 of file system_projection.C.

UniquePtr<FunctionBase<Number> > libMesh::BoundaryProjectSolution::f
private

Definition at line 527 of file system_projection.C.

UniquePtr<FunctionBase<Gradient> > libMesh::BoundaryProjectSolution::g
private

Definition at line 528 of file system_projection.C.

NumericVector<Number>& libMesh::BoundaryProjectSolution::new_vector
private

Definition at line 530 of file system_projection.C.

const Parameters& libMesh::BoundaryProjectSolution::parameters
private

Definition at line 529 of file system_projection.C.

const System& libMesh::BoundaryProjectSolution::system
private

Definition at line 526 of file system_projection.C.

const std::vector<unsigned int>& libMesh::BoundaryProjectSolution::variables
private

Definition at line 525 of file system_projection.C.


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