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
 
std::unique_ptr< FunctionBase< Number > > f
 
std::unique_ptr< 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 147 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 159 of file system_projection.C.

References libMesh::NumericVector< T >::init().

165  :
166  b(b_in),
167  variables(variables_in),
168  system(system_in),
169  f(f_in ? f_in->clone() : std::unique_ptr<FunctionBase<Number>>()),
170  g(g_in ? g_in->clone() : std::unique_ptr<FunctionBase<Gradient>>()),
171  parameters(parameters_in),
172  new_vector(new_v_in)
173  {
174  libmesh_assert(f.get());
175  f->init();
176  if (g.get())
177  g->init();
178  }
std::unique_ptr< FunctionBase< Number > > f
const std::vector< unsigned int > & variables
NumericVector< Number > & new_vector
virtual std::unique_ptr< FunctionBase< Output > > clone() const =0
const std::set< boundary_id_type > & b
std::unique_ptr< FunctionBase< Gradient > > g
libMesh::BoundaryProjectSolution::BoundaryProjectSolution ( const BoundaryProjectSolution in)
inline

Definition at line 180 of file system_projection.C.

References libMesh::NumericVector< T >::init().

180  :
181  b(in.b),
182  variables(in.variables),
183  system(in.system),
184  f(in.f.get() ? in.f->clone() : std::unique_ptr<FunctionBase<Number>>()),
185  g(in.g.get() ? in.g->clone() : std::unique_ptr<FunctionBase<Gradient>>()),
186  parameters(in.parameters),
187  new_vector(in.new_vector)
188  {
189  libmesh_assert(f.get());
190  f->init();
191  if (g.get())
192  g->init();
193  }
std::unique_ptr< FunctionBase< Number > > f
const std::vector< unsigned int > & variables
NumericVector< Number > & new_vector
const std::set< boundary_id_type > & b
std::unique_ptr< FunctionBase< Gradient > > g

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 1178 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::HERMITE, libMesh::Elem::is_edge_on_side(), libMesh::Elem::is_node_on_side(), libMesh::Elem::is_vertex(), libmesh_nullptr, 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::TOLERANCE, libMesh::Variable::type(), libMesh::Elem::type(), libMesh::DofMap::variable(), libMesh::DenseMatrix< T >::zero(), and libMesh::DenseVector< T >::zero().

1179 {
1180  // We need data to project
1181  libmesh_assert(f.get());
1182 
1190  // The dimensionality of the current mesh
1191  const unsigned int dim = system.get_mesh().mesh_dimension();
1192 
1193  // The DofMap for this system
1194  const DofMap & dof_map = system.get_dof_map();
1195 
1196  // Boundary info for the current mesh
1197  const BoundaryInfo & boundary_info =
1199 
1200  // The element matrix and RHS for projections.
1201  // Note that Ke is always real-valued, whereas
1202  // Fe may be complex valued if complex number
1203  // support is enabled
1204  DenseMatrix<Real> Ke;
1205  DenseVector<Number> Fe;
1206  // The new element coefficients
1207  DenseVector<Number> Ue;
1208 
1209 
1210  // Loop over all the variables we've been requested to project
1211  for (std::size_t v=0; v!=variables.size(); v++)
1212  {
1213  const unsigned int var = variables[v];
1214 
1215  const Variable & variable = dof_map.variable(var);
1216 
1217  const FEType & fe_type = variable.type();
1218 
1219  if (fe_type.family == SCALAR)
1220  continue;
1221 
1222  const unsigned int var_component =
1224 
1225  // Get FE objects of the appropriate type
1226  std::unique_ptr<FEBase> fe (FEBase::build(dim, fe_type));
1227 
1228  // Prepare variables for projection
1229  std::unique_ptr<QBase> qedgerule (fe_type.default_quadrature_rule(1));
1230  std::unique_ptr<QBase> qsiderule (fe_type.default_quadrature_rule(dim-1));
1231 
1232  // The values of the shape functions at the quadrature
1233  // points
1234  const std::vector<std::vector<Real>> & phi = fe->get_phi();
1235 
1236  // The gradients of the shape functions at the quadrature
1237  // points on the child element.
1238  const std::vector<std::vector<RealGradient>> * dphi = libmesh_nullptr;
1239 
1240  const FEContinuity cont = fe->get_continuity();
1241 
1242  if (cont == C_ONE)
1243  {
1244  // We'll need gradient data for a C1 projection
1245  libmesh_assert(g.get());
1246 
1247  const std::vector<std::vector<RealGradient>> &
1248  ref_dphi = fe->get_dphi();
1249  dphi = &ref_dphi;
1250  }
1251 
1252  // The Jacobian * quadrature weight at the quadrature points
1253  const std::vector<Real> & JxW =
1254  fe->get_JxW();
1255 
1256  // The XYZ locations of the quadrature points
1257  const std::vector<Point> & xyz_values =
1258  fe->get_xyz();
1259 
1260  // The global DOF indices
1261  std::vector<dof_id_type> dof_indices;
1262  // Side/edge DOF indices
1263  std::vector<unsigned int> side_dofs;
1264 
1265  // Container to catch IDs passed back from BoundaryInfo.
1266  std::vector<boundary_id_type> bc_ids;
1267 
1268  // Iterate over all the elements in the range
1269  for (ConstElemRange::const_iterator elem_it=range.begin(); elem_it != range.end(); ++elem_it)
1270  {
1271  const Elem * elem = *elem_it;
1272 
1273  // Per-subdomain variables don't need to be projected on
1274  // elements where they're not active
1275  if (!variable.active_on_subdomain(elem->subdomain_id()))
1276  continue;
1277 
1278  const unsigned short n_nodes = elem->n_nodes();
1279  const unsigned short n_edges = elem->n_edges();
1280  const unsigned short n_sides = elem->n_sides();
1281 
1282  // Find out which nodes, edges and sides are on a requested
1283  // boundary:
1284  std::vector<bool> is_boundary_node(n_nodes, false),
1285  is_boundary_edge(n_edges, false),
1286  is_boundary_side(n_sides, false);
1287  for (unsigned char s=0; s != n_sides; ++s)
1288  {
1289  // First see if this side has been requested
1290  boundary_info.boundary_ids (elem, s, bc_ids);
1291  bool do_this_side = false;
1292  for (std::vector<boundary_id_type>::iterator i=bc_ids.begin();
1293  i!=bc_ids.end(); ++i)
1294  if (b.count(*i))
1295  {
1296  do_this_side = true;
1297  break;
1298  }
1299  if (!do_this_side)
1300  continue;
1301 
1302  is_boundary_side[s] = true;
1303 
1304  // Then see what nodes and what edges are on it
1305  for (unsigned int n=0; n != n_nodes; ++n)
1306  if (elem->is_node_on_side(n,s))
1307  is_boundary_node[n] = true;
1308  for (unsigned int e=0; e != n_edges; ++e)
1309  if (elem->is_edge_on_side(e,s))
1310  is_boundary_edge[e] = true;
1311  }
1312 
1313  // Update the DOF indices for this element based on
1314  // the current mesh
1315  dof_map.dof_indices (elem, dof_indices, var);
1316 
1317  // The number of DOFs on the element
1318  const unsigned int n_dofs =
1319  cast_int<unsigned int>(dof_indices.size());
1320 
1321  // Fixed vs. free DoFs on edge/face projections
1322  std::vector<char> dof_is_fixed(n_dofs, false); // bools
1323  std::vector<int> free_dof(n_dofs, 0);
1324 
1325  // The element type
1326  const ElemType elem_type = elem->type();
1327 
1328  // Zero the interpolated values
1329  Ue.resize (n_dofs); Ue.zero();
1330 
1331  // In general, we need a series of
1332  // projections to ensure a unique and continuous
1333  // solution. We start by interpolating boundary nodes, then
1334  // hold those fixed and project boundary edges, then hold
1335  // those fixed and project boundary faces,
1336 
1337  // Interpolate node values first
1338  unsigned int current_dof = 0;
1339  for (unsigned short n = 0; n != n_nodes; ++n)
1340  {
1341  // FIXME: this should go through the DofMap,
1342  // not duplicate dof_indices code badly!
1343  const unsigned int nc =
1344  FEInterface::n_dofs_at_node (dim, fe_type, elem_type,
1345  n);
1346  if (!elem->is_vertex(n) || !is_boundary_node[n])
1347  {
1348  current_dof += nc;
1349  continue;
1350  }
1351  if (cont == DISCONTINUOUS)
1352  {
1353  libmesh_assert_equal_to (nc, 0);
1354  }
1355  // Assume that C_ZERO elements have a single nodal
1356  // value shape function
1357  else if (cont == C_ZERO)
1358  {
1359  libmesh_assert_equal_to (nc, 1);
1360  Ue(current_dof) = f->component(var_component,
1361  elem->point(n),
1362  system.time);
1363  dof_is_fixed[current_dof] = true;
1364  current_dof++;
1365  }
1366  // The hermite element vertex shape functions are weird
1367  else if (fe_type.family == HERMITE)
1368  {
1369  Ue(current_dof) = f->component(var_component,
1370  elem->point(n),
1371  system.time);
1372  dof_is_fixed[current_dof] = true;
1373  current_dof++;
1374  Gradient grad = g->component(var_component,
1375  elem->point(n),
1376  system.time);
1377  // x derivative
1378  Ue(current_dof) = grad(0);
1379  dof_is_fixed[current_dof] = true;
1380  current_dof++;
1381  if (dim > 1)
1382  {
1383  // We'll finite difference mixed derivatives
1384  Point nxminus = elem->point(n),
1385  nxplus = elem->point(n);
1386  nxminus(0) -= TOLERANCE;
1387  nxplus(0) += TOLERANCE;
1388  Gradient gxminus = g->component(var_component,
1389  nxminus,
1390  system.time);
1391  Gradient gxplus = g->component(var_component,
1392  nxplus,
1393  system.time);
1394  // y derivative
1395  Ue(current_dof) = grad(1);
1396  dof_is_fixed[current_dof] = true;
1397  current_dof++;
1398  // xy derivative
1399  Ue(current_dof) = (gxplus(1) - gxminus(1))
1400  / 2. / TOLERANCE;
1401  dof_is_fixed[current_dof] = true;
1402  current_dof++;
1403 
1404  if (dim > 2)
1405  {
1406  // z derivative
1407  Ue(current_dof) = grad(2);
1408  dof_is_fixed[current_dof] = true;
1409  current_dof++;
1410  // xz derivative
1411  Ue(current_dof) = (gxplus(2) - gxminus(2))
1412  / 2. / TOLERANCE;
1413  dof_is_fixed[current_dof] = true;
1414  current_dof++;
1415  // We need new points for yz
1416  Point nyminus = elem->point(n),
1417  nyplus = elem->point(n);
1418  nyminus(1) -= TOLERANCE;
1419  nyplus(1) += TOLERANCE;
1420  Gradient gyminus = g->component(var_component,
1421  nyminus,
1422  system.time);
1423  Gradient gyplus = g->component(var_component,
1424  nyplus,
1425  system.time);
1426  // xz derivative
1427  Ue(current_dof) = (gyplus(2) - gyminus(2))
1428  / 2. / TOLERANCE;
1429  dof_is_fixed[current_dof] = true;
1430  current_dof++;
1431  // Getting a 2nd order xyz is more tedious
1432  Point nxmym = elem->point(n),
1433  nxmyp = elem->point(n),
1434  nxpym = elem->point(n),
1435  nxpyp = elem->point(n);
1436  nxmym(0) -= TOLERANCE;
1437  nxmym(1) -= TOLERANCE;
1438  nxmyp(0) -= TOLERANCE;
1439  nxmyp(1) += TOLERANCE;
1440  nxpym(0) += TOLERANCE;
1441  nxpym(1) -= TOLERANCE;
1442  nxpyp(0) += TOLERANCE;
1443  nxpyp(1) += TOLERANCE;
1444  Gradient gxmym = g->component(var_component,
1445  nxmym,
1446  system.time);
1447  Gradient gxmyp = g->component(var_component,
1448  nxmyp,
1449  system.time);
1450  Gradient gxpym = g->component(var_component,
1451  nxpym,
1452  system.time);
1453  Gradient gxpyp = g->component(var_component,
1454  nxpyp,
1455  system.time);
1456  Number gxzplus = (gxpyp(2) - gxmyp(2))
1457  / 2. / TOLERANCE;
1458  Number gxzminus = (gxpym(2) - gxmym(2))
1459  / 2. / TOLERANCE;
1460  // xyz derivative
1461  Ue(current_dof) = (gxzplus - gxzminus)
1462  / 2. / TOLERANCE;
1463  dof_is_fixed[current_dof] = true;
1464  current_dof++;
1465  }
1466  }
1467  }
1468  // Assume that other C_ONE elements have a single nodal
1469  // value shape function and nodal gradient component
1470  // shape functions
1471  else if (cont == C_ONE)
1472  {
1473  libmesh_assert_equal_to (nc, 1 + dim);
1474  Ue(current_dof) = f->component(var_component,
1475  elem->point(n),
1476  system.time);
1477  dof_is_fixed[current_dof] = true;
1478  current_dof++;
1479  Gradient grad = g->component(var_component,
1480  elem->point(n),
1481  system.time);
1482  for (unsigned int i=0; i!= dim; ++i)
1483  {
1484  Ue(current_dof) = grad(i);
1485  dof_is_fixed[current_dof] = true;
1486  current_dof++;
1487  }
1488  }
1489  else
1490  libmesh_error_msg("Unknown continuity " << cont);
1491  }
1492 
1493  // In 3D, project any edge values next
1494  if (dim > 2 && cont != DISCONTINUOUS)
1495  for (unsigned short e = 0; e != n_edges; ++e)
1496  {
1497  if (!is_boundary_edge[e])
1498  continue;
1499 
1500  FEInterface::dofs_on_edge(elem, dim, fe_type, e,
1501  side_dofs);
1502 
1503  // Some edge dofs are on nodes and already
1504  // fixed, others are free to calculate
1505  unsigned int free_dofs = 0;
1506  for (std::size_t i=0; i != side_dofs.size(); ++i)
1507  if (!dof_is_fixed[side_dofs[i]])
1508  free_dof[free_dofs++] = i;
1509 
1510  // There may be nothing to project
1511  if (!free_dofs)
1512  continue;
1513 
1514  Ke.resize (free_dofs, free_dofs); Ke.zero();
1515  Fe.resize (free_dofs); Fe.zero();
1516  // The new edge coefficients
1517  DenseVector<Number> Uedge(free_dofs);
1518 
1519  // Initialize FE data on the edge
1520  fe->attach_quadrature_rule (qedgerule.get());
1521  fe->edge_reinit (elem, e);
1522  const unsigned int n_qp = qedgerule->n_points();
1523 
1524  // Loop over the quadrature points
1525  for (unsigned int qp=0; qp<n_qp; qp++)
1526  {
1527  // solution at the quadrature point
1528  Number fineval = f->component(var_component,
1529  xyz_values[qp],
1530  system.time);
1531  // solution grad at the quadrature point
1532  Gradient finegrad;
1533  if (cont == C_ONE)
1534  finegrad = g->component(var_component,
1535  xyz_values[qp],
1536  system.time);
1537 
1538  // Form edge projection matrix
1539  for (std::size_t sidei=0, freei=0;
1540  sidei != side_dofs.size(); ++sidei)
1541  {
1542  unsigned int i = side_dofs[sidei];
1543  // fixed DoFs aren't test functions
1544  if (dof_is_fixed[i])
1545  continue;
1546  for (std::size_t sidej=0, freej=0;
1547  sidej != side_dofs.size(); ++sidej)
1548  {
1549  unsigned int j = side_dofs[sidej];
1550  if (dof_is_fixed[j])
1551  Fe(freei) -= phi[i][qp] * phi[j][qp] *
1552  JxW[qp] * Ue(j);
1553  else
1554  Ke(freei,freej) += phi[i][qp] *
1555  phi[j][qp] * JxW[qp];
1556  if (cont == C_ONE)
1557  {
1558  if (dof_is_fixed[j])
1559  Fe(freei) -= ((*dphi)[i][qp] *
1560  (*dphi)[j][qp]) *
1561  JxW[qp] * Ue(j);
1562  else
1563  Ke(freei,freej) += ((*dphi)[i][qp] *
1564  (*dphi)[j][qp])
1565  * JxW[qp];
1566  }
1567  if (!dof_is_fixed[j])
1568  freej++;
1569  }
1570  Fe(freei) += phi[i][qp] * fineval * JxW[qp];
1571  if (cont == C_ONE)
1572  Fe(freei) += (finegrad * (*dphi)[i][qp]) *
1573  JxW[qp];
1574  freei++;
1575  }
1576  }
1577 
1578  Ke.cholesky_solve(Fe, Uedge);
1579 
1580  // Transfer new edge solutions to element
1581  for (unsigned int i=0; i != free_dofs; ++i)
1582  {
1583  Number & ui = Ue(side_dofs[free_dof[i]]);
1584  libmesh_assert(std::abs(ui) < TOLERANCE ||
1585  std::abs(ui - Uedge(i)) < TOLERANCE);
1586  ui = Uedge(i);
1587  dof_is_fixed[side_dofs[free_dof[i]]] = true;
1588  }
1589  }
1590 
1591  // Project any side values (edges in 2D, faces in 3D)
1592  if (dim > 1 && cont != DISCONTINUOUS)
1593  for (unsigned short s = 0; s != n_sides; ++s)
1594  {
1595  if (!is_boundary_side[s])
1596  continue;
1597 
1598  FEInterface::dofs_on_side(elem, dim, fe_type, s,
1599  side_dofs);
1600 
1601  // Some side dofs are on nodes/edges and already
1602  // fixed, others are free to calculate
1603  unsigned int free_dofs = 0;
1604  for (std::size_t i=0; i != side_dofs.size(); ++i)
1605  if (!dof_is_fixed[side_dofs[i]])
1606  free_dof[free_dofs++] = i;
1607 
1608  // There may be nothing to project
1609  if (!free_dofs)
1610  continue;
1611 
1612  Ke.resize (free_dofs, free_dofs); Ke.zero();
1613  Fe.resize (free_dofs); Fe.zero();
1614  // The new side coefficients
1615  DenseVector<Number> Uside(free_dofs);
1616 
1617  // Initialize FE data on the side
1618  fe->attach_quadrature_rule (qsiderule.get());
1619  fe->reinit (elem, s);
1620  const unsigned int n_qp = qsiderule->n_points();
1621 
1622  // Loop over the quadrature points
1623  for (unsigned int qp=0; qp<n_qp; qp++)
1624  {
1625  // solution at the quadrature point
1626  Number fineval = f->component(var_component,
1627  xyz_values[qp],
1628  system.time);
1629  // solution grad at the quadrature point
1630  Gradient finegrad;
1631  if (cont == C_ONE)
1632  finegrad = g->component(var_component,
1633  xyz_values[qp],
1634  system.time);
1635 
1636  // Form side projection matrix
1637  for (std::size_t sidei=0, freei=0;
1638  sidei != side_dofs.size(); ++sidei)
1639  {
1640  unsigned int i = side_dofs[sidei];
1641  // fixed DoFs aren't test functions
1642  if (dof_is_fixed[i])
1643  continue;
1644  for (std::size_t sidej=0, freej=0;
1645  sidej != side_dofs.size(); ++sidej)
1646  {
1647  unsigned int j = side_dofs[sidej];
1648  if (dof_is_fixed[j])
1649  Fe(freei) -= phi[i][qp] * phi[j][qp] *
1650  JxW[qp] * Ue(j);
1651  else
1652  Ke(freei,freej) += phi[i][qp] *
1653  phi[j][qp] * JxW[qp];
1654  if (cont == C_ONE)
1655  {
1656  if (dof_is_fixed[j])
1657  Fe(freei) -= ((*dphi)[i][qp] *
1658  (*dphi)[j][qp]) *
1659  JxW[qp] * Ue(j);
1660  else
1661  Ke(freei,freej) += ((*dphi)[i][qp] *
1662  (*dphi)[j][qp])
1663  * JxW[qp];
1664  }
1665  if (!dof_is_fixed[j])
1666  freej++;
1667  }
1668  Fe(freei) += (fineval * phi[i][qp]) * JxW[qp];
1669  if (cont == C_ONE)
1670  Fe(freei) += (finegrad * (*dphi)[i][qp]) *
1671  JxW[qp];
1672  freei++;
1673  }
1674  }
1675 
1676  Ke.cholesky_solve(Fe, Uside);
1677 
1678  // Transfer new side solutions to element
1679  for (unsigned int i=0; i != free_dofs; ++i)
1680  {
1681  Number & ui = Ue(side_dofs[free_dof[i]]);
1682  libmesh_assert(std::abs(ui) < TOLERANCE ||
1683  std::abs(ui - Uside(i)) < TOLERANCE);
1684  ui = Uside(i);
1685  dof_is_fixed[side_dofs[free_dof[i]]] = true;
1686  }
1687  }
1688 
1689  const dof_id_type
1690  first = new_vector.first_local_index(),
1691  last = new_vector.last_local_index();
1692 
1693  // Lock the new_vector since it is shared among threads.
1694  {
1695  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1696 
1697  for (unsigned int i = 0; i < n_dofs; i++)
1698  if (dof_is_fixed[i] &&
1699  (dof_indices[i] >= first) &&
1700  (dof_indices[i] < last))
1701  {
1702  new_vector.set(dof_indices[i], Ue(i));
1703  }
1704  }
1705  } // end elem loop
1706  } // end variables loop
1707 }
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:493
const BoundaryInfo & get_boundary_info() const
Definition: mesh_base.h:118
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:480
std::unique_ptr< FunctionBase< Number > > f
const class libmesh_nullptr_t libmesh_nullptr
static const Real TOLERANCE
const std::vector< unsigned int > & variables
const dof_id_type n_nodes
Definition: tecplot_io.C:67
const MeshBase & get_mesh() const
Definition: system.h:2014
spin_mutex spin_mtx
Definition: threads.C:29
const DofMap & get_dof_map() const
Definition: system.h:2030
NumericVector< Number > & new_vector
NumberVectorValue Gradient
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
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:430
unsigned int mesh_dimension() const
Definition: mesh_base.C:150
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
const std::set< boundary_id_type > & b
std::unique_ptr< FunctionBase< Gradient > > g
uint8_t dof_id_type
Definition: id_types.h:64

Member Data Documentation

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

Definition at line 150 of file system_projection.C.

std::unique_ptr<FunctionBase<Number> > libMesh::BoundaryProjectSolution::f
private

Definition at line 153 of file system_projection.C.

std::unique_ptr<FunctionBase<Gradient> > libMesh::BoundaryProjectSolution::g
private

Definition at line 154 of file system_projection.C.

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

Definition at line 156 of file system_projection.C.

const Parameters& libMesh::BoundaryProjectSolution::parameters
private

Definition at line 155 of file system_projection.C.

const System& libMesh::BoundaryProjectSolution::system
private

Definition at line 152 of file system_projection.C.

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

Definition at line 151 of file system_projection.C.


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