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 exectued 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 2130 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().

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