libMesh
element_postprocess.C
Go to the documentation of this file.
1 // General libMesh includes
2 #include "libmesh/libmesh_common.h"
3 #include "libmesh/elem.h"
4 #include "libmesh/fe_base.h"
5 #include "libmesh/fem_context.h"
6 #include "libmesh/point.h"
7 #include "libmesh/quadrature.h"
8 
9 // Local includes
10 #include "poisson.h"
11 
12 // Bring in everything from the libMesh namespace
13 using namespace libMesh;
14 
15 // Define the postprocess function to compute QoI 0, the weighted flux
16 
18 {
19  FEMContext & c = cast_ref<FEMContext &>(context);
20 
21  FEBase * elem_fe = NULL;
22  c.get_element_fe( 0, elem_fe );
23 
24  // Element Jacobian * quadrature weights for interior integration
25  const std::vector<Real> & JxW = elem_fe->get_JxW();
26 
27  const std::vector<Point> & xyz = elem_fe->get_xyz();
28 
29  // The number of local degrees of freedom in each variable
30 
31  unsigned int n_qpoints = c.get_element_qrule().n_points();
32 
33  Number dQoI_0 = 0.;
34 
35  // Loop over quadrature points
36 
37  for (unsigned int qp = 0; qp != n_qpoints; qp++)
38  {
39  // Get co-ordinate locations of the current quadrature point
40  const Real x = xyz[qp](0);
41  const Real y = xyz[qp](1);
42 
43  Real f = -alpha * ( ( (- 4 * alpha * alpha) * exp(-alpha*x) * y * (1 - y) ) + ( -8 + ( 8 * exp(-alpha*x) ) + ( 8 * ( 1 - exp(-alpha) )* x) ) );
44 
45  Real s = x * (1 - x) * (1 - y);
46  RealVectorValue U ((1-y)*(1-(2*x)),-x*(1-x));
47 
48  Gradient grad_u = c.interior_gradient(0, qp);
49 
50  // Flux with weight s = R(u^h, s) = int ( f*s - alpha*(grad_u*grad_s) ) dx
51  dQoI_0 += JxW[qp] * ( (f * s) - (alpha * (U * grad_u)) );
52  }
53 
54  // Update the computed value of the global functional R, by adding the contribution from this element
55 
56  computed_QoI[0] = computed_QoI[0] + dQoI_0;
57 
58 }
This class provides all data required for a physics package (e.g.
Definition: diff_context.h:54
void get_element_fe(unsigned int var, FEGenericBase< OutputShape > *&fe) const
Accessor for interior finite element object for variable var for the largest dimension in the mesh...
Definition: fem_context.h:262
The libMesh namespace provides an interface to certain functionality in the library.
PetscErrorCode Vec x
virtual void element_postprocess(DiffContext &context)
Does any work that needs to be done on elem in a postprocessing loop.
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:61
Gradient interior_gradient(unsigned int var, unsigned int qp) const
Definition: fem_context.C:381
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const QBase & get_element_qrule() const
Accessor for element interior quadrature rule for the dimension of the current _elem.
Definition: fem_context.h:765
unsigned int n_points() const
Definition: quadrature.h:113
This class forms the foundation from which generic finite elements may be derived.