libMesh
Public Member Functions | Private Attributes | List of all members
libMesh::BoundaryProjectSolution Class Reference

This class implements projecting an arbitrary boundary function to the current mesh. More...

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 
)

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)

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 2126 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(), dim, 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().

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