2
#include "rheolef/catchmark.h"
3
using namespace rheolef;
5
int main(int argc, char**argv) {
6
environment rheolef (argc, argv);
7
const Float epsilon = 0.01;
9
space Xh (omega, "P1");
12
space Qh (omega, "P0");
15
eta ["west"] = epsilon;
16
field lh = riesz (Xh, field(Xh,1));
17
form b = -form(Xh, Qh, "grad");
18
form inv_m (Qh, Qh, "inv_mass", 1/eta);
19
form a = trans(b)*inv_m*b;
21
uh.u = fact.solve (lh.u - a.ub*uh.b);
22
cout << catchmark("epsilon") << epsilon << endl
23
<< catchmark("u") << uh;