2
// approx P1 -> P0 and P2 -> P1d
3
// option -x i ==> "d_dxi" i=1,2,3
6
NAME: Computing the Gradient
8
\cindex{discontinuous approximation}
9
\cindex{block-diagonal matrix}
12
Remark that the gradient ${\bf q}_h = \nabla u_h$ of a continuous picewise linear function $u_h$
13
is a discontinuous piecewise constant function.
14
Conversely, the gradient of a continuous picewise quadratic
15
function is a discontinuous piecewise linear function.
16
\subsection{Formulation}
17
For all $i=1\ldots N$, the $i$-th component $q_{h,i}$ of the gradient belongs to the space:
19
T_h \ \{ \phi_h \in L^2(\Omega); \
20
\phi_{h/K} \in P_{k-1}, \
21
\forall K \in {\cal T}_h
24
By introducing the two following bilinear forms:
26
m_T(\psi_h,\phi_h) = \int_\Omega \psi_h \phi_h \, dx
29
b_i(u_h,\phi_h) = \int_\Omega \frac{\partial u_h}{\partial x_i} \phi_h \, dx
31
the gradient satisfies the following variational formulation:
33
m_T(q_{h,i},\phi_h) = b_i(u_h, \phi_h), \ \forall \psi_h \in Th
35
Remark that the matrix associated to the $m_T(.,,)$ bilinear form
36
is block-diagonal. Each block is associated to an element, and can
37
be inverted at the element level.
38
The following code uses this property on step (e),
39
as the \code{"inv\\_mass"} form.
43
#include "rheolef/rheolef.h"
46
string get_approx_grad (const space &Vh)
48
string Pk = Vh.get_approx();
49
if (Pk == "P1") return "P0";
50
if (Pk == "P2") return "P1d";
51
cerr << "unexpected approximation " << Pk << endl;
54
int main(int argc, char**argv)
58
space Vh = uh.get_space();
59
geo omega_h = Vh.get_geo();
60
string approx_grad = get_approx_grad(Vh);
61
space Th (omega_h, approx_grad);
62
form b(Vh, Th, "d_dx0");
63
form inv_mass (Th, Th, "inv_mass");
64
field du_dx0 = inv_mass*(b*uh);
65
int digits10 = numeric_limits<Float>::digits10;
66
cout << setprecision(digits10) << du_dx0;
added
removed
doc/usrman/d-dx.cc
