libMesh::ProjectSolution Class Reference

Public Member Functions

 ProjectSolution (const System &system_in, FunctionBase< Number > *f_in, FunctionBase< Gradient > *g_in, const Parameters &parameters_in, NumericVector< Number > &new_v_in)
 
 ProjectSolution (const ProjectSolution &in)
 
void operator() (const ConstElemRange &range) const
 

Private Attributes

const Systemsystem
 
AutoPtr< FunctionBase< Number > > f
 
AutoPtr< FunctionBase< Gradient > > g
 
const Parametersparameters
 
NumericVector< Number > & new_vector
 

Detailed Description

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

Definition at line 110 of file system_projection.C.

Constructor & Destructor Documentation

libMesh::ProjectSolution::ProjectSolution ( 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 121 of file system_projection.C.

References f, g, and libMesh::libmesh_assert().

125  :
126  system(system_in),
127  f(f_in ? f_in->clone() : AutoPtr<FunctionBase<Number> >(NULL)),
128  g(g_in ? g_in->clone() : AutoPtr<FunctionBase<Gradient> >(NULL)),
129  parameters(parameters_in),
130  new_vector(new_v_in)
131  {
132  libmesh_assert(f.get());
133  f->init();
134  if (g.get())
135  g->init();
136  }
libMesh::ProjectSolution::ProjectSolution ( const ProjectSolution in)
inline

Definition at line 138 of file system_projection.C.

References f, g, and libMesh::libmesh_assert().

138  :
139  system(in.system),
140  f(in.f.get() ? in.f->clone() : AutoPtr<FunctionBase<Number> >(NULL)),
141  g(in.g.get() ? in.g->clone() : AutoPtr<FunctionBase<Gradient> >(NULL)),
142  parameters(in.parameters),
143  new_vector(in.new_vector)
144  {
145  libmesh_assert(f.get());
146  f->init();
147  if (g.get())
148  g->init();
149  }

Member Function Documentation

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

This method projects an arbitrary 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 1252 of file system_projection.C.

References std::abs(), libMesh::Variable::active_on_subdomain(), libMesh::StoredRange< iterator_type, object_type >::begin(), libMesh::FEGenericBase< T >::build(), libMeshEnums::C_ONE, libMeshEnums::C_ZERO, libMesh::DenseMatrix< T >::cholesky_solve(), libMesh::FunctionBase< Output >::component(), libMesh::FEType::default_quadrature_rule(), libMesh::dim, libMeshEnums::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::NumericVector< T >::first_local_index(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMeshEnums::HERMITE, libMesh::Elem::is_vertex(), libMesh::NumericVector< T >::last_local_index(), libMesh::libmesh_assert(), 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::System::n_vars(), libMesh::Elem::point(), libMesh::DenseVector< T >::resize(), libMesh::DenseMatrix< T >::resize(), libMeshEnums::SCALAR, libMesh::NumericVector< T >::set(), libMesh::Threads::spin_mtx, libMesh::Elem::subdomain_id(), system, libMesh::System::time, libMesh::TOLERANCE, libMesh::Variable::type(), libMesh::Elem::type(), libMesh::DofMap::variable(), libMesh::System::variable_scalar_number(), libMesh::DenseMatrix< T >::zero(), and libMesh::DenseVector< T >::zero().

1253 {
1254  // We need data to project
1255  libmesh_assert(f.get());
1256 
1264  // The number of variables in this system
1265  const unsigned int n_variables = system.n_vars();
1266 
1267  // The dimensionality of the current mesh
1268  const unsigned int dim = system.get_mesh().mesh_dimension();
1269 
1270  // The DofMap for this system
1271  const DofMap& dof_map = system.get_dof_map();
1272 
1273  // The element matrix and RHS for projections.
1274  // Note that Ke is always real-valued, whereas
1275  // Fe may be complex valued if complex number
1276  // support is enabled
1277  DenseMatrix<Real> Ke;
1278  DenseVector<Number> Fe;
1279  // The new element coefficients
1280  DenseVector<Number> Ue;
1281 
1282 
1283  // Loop over all the variables in the system
1284  for (unsigned int var=0; var<n_variables; var++)
1285  {
1286  const Variable& variable = dof_map.variable(var);
1287 
1288  const FEType& fe_type = variable.type();
1289 
1290  if (fe_type.family == SCALAR)
1291  continue;
1292 
1293  const unsigned int var_component =
1295 
1296  // Get FE objects of the appropriate type
1297  AutoPtr<FEBase> fe (FEBase::build(dim, fe_type));
1298 
1299  // Prepare variables for projection
1300  AutoPtr<QBase> qrule (fe_type.default_quadrature_rule(dim));
1301  AutoPtr<QBase> qedgerule (fe_type.default_quadrature_rule(1));
1302  AutoPtr<QBase> qsiderule (fe_type.default_quadrature_rule(dim-1));
1303 
1304  // The values of the shape functions at the quadrature
1305  // points
1306  const std::vector<std::vector<Real> >& phi = fe->get_phi();
1307 
1308  // The gradients of the shape functions at the quadrature
1309  // points on the child element.
1310  const std::vector<std::vector<RealGradient> > *dphi = NULL;
1311 
1312  const FEContinuity cont = fe->get_continuity();
1313 
1314  if (cont == C_ONE)
1315  {
1316  // We'll need gradient data for a C1 projection
1317  libmesh_assert(g.get());
1318 
1319  const std::vector<std::vector<RealGradient> >&
1320  ref_dphi = fe->get_dphi();
1321  dphi = &ref_dphi;
1322  }
1323 
1324  // The Jacobian * quadrature weight at the quadrature points
1325  const std::vector<Real>& JxW =
1326  fe->get_JxW();
1327 
1328  // The XYZ locations of the quadrature points
1329  const std::vector<Point>& xyz_values =
1330  fe->get_xyz();
1331 
1332  // The global DOF indices
1333  std::vector<dof_id_type> dof_indices;
1334  // Side/edge DOF indices
1335  std::vector<unsigned int> side_dofs;
1336 
1337  // Iterate over all the elements in the range
1338  for (ConstElemRange::const_iterator elem_it=range.begin(); elem_it != range.end(); ++elem_it)
1339  {
1340  const Elem* elem = *elem_it;
1341 
1342  // Per-subdomain variables don't need to be projected on
1343  // elements where they're not active
1344  if (!variable.active_on_subdomain(elem->subdomain_id()))
1345  continue;
1346 
1347  // Update the DOF indices for this element based on
1348  // the current mesh
1349  dof_map.dof_indices (elem, dof_indices, var);
1350 
1351  // The number of DOFs on the element
1352  const unsigned int n_dofs =
1353  libmesh_cast_int<unsigned int>(dof_indices.size());
1354 
1355  // Fixed vs. free DoFs on edge/face projections
1356  std::vector<char> dof_is_fixed(n_dofs, false); // bools
1357  std::vector<int> free_dof(n_dofs, 0);
1358 
1359  // The element type
1360  const ElemType elem_type = elem->type();
1361 
1362  // The number of nodes on the new element
1363  const unsigned int n_nodes = elem->n_nodes();
1364 
1365  // Zero the interpolated values
1366  Ue.resize (n_dofs); Ue.zero();
1367 
1368  // In general, we need a series of
1369  // projections to ensure a unique and continuous
1370  // solution. We start by interpolating nodes, then
1371  // hold those fixed and project edges, then
1372  // hold those fixed and project faces, then
1373  // hold those fixed and project interiors
1374 
1375  // Interpolate node values first
1376  unsigned int current_dof = 0;
1377  for (unsigned int n=0; n!= n_nodes; ++n)
1378  {
1379  // FIXME: this should go through the DofMap,
1380  // not duplicate dof_indices code badly!
1381  const unsigned int nc =
1382  FEInterface::n_dofs_at_node (dim, fe_type, elem_type,
1383  n);
1384  if (!elem->is_vertex(n))
1385  {
1386  current_dof += nc;
1387  continue;
1388  }
1389  if (cont == DISCONTINUOUS)
1390  {
1391  libmesh_assert_equal_to (nc, 0);
1392  }
1393  // Assume that C_ZERO elements have a single nodal
1394  // value shape function
1395  else if (cont == C_ZERO)
1396  {
1397  libmesh_assert_equal_to (nc, 1);
1398  Ue(current_dof) = f->component(var_component,
1399  elem->point(n),
1400  system.time);
1401  dof_is_fixed[current_dof] = true;
1402  current_dof++;
1403  }
1404  // The hermite element vertex shape functions are weird
1405  else if (fe_type.family == HERMITE)
1406  {
1407  Ue(current_dof) = f->component(var_component,
1408  elem->point(n),
1409  system.time);
1410  dof_is_fixed[current_dof] = true;
1411  current_dof++;
1412  Gradient grad = g->component(var_component,
1413  elem->point(n),
1414  system.time);
1415  // x derivative
1416  Ue(current_dof) = grad(0);
1417  dof_is_fixed[current_dof] = true;
1418  current_dof++;
1419  if (dim > 1)
1420  {
1421  // We'll finite difference mixed derivatives
1422  Point nxminus = elem->point(n),
1423  nxplus = elem->point(n);
1424  nxminus(0) -= TOLERANCE;
1425  nxplus(0) += TOLERANCE;
1426  Gradient gxminus = g->component(var_component,
1427  nxminus,
1428  system.time);
1429  Gradient gxplus = g->component(var_component,
1430  nxplus,
1431  system.time);
1432  // y derivative
1433  Ue(current_dof) = grad(1);
1434  dof_is_fixed[current_dof] = true;
1435  current_dof++;
1436  // xy derivative
1437  Ue(current_dof) = (gxplus(1) - gxminus(1))
1438  / 2. / TOLERANCE;
1439  dof_is_fixed[current_dof] = true;
1440  current_dof++;
1441 
1442  if (dim > 2)
1443  {
1444  // z derivative
1445  Ue(current_dof) = grad(2);
1446  dof_is_fixed[current_dof] = true;
1447  current_dof++;
1448  // xz derivative
1449  Ue(current_dof) = (gxplus(2) - gxminus(2))
1450  / 2. / TOLERANCE;
1451  dof_is_fixed[current_dof] = true;
1452  current_dof++;
1453  // We need new points for yz
1454  Point nyminus = elem->point(n),
1455  nyplus = elem->point(n);
1456  nyminus(1) -= TOLERANCE;
1457  nyplus(1) += TOLERANCE;
1458  Gradient gyminus = g->component(var_component,
1459  nyminus,
1460  system.time);
1461  Gradient gyplus = g->component(var_component,
1462  nyplus,
1463  system.time);
1464  // xz derivative
1465  Ue(current_dof) = (gyplus(2) - gyminus(2))
1466  / 2. / TOLERANCE;
1467  dof_is_fixed[current_dof] = true;
1468  current_dof++;
1469  // Getting a 2nd order xyz is more tedious
1470  Point nxmym = elem->point(n),
1471  nxmyp = elem->point(n),
1472  nxpym = elem->point(n),
1473  nxpyp = elem->point(n);
1474  nxmym(0) -= TOLERANCE;
1475  nxmym(1) -= TOLERANCE;
1476  nxmyp(0) -= TOLERANCE;
1477  nxmyp(1) += TOLERANCE;
1478  nxpym(0) += TOLERANCE;
1479  nxpym(1) -= TOLERANCE;
1480  nxpyp(0) += TOLERANCE;
1481  nxpyp(1) += TOLERANCE;
1482  Gradient gxmym = g->component(var_component,
1483  nxmym,
1484  system.time);
1485  Gradient gxmyp = g->component(var_component,
1486  nxmyp,
1487  system.time);
1488  Gradient gxpym = g->component(var_component,
1489  nxpym,
1490  system.time);
1491  Gradient gxpyp = g->component(var_component,
1492  nxpyp,
1493  system.time);
1494  Number gxzplus = (gxpyp(2) - gxmyp(2))
1495  / 2. / TOLERANCE;
1496  Number gxzminus = (gxpym(2) - gxmym(2))
1497  / 2. / TOLERANCE;
1498  // xyz derivative
1499  Ue(current_dof) = (gxzplus - gxzminus)
1500  / 2. / TOLERANCE;
1501  dof_is_fixed[current_dof] = true;
1502  current_dof++;
1503  }
1504  }
1505  }
1506  // Assume that other C_ONE elements have a single nodal
1507  // value shape function and nodal gradient component
1508  // shape functions
1509  else if (cont == C_ONE)
1510  {
1511  libmesh_assert_equal_to (nc, 1 + dim);
1512  Ue(current_dof) = f->component(var_component,
1513  elem->point(n),
1514  system.time);
1515  dof_is_fixed[current_dof] = true;
1516  current_dof++;
1517  Gradient grad = g->component(var_component,
1518  elem->point(n),
1519  system.time);
1520  for (unsigned int i=0; i!= dim; ++i)
1521  {
1522  Ue(current_dof) = grad(i);
1523  dof_is_fixed[current_dof] = true;
1524  current_dof++;
1525  }
1526  }
1527  else
1528  libmesh_error();
1529  }
1530 
1531  // In 3D, project any edge values next
1532  if (dim > 2 && cont != DISCONTINUOUS)
1533  for (unsigned int e=0; e != elem->n_edges(); ++e)
1534  {
1535  FEInterface::dofs_on_edge(elem, dim, fe_type, e,
1536  side_dofs);
1537 
1538  // Some edge dofs are on nodes and already
1539  // fixed, others are free to calculate
1540  unsigned int free_dofs = 0;
1541  for (unsigned int i=0; i != side_dofs.size(); ++i)
1542  if (!dof_is_fixed[side_dofs[i]])
1543  free_dof[free_dofs++] = i;
1544 
1545  // There may be nothing to project
1546  if (!free_dofs)
1547  continue;
1548 
1549  Ke.resize (free_dofs, free_dofs); Ke.zero();
1550  Fe.resize (free_dofs); Fe.zero();
1551  // The new edge coefficients
1552  DenseVector<Number> Uedge(free_dofs);
1553 
1554  // Initialize FE data on the edge
1555  fe->attach_quadrature_rule (qedgerule.get());
1556  fe->edge_reinit (elem, e);
1557  const unsigned int n_qp = qedgerule->n_points();
1558 
1559  // Loop over the quadrature points
1560  for (unsigned int qp=0; qp<n_qp; qp++)
1561  {
1562  // solution at the quadrature point
1563  Number fineval = f->component(var_component,
1564  xyz_values[qp],
1565  system.time);
1566  // solution grad at the quadrature point
1567  Gradient finegrad;
1568  if (cont == C_ONE)
1569  finegrad = g->component(var_component,
1570  xyz_values[qp],
1571  system.time);
1572 
1573  // Form edge projection matrix
1574  for (unsigned int sidei=0, freei=0;
1575  sidei != side_dofs.size(); ++sidei)
1576  {
1577  unsigned int i = side_dofs[sidei];
1578  // fixed DoFs aren't test functions
1579  if (dof_is_fixed[i])
1580  continue;
1581  for (unsigned int sidej=0, freej=0;
1582  sidej != side_dofs.size(); ++sidej)
1583  {
1584  unsigned int j = side_dofs[sidej];
1585  if (dof_is_fixed[j])
1586  Fe(freei) -= phi[i][qp] * phi[j][qp] *
1587  JxW[qp] * Ue(j);
1588  else
1589  Ke(freei,freej) += phi[i][qp] *
1590  phi[j][qp] * JxW[qp];
1591  if (cont == C_ONE)
1592  {
1593  if (dof_is_fixed[j])
1594  Fe(freei) -= ((*dphi)[i][qp] *
1595  (*dphi)[j][qp]) *
1596  JxW[qp] * Ue(j);
1597  else
1598  Ke(freei,freej) += ((*dphi)[i][qp] *
1599  (*dphi)[j][qp])
1600  * JxW[qp];
1601  }
1602  if (!dof_is_fixed[j])
1603  freej++;
1604  }
1605  Fe(freei) += phi[i][qp] * fineval * JxW[qp];
1606  if (cont == C_ONE)
1607  Fe(freei) += (finegrad * (*dphi)[i][qp]) *
1608  JxW[qp];
1609  freei++;
1610  }
1611  }
1612 
1613  Ke.cholesky_solve(Fe, Uedge);
1614 
1615  // Transfer new edge solutions to element
1616  for (unsigned int i=0; i != free_dofs; ++i)
1617  {
1618  Number &ui = Ue(side_dofs[free_dof[i]]);
1620  std::abs(ui - Uedge(i)) < TOLERANCE);
1621  ui = Uedge(i);
1622  dof_is_fixed[side_dofs[free_dof[i]]] = true;
1623  }
1624  }
1625 
1626  // Project any side values (edges in 2D, faces in 3D)
1627  if (dim > 1 && cont != DISCONTINUOUS)
1628  for (unsigned int s=0; s != elem->n_sides(); ++s)
1629  {
1630  FEInterface::dofs_on_side(elem, dim, fe_type, s,
1631  side_dofs);
1632 
1633  // Some side dofs are on nodes/edges and already
1634  // fixed, others are free to calculate
1635  unsigned int free_dofs = 0;
1636  for (unsigned int i=0; i != side_dofs.size(); ++i)
1637  if (!dof_is_fixed[side_dofs[i]])
1638  free_dof[free_dofs++] = i;
1639 
1640  // There may be nothing to project
1641  if (!free_dofs)
1642  continue;
1643 
1644  Ke.resize (free_dofs, free_dofs); Ke.zero();
1645  Fe.resize (free_dofs); Fe.zero();
1646  // The new side coefficients
1647  DenseVector<Number> Uside(free_dofs);
1648 
1649  // Initialize FE data on the side
1650  fe->attach_quadrature_rule (qsiderule.get());
1651  fe->reinit (elem, s);
1652  const unsigned int n_qp = qsiderule->n_points();
1653 
1654  // Loop over the quadrature points
1655  for (unsigned int qp=0; qp<n_qp; qp++)
1656  {
1657  // solution at the quadrature point
1658  Number fineval = f->component(var_component,
1659  xyz_values[qp],
1660  system.time);
1661  // solution grad at the quadrature point
1662  Gradient finegrad;
1663  if (cont == C_ONE)
1664  finegrad = g->component(var_component,
1665  xyz_values[qp],
1666  system.time);
1667 
1668  // Form side projection matrix
1669  for (unsigned int sidei=0, freei=0;
1670  sidei != side_dofs.size(); ++sidei)
1671  {
1672  unsigned int i = side_dofs[sidei];
1673  // fixed DoFs aren't test functions
1674  if (dof_is_fixed[i])
1675  continue;
1676  for (unsigned int sidej=0, freej=0;
1677  sidej != side_dofs.size(); ++sidej)
1678  {
1679  unsigned int j = side_dofs[sidej];
1680  if (dof_is_fixed[j])
1681  Fe(freei) -= phi[i][qp] * phi[j][qp] *
1682  JxW[qp] * Ue(j);
1683  else
1684  Ke(freei,freej) += phi[i][qp] *
1685  phi[j][qp] * JxW[qp];
1686  if (cont == C_ONE)
1687  {
1688  if (dof_is_fixed[j])
1689  Fe(freei) -= ((*dphi)[i][qp] *
1690  (*dphi)[j][qp]) *
1691  JxW[qp] * Ue(j);
1692  else
1693  Ke(freei,freej) += ((*dphi)[i][qp] *
1694  (*dphi)[j][qp])
1695  * JxW[qp];
1696  }
1697  if (!dof_is_fixed[j])
1698  freej++;
1699  }
1700  Fe(freei) += (fineval * phi[i][qp]) * JxW[qp];
1701  if (cont == C_ONE)
1702  Fe(freei) += (finegrad * (*dphi)[i][qp]) *
1703  JxW[qp];
1704  freei++;
1705  }
1706  }
1707 
1708  Ke.cholesky_solve(Fe, Uside);
1709 
1710  // Transfer new side solutions to element
1711  for (unsigned int i=0; i != free_dofs; ++i)
1712  {
1713  Number &ui = Ue(side_dofs[free_dof[i]]);
1715  std::abs(ui - Uside(i)) < TOLERANCE);
1716  ui = Uside(i);
1717  dof_is_fixed[side_dofs[free_dof[i]]] = true;
1718  }
1719  }
1720 
1721  // Project the interior values, finally
1722 
1723  // Some interior dofs are on nodes/edges/sides and
1724  // already fixed, others are free to calculate
1725  unsigned int free_dofs = 0;
1726  for (unsigned int i=0; i != n_dofs; ++i)
1727  if (!dof_is_fixed[i])
1728  free_dof[free_dofs++] = i;
1729 
1730  // There may be nothing to project
1731  if (free_dofs)
1732  {
1733 
1734  Ke.resize (free_dofs, free_dofs); Ke.zero();
1735  Fe.resize (free_dofs); Fe.zero();
1736  // The new interior coefficients
1737  DenseVector<Number> Uint(free_dofs);
1738 
1739  // Initialize FE data
1740  fe->attach_quadrature_rule (qrule.get());
1741  fe->reinit (elem);
1742  const unsigned int n_qp = qrule->n_points();
1743 
1744  // Loop over the quadrature points
1745  for (unsigned int qp=0; qp<n_qp; qp++)
1746  {
1747  // solution at the quadrature point
1748  Number fineval = f->component(var_component,
1749  xyz_values[qp],
1750  system.time);
1751  // solution grad at the quadrature point
1752  Gradient finegrad;
1753  if (cont == C_ONE)
1754  finegrad = g->component(var_component,
1755  xyz_values[qp],
1756  system.time);
1757 
1758  // Form interior projection matrix
1759  for (unsigned int i=0, freei=0; i != n_dofs; ++i)
1760  {
1761  // fixed DoFs aren't test functions
1762  if (dof_is_fixed[i])
1763  continue;
1764  for (unsigned int j=0, freej=0; j != n_dofs; ++j)
1765  {
1766  if (dof_is_fixed[j])
1767  Fe(freei) -= phi[i][qp] * phi[j][qp] * JxW[qp]
1768  * Ue(j);
1769  else
1770  Ke(freei,freej) += phi[i][qp] * phi[j][qp] *
1771  JxW[qp];
1772  if (cont == C_ONE)
1773  {
1774  if (dof_is_fixed[j])
1775  Fe(freei) -= ((*dphi)[i][qp] *
1776  (*dphi)[j][qp]) * JxW[qp] *
1777  Ue(j);
1778  else
1779  Ke(freei,freej) += ((*dphi)[i][qp] *
1780  (*dphi)[j][qp]) *
1781  JxW[qp];
1782  }
1783  if (!dof_is_fixed[j])
1784  freej++;
1785  }
1786  Fe(freei) += phi[i][qp] * fineval * JxW[qp];
1787  if (cont == C_ONE)
1788  Fe(freei) += (finegrad * (*dphi)[i][qp]) * JxW[qp];
1789  freei++;
1790  }
1791  }
1792  Ke.cholesky_solve(Fe, Uint);
1793 
1794  // Transfer new interior solutions to element
1795  for (unsigned int i=0; i != free_dofs; ++i)
1796  {
1797  Number &ui = Ue(free_dof[i]);
1799  std::abs(ui - Uint(i)) < TOLERANCE);
1800  ui = Uint(i);
1801  dof_is_fixed[free_dof[i]] = true;
1802  }
1803 
1804  } // if there are free interior dofs
1805 
1806  // Make sure every DoF got reached!
1807  for (unsigned int i=0; i != n_dofs; ++i)
1808  libmesh_assert(dof_is_fixed[i]);
1809 
1810  const dof_id_type
1811  first = new_vector.first_local_index(),
1812  last = new_vector.last_local_index();
1813 
1814  // Lock the new_vector since it is shared among threads.
1815  {
1816  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1817 
1818  for (unsigned int i = 0; i < n_dofs; i++)
1819  // We may be projecting a new zero value onto
1820  // an old nonzero approximation - RHS
1821  // if (Ue(i) != 0.)
1822  if ((dof_indices[i] >= first) &&
1823  (dof_indices[i] < last))
1824  {
1825  new_vector.set(dof_indices[i], Ue(i));
1826  }
1827  }
1828  } // end elem loop
1829  } // end variables loop
1830 }

Member Data Documentation

AutoPtr<FunctionBase<Number> > libMesh::ProjectSolution::f
private

Definition at line 115 of file system_projection.C.

Referenced by ProjectSolution().

AutoPtr<FunctionBase<Gradient> > libMesh::ProjectSolution::g
private

Definition at line 116 of file system_projection.C.

Referenced by ProjectSolution().

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

Definition at line 118 of file system_projection.C.

const Parameters& libMesh::ProjectSolution::parameters
private

Definition at line 117 of file system_projection.C.

const System& libMesh::ProjectSolution::system
private

Definition at line 113 of file system_projection.C.

Referenced by operator()().


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

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:08 UTC

Hosted By:
SourceForge.net Logo