/usr/share/doc/rheolef-doc/examples/navier_stokes_dg1.icc is in rheolef-doc 6.5-1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 | navier_stokes_dg::navier_stokes_dg (
Float Re1, const geo& omega, string approx)
: Re(Re1), Xh(), Qh(), qopt(), a0(), b(), c(), mu(), mp(), lh0(), lh(), kh(),
smu(), smp(), a1(), stokes1()
{
Xh = space (omega, approx, "vector");
Qh = space (omega, approx);
qopt.set_family(quadrature_option_type::gauss);
qopt.set_order(2*Xh.degree()+1);
stokes_dirichlet_dg (Xh, Qh, a0, b, c, mp, lh0, kh, qopt);
trial u (Xh); test v (Xh);
lh = lh0 + Re*inertia_fix_rhs (v, qopt);
mu = integrate (dot(u,v), qopt);
smu = solver(mu.uu());
smp = solver(mp.uu());
}
navier_stokes_dg::value_type
navier_stokes_dg::initial (string restart) const {
value_type xh(2);
xh[0] = field (Xh, 0);
xh[1] = field (Qh, 0);
Float Re0 = 0;
if (restart == "") {
solver_abtb stokes0 (a0.uu(), b.uu(), c.uu(), mp.uu());
stokes0.solve (lh0.u(), kh.u(), xh[0].set_u(), xh[1].set_u());
} else {
idiststream in (restart);
in >> catchmark("Re") >> Re0
>> catchmark("u") >> xh[0]
>> catchmark("p") >> xh[1];
check_macro (xh[1].get_space() == Qh, "unexpected " << xh[0].get_space().stamp()
<< " approximation in file \"" << restart << "\" (" << Xh.stamp() << " expected)");
}
derr << "# continuation: from Re=" << Re0 << " to " << Re << endl;
return xh;
}
navier_stokes_dg::value_type
navier_stokes_dg::residue (const value_type& xh) const {
trial u (Xh); test v (Xh);
form a = a0 + Re*inertia(xh[0], u, v, qopt);
value_type mrh(2);
mrh[0] = a*xh[0] + b.trans_mult(xh[1]) - lh;
mrh[1] = b*xh[0] - c*xh[1] - kh;
return mrh;
}
void navier_stokes_dg::update_derivative (const value_type& xh) const {
trial u (Xh); test v (Xh);
a1 = a0 + Re*(inertia(xh[0], u, v, qopt) + inertia(u, xh[0], v, qopt));
stokes1 = solver_abtb (a1.uu(), b.uu(), c.uu(), mp.uu());
}
navier_stokes_dg::value_type
navier_stokes_dg::derivative_solve (const value_type& mrh) const {
value_type delta_xh(2);
delta_xh[0] = field (Xh, 0);
delta_xh[1] = field (Qh, 0);
stokes1.solve (mrh[0].u(), mrh[1].u(),
delta_xh[0].set_u(), delta_xh[1].set_u());
return delta_xh;
}
navier_stokes_dg::value_type
navier_stokes_dg::derivative_trans_mult (const value_type& mrh) const {
value_type rh(2);
rh[0] = field (Xh);
rh[1] = field (Qh);
rh[0].set_u() = smu.solve(mrh[0].u());
rh[1].set_u() = smp.solve(mrh[1].u());
value_type mgh(2);
mgh[0] = a1.trans_mult(rh[0]) + b.trans_mult(rh[1]);
mgh[1] = b*rh[0] - c*rh[1];
return mgh;
}
|