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
 
AutoPtr< FunctionBase< Number > > f
 
AutoPtr< 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 200 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 212 of file system_projection.C.

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

218  :
219  b(b_in),
220  variables(variables_in),
221  system(system_in),
222  f(f_in ? f_in->clone() : AutoPtr<FunctionBase<Number> >(NULL)),
223  g(g_in ? g_in->clone() : AutoPtr<FunctionBase<Gradient> >(NULL)),
224  parameters(parameters_in),
225  new_vector(new_v_in)
226  {
227  libmesh_assert(f.get());
228  f->init();
229  if (g.get())
230  g->init();
231  }
libMesh::BoundaryProjectSolution::BoundaryProjectSolution ( const BoundaryProjectSolution in)
inline

Definition at line 233 of file system_projection.C.

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

233  :
234  b(in.b),
235  variables(in.variables),
236  system(in.system),
237  f(in.f.get() ? in.f->clone() : AutoPtr<FunctionBase<Number> >(NULL)),
238  g(in.g.get() ? in.g->clone() : AutoPtr<FunctionBase<Gradient> >(NULL)),
239  parameters(in.parameters),
240  new_vector(in.new_vector)
241  {
242  libmesh_assert(f.get());
243  f->init();
244  if (g.get())
245  g->init();
246  }

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 2467 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::MeshBase::boundary_info, 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_edge_on_side(), libMesh::Elem::is_node_on_side(), 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::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().

2468 {
2469  // We need data to project
2470  libmesh_assert(f.get());
2471 
2479  // The dimensionality of the current mesh
2480  const unsigned int dim = system.get_mesh().mesh_dimension();
2481 
2482  // The DofMap for this system
2483  const DofMap& dof_map = system.get_dof_map();
2484 
2485  // Boundary info for the current mesh
2486  const BoundaryInfo& boundary_info = *system.get_mesh().boundary_info;
2487 
2488  // The element matrix and RHS for projections.
2489  // Note that Ke is always real-valued, whereas
2490  // Fe may be complex valued if complex number
2491  // support is enabled
2492  DenseMatrix<Real> Ke;
2493  DenseVector<Number> Fe;
2494  // The new element coefficients
2495  DenseVector<Number> Ue;
2496 
2497 
2498  // Loop over all the variables we've been requested to project
2499  for (unsigned int v=0; v!=variables.size(); v++)
2500  {
2501  const unsigned int var = variables[v];
2502 
2503  const Variable& variable = dof_map.variable(var);
2504 
2505  const FEType& fe_type = variable.type();
2506 
2507  if (fe_type.family == SCALAR)
2508  continue;
2509 
2510  const unsigned int var_component =
2512 
2513  // Get FE objects of the appropriate type
2514  AutoPtr<FEBase> fe (FEBase::build(dim, fe_type));
2515 
2516  // Prepare variables for projection
2517  AutoPtr<QBase> qedgerule (fe_type.default_quadrature_rule(1));
2518  AutoPtr<QBase> qsiderule (fe_type.default_quadrature_rule(dim-1));
2519 
2520  // The values of the shape functions at the quadrature
2521  // points
2522  const std::vector<std::vector<Real> >& phi = fe->get_phi();
2523 
2524  // The gradients of the shape functions at the quadrature
2525  // points on the child element.
2526  const std::vector<std::vector<RealGradient> > *dphi = NULL;
2527 
2528  const FEContinuity cont = fe->get_continuity();
2529 
2530  if (cont == C_ONE)
2531  {
2532  // We'll need gradient data for a C1 projection
2533  libmesh_assert(g.get());
2534 
2535  const std::vector<std::vector<RealGradient> >&
2536  ref_dphi = fe->get_dphi();
2537  dphi = &ref_dphi;
2538  }
2539 
2540  // The Jacobian * quadrature weight at the quadrature points
2541  const std::vector<Real>& JxW =
2542  fe->get_JxW();
2543 
2544  // The XYZ locations of the quadrature points
2545  const std::vector<Point>& xyz_values =
2546  fe->get_xyz();
2547 
2548  // The global DOF indices
2549  std::vector<dof_id_type> dof_indices;
2550  // Side/edge DOF indices
2551  std::vector<unsigned int> side_dofs;
2552 
2553  // Iterate over all the elements in the range
2554  for (ConstElemRange::const_iterator elem_it=range.begin(); elem_it != range.end(); ++elem_it)
2555  {
2556  const Elem* elem = *elem_it;
2557 
2558  // Per-subdomain variables don't need to be projected on
2559  // elements where they're not active
2560  if (!variable.active_on_subdomain(elem->subdomain_id()))
2561  continue;
2562 
2563  // Find out which nodes, edges and sides are on a requested
2564  // boundary:
2565  std::vector<bool> is_boundary_node(elem->n_nodes(), false),
2566  is_boundary_edge(elem->n_edges(), false),
2567  is_boundary_side(elem->n_sides(), false);
2568  for (unsigned char s=0; s != elem->n_sides(); ++s)
2569  {
2570  // First see if this side has been requested
2571  const std::vector<boundary_id_type>& bc_ids =
2572  boundary_info.boundary_ids (elem, s);
2573  bool do_this_side = false;
2574  for (unsigned int i=0; i != bc_ids.size(); ++i)
2575  if (b.count(bc_ids[i]))
2576  {
2577  do_this_side = true;
2578  break;
2579  }
2580  if (!do_this_side)
2581  continue;
2582 
2583  is_boundary_side[s] = true;
2584 
2585  // Then see what nodes and what edges are on it
2586  for (unsigned int n=0; n != elem->n_nodes(); ++n)
2587  if (elem->is_node_on_side(n,s))
2588  is_boundary_node[n] = true;
2589  for (unsigned int e=0; e != elem->n_edges(); ++e)
2590  if (elem->is_edge_on_side(e,s))
2591  is_boundary_edge[e] = true;
2592  }
2593 
2594  // Update the DOF indices for this element based on
2595  // the current mesh
2596  dof_map.dof_indices (elem, dof_indices, var);
2597 
2598  // The number of DOFs on the element
2599  const unsigned int n_dofs =
2600  libmesh_cast_int<unsigned int>(dof_indices.size());
2601 
2602  // Fixed vs. free DoFs on edge/face projections
2603  std::vector<char> dof_is_fixed(n_dofs, false); // bools
2604  std::vector<int> free_dof(n_dofs, 0);
2605 
2606  // The element type
2607  const ElemType elem_type = elem->type();
2608 
2609  // The number of nodes on the new element
2610  const unsigned int n_nodes = elem->n_nodes();
2611 
2612  // Zero the interpolated values
2613  Ue.resize (n_dofs); Ue.zero();
2614 
2615  // In general, we need a series of
2616  // projections to ensure a unique and continuous
2617  // solution. We start by interpolating boundary nodes, then
2618  // hold those fixed and project boundary edges, then hold
2619  // those fixed and project boundary faces,
2620 
2621  // Interpolate node values first
2622  unsigned int current_dof = 0;
2623  for (unsigned int n=0; n!= n_nodes; ++n)
2624  {
2625  // FIXME: this should go through the DofMap,
2626  // not duplicate dof_indices code badly!
2627  const unsigned int nc =
2628  FEInterface::n_dofs_at_node (dim, fe_type, elem_type,
2629  n);
2630  if (!elem->is_vertex(n) || !is_boundary_node[n])
2631  {
2632  current_dof += nc;
2633  continue;
2634  }
2635  if (cont == DISCONTINUOUS)
2636  {
2637  libmesh_assert_equal_to (nc, 0);
2638  }
2639  // Assume that C_ZERO elements have a single nodal
2640  // value shape function
2641  else if (cont == C_ZERO)
2642  {
2643  libmesh_assert_equal_to (nc, 1);
2644  Ue(current_dof) = f->component(var_component,
2645  elem->point(n),
2646  system.time);
2647  dof_is_fixed[current_dof] = true;
2648  current_dof++;
2649  }
2650  // The hermite element vertex shape functions are weird
2651  else if (fe_type.family == HERMITE)
2652  {
2653  Ue(current_dof) = f->component(var_component,
2654  elem->point(n),
2655  system.time);
2656  dof_is_fixed[current_dof] = true;
2657  current_dof++;
2658  Gradient grad = g->component(var_component,
2659  elem->point(n),
2660  system.time);
2661  // x derivative
2662  Ue(current_dof) = grad(0);
2663  dof_is_fixed[current_dof] = true;
2664  current_dof++;
2665  if (dim > 1)
2666  {
2667  // We'll finite difference mixed derivatives
2668  Point nxminus = elem->point(n),
2669  nxplus = elem->point(n);
2670  nxminus(0) -= TOLERANCE;
2671  nxplus(0) += TOLERANCE;
2672  Gradient gxminus = g->component(var_component,
2673  nxminus,
2674  system.time);
2675  Gradient gxplus = g->component(var_component,
2676  nxplus,
2677  system.time);
2678  // y derivative
2679  Ue(current_dof) = grad(1);
2680  dof_is_fixed[current_dof] = true;
2681  current_dof++;
2682  // xy derivative
2683  Ue(current_dof) = (gxplus(1) - gxminus(1))
2684  / 2. / TOLERANCE;
2685  dof_is_fixed[current_dof] = true;
2686  current_dof++;
2687 
2688  if (dim > 2)
2689  {
2690  // z derivative
2691  Ue(current_dof) = grad(2);
2692  dof_is_fixed[current_dof] = true;
2693  current_dof++;
2694  // xz derivative
2695  Ue(current_dof) = (gxplus(2) - gxminus(2))
2696  / 2. / TOLERANCE;
2697  dof_is_fixed[current_dof] = true;
2698  current_dof++;
2699  // We need new points for yz
2700  Point nyminus = elem->point(n),
2701  nyplus = elem->point(n);
2702  nyminus(1) -= TOLERANCE;
2703  nyplus(1) += TOLERANCE;
2704  Gradient gyminus = g->component(var_component,
2705  nyminus,
2706  system.time);
2707  Gradient gyplus = g->component(var_component,
2708  nyplus,
2709  system.time);
2710  // xz derivative
2711  Ue(current_dof) = (gyplus(2) - gyminus(2))
2712  / 2. / TOLERANCE;
2713  dof_is_fixed[current_dof] = true;
2714  current_dof++;
2715  // Getting a 2nd order xyz is more tedious
2716  Point nxmym = elem->point(n),
2717  nxmyp = elem->point(n),
2718  nxpym = elem->point(n),
2719  nxpyp = elem->point(n);
2720  nxmym(0) -= TOLERANCE;
2721  nxmym(1) -= TOLERANCE;
2722  nxmyp(0) -= TOLERANCE;
2723  nxmyp(1) += TOLERANCE;
2724  nxpym(0) += TOLERANCE;
2725  nxpym(1) -= TOLERANCE;
2726  nxpyp(0) += TOLERANCE;
2727  nxpyp(1) += TOLERANCE;
2728  Gradient gxmym = g->component(var_component,
2729  nxmym,
2730  system.time);
2731  Gradient gxmyp = g->component(var_component,
2732  nxmyp,
2733  system.time);
2734  Gradient gxpym = g->component(var_component,
2735  nxpym,
2736  system.time);
2737  Gradient gxpyp = g->component(var_component,
2738  nxpyp,
2739  system.time);
2740  Number gxzplus = (gxpyp(2) - gxmyp(2))
2741  / 2. / TOLERANCE;
2742  Number gxzminus = (gxpym(2) - gxmym(2))
2743  / 2. / TOLERANCE;
2744  // xyz derivative
2745  Ue(current_dof) = (gxzplus - gxzminus)
2746  / 2. / TOLERANCE;
2747  dof_is_fixed[current_dof] = true;
2748  current_dof++;
2749  }
2750  }
2751  }
2752  // Assume that other C_ONE elements have a single nodal
2753  // value shape function and nodal gradient component
2754  // shape functions
2755  else if (cont == C_ONE)
2756  {
2757  libmesh_assert_equal_to (nc, 1 + dim);
2758  Ue(current_dof) = f->component(var_component,
2759  elem->point(n),
2760  system.time);
2761  dof_is_fixed[current_dof] = true;
2762  current_dof++;
2763  Gradient grad = g->component(var_component,
2764  elem->point(n),
2765  system.time);
2766  for (unsigned int i=0; i!= dim; ++i)
2767  {
2768  Ue(current_dof) = grad(i);
2769  dof_is_fixed[current_dof] = true;
2770  current_dof++;
2771  }
2772  }
2773  else
2774  libmesh_error();
2775  }
2776 
2777  // In 3D, project any edge values next
2778  if (dim > 2 && cont != DISCONTINUOUS)
2779  for (unsigned int e=0; e != elem->n_edges(); ++e)
2780  {
2781  if (!is_boundary_edge[e])
2782  continue;
2783 
2784  FEInterface::dofs_on_edge(elem, dim, fe_type, e,
2785  side_dofs);
2786 
2787  // Some edge dofs are on nodes and already
2788  // fixed, others are free to calculate
2789  unsigned int free_dofs = 0;
2790  for (unsigned int i=0; i != side_dofs.size(); ++i)
2791  if (!dof_is_fixed[side_dofs[i]])
2792  free_dof[free_dofs++] = i;
2793 
2794  // There may be nothing to project
2795  if (!free_dofs)
2796  continue;
2797 
2798  Ke.resize (free_dofs, free_dofs); Ke.zero();
2799  Fe.resize (free_dofs); Fe.zero();
2800  // The new edge coefficients
2801  DenseVector<Number> Uedge(free_dofs);
2802 
2803  // Initialize FE data on the edge
2804  fe->attach_quadrature_rule (qedgerule.get());
2805  fe->edge_reinit (elem, e);
2806  const unsigned int n_qp = qedgerule->n_points();
2807 
2808  // Loop over the quadrature points
2809  for (unsigned int qp=0; qp<n_qp; qp++)
2810  {
2811  // solution at the quadrature point
2812  Number fineval = f->component(var_component,
2813  xyz_values[qp],
2814  system.time);
2815  // solution grad at the quadrature point
2816  Gradient finegrad;
2817  if (cont == C_ONE)
2818  finegrad = g->component(var_component,
2819  xyz_values[qp],
2820  system.time);
2821 
2822  // Form edge projection matrix
2823  for (unsigned int sidei=0, freei=0;
2824  sidei != side_dofs.size(); ++sidei)
2825  {
2826  unsigned int i = side_dofs[sidei];
2827  // fixed DoFs aren't test functions
2828  if (dof_is_fixed[i])
2829  continue;
2830  for (unsigned int sidej=0, freej=0;
2831  sidej != side_dofs.size(); ++sidej)
2832  {
2833  unsigned int j = side_dofs[sidej];
2834  if (dof_is_fixed[j])
2835  Fe(freei) -= phi[i][qp] * phi[j][qp] *
2836  JxW[qp] * Ue(j);
2837  else
2838  Ke(freei,freej) += phi[i][qp] *
2839  phi[j][qp] * JxW[qp];
2840  if (cont == C_ONE)
2841  {
2842  if (dof_is_fixed[j])
2843  Fe(freei) -= ((*dphi)[i][qp] *
2844  (*dphi)[j][qp]) *
2845  JxW[qp] * Ue(j);
2846  else
2847  Ke(freei,freej) += ((*dphi)[i][qp] *
2848  (*dphi)[j][qp])
2849  * JxW[qp];
2850  }
2851  if (!dof_is_fixed[j])
2852  freej++;
2853  }
2854  Fe(freei) += phi[i][qp] * fineval * JxW[qp];
2855  if (cont == C_ONE)
2856  Fe(freei) += (finegrad * (*dphi)[i][qp]) *
2857  JxW[qp];
2858  freei++;
2859  }
2860  }
2861 
2862  Ke.cholesky_solve(Fe, Uedge);
2863 
2864  // Transfer new edge solutions to element
2865  for (unsigned int i=0; i != free_dofs; ++i)
2866  {
2867  Number &ui = Ue(side_dofs[free_dof[i]]);
2869  std::abs(ui - Uedge(i)) < TOLERANCE);
2870  ui = Uedge(i);
2871  dof_is_fixed[side_dofs[free_dof[i]]] = true;
2872  }
2873  }
2874 
2875  // Project any side values (edges in 2D, faces in 3D)
2876  if (dim > 1 && cont != DISCONTINUOUS)
2877  for (unsigned int s=0; s != elem->n_sides(); ++s)
2878  {
2879  if (!is_boundary_side[s])
2880  continue;
2881 
2882  FEInterface::dofs_on_side(elem, dim, fe_type, s,
2883  side_dofs);
2884 
2885  // Some side dofs are on nodes/edges and already
2886  // fixed, others are free to calculate
2887  unsigned int free_dofs = 0;
2888  for (unsigned int i=0; i != side_dofs.size(); ++i)
2889  if (!dof_is_fixed[side_dofs[i]])
2890  free_dof[free_dofs++] = i;
2891 
2892  // There may be nothing to project
2893  if (!free_dofs)
2894  continue;
2895 
2896  Ke.resize (free_dofs, free_dofs); Ke.zero();
2897  Fe.resize (free_dofs); Fe.zero();
2898  // The new side coefficients
2899  DenseVector<Number> Uside(free_dofs);
2900 
2901  // Initialize FE data on the side
2902  fe->attach_quadrature_rule (qsiderule.get());
2903  fe->reinit (elem, s);
2904  const unsigned int n_qp = qsiderule->n_points();
2905 
2906  // Loop over the quadrature points
2907  for (unsigned int qp=0; qp<n_qp; qp++)
2908  {
2909  // solution at the quadrature point
2910  Number fineval = f->component(var_component,
2911  xyz_values[qp],
2912  system.time);
2913  // solution grad at the quadrature point
2914  Gradient finegrad;
2915  if (cont == C_ONE)
2916  finegrad = g->component(var_component,
2917  xyz_values[qp],
2918  system.time);
2919 
2920  // Form side projection matrix
2921  for (unsigned int sidei=0, freei=0;
2922  sidei != side_dofs.size(); ++sidei)
2923  {
2924  unsigned int i = side_dofs[sidei];
2925  // fixed DoFs aren't test functions
2926  if (dof_is_fixed[i])
2927  continue;
2928  for (unsigned int sidej=0, freej=0;
2929  sidej != side_dofs.size(); ++sidej)
2930  {
2931  unsigned int j = side_dofs[sidej];
2932  if (dof_is_fixed[j])
2933  Fe(freei) -= phi[i][qp] * phi[j][qp] *
2934  JxW[qp] * Ue(j);
2935  else
2936  Ke(freei,freej) += phi[i][qp] *
2937  phi[j][qp] * JxW[qp];
2938  if (cont == C_ONE)
2939  {
2940  if (dof_is_fixed[j])
2941  Fe(freei) -= ((*dphi)[i][qp] *
2942  (*dphi)[j][qp]) *
2943  JxW[qp] * Ue(j);
2944  else
2945  Ke(freei,freej) += ((*dphi)[i][qp] *
2946  (*dphi)[j][qp])
2947  * JxW[qp];
2948  }
2949  if (!dof_is_fixed[j])
2950  freej++;
2951  }
2952  Fe(freei) += (fineval * phi[i][qp]) * JxW[qp];
2953  if (cont == C_ONE)
2954  Fe(freei) += (finegrad * (*dphi)[i][qp]) *
2955  JxW[qp];
2956  freei++;
2957  }
2958  }
2959 
2960  Ke.cholesky_solve(Fe, Uside);
2961 
2962  // Transfer new side solutions to element
2963  for (unsigned int i=0; i != free_dofs; ++i)
2964  {
2965  Number &ui = Ue(side_dofs[free_dof[i]]);
2967  std::abs(ui - Uside(i)) < TOLERANCE);
2968  ui = Uside(i);
2969  dof_is_fixed[side_dofs[free_dof[i]]] = true;
2970  }
2971  }
2972 
2973  const dof_id_type
2974  first = new_vector.first_local_index(),
2975  last = new_vector.last_local_index();
2976 
2977  // Lock the new_vector since it is shared among threads.
2978  {
2979  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
2980 
2981  for (unsigned int i = 0; i < n_dofs; i++)
2982  if (dof_is_fixed[i] &&
2983  (dof_indices[i] >= first) &&
2984  (dof_indices[i] < last))
2985  {
2986  new_vector.set(dof_indices[i], Ue(i));
2987  }
2988  }
2989  } // end elem loop
2990  } // end variables loop
2991 }

Member Data Documentation

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

Definition at line 203 of file system_projection.C.

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

Definition at line 206 of file system_projection.C.

Referenced by BoundaryProjectSolution().

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

Definition at line 207 of file system_projection.C.

Referenced by BoundaryProjectSolution().

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

Definition at line 209 of file system_projection.C.

const Parameters& libMesh::BoundaryProjectSolution::parameters
private

Definition at line 208 of file system_projection.C.

const System& libMesh::BoundaryProjectSolution::system
private

Definition at line 205 of file system_projection.C.

Referenced by operator()().

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

Definition at line 204 of file system_projection.C.


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