/usr/share/doc/rheolef-doc/examples/laplace_band.cc 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 | #include "rheolef.h"
using namespace std;
using namespace rheolef;
#include "torus.icc"
int main (int argc, char**argv) {
environment rheolef(argc, argv);
geo lambda (argv[1]);
size_t d = lambda.dimension();
space Xh (lambda, "P1");
field phi_h = interpolate(Xh, phi);
band gamma_h (phi_h);
field phi_h_band = phi_h [gamma_h.band()];
space Bh (gamma_h.band(), "P1");
Bh.block ("isolated");
Bh.unblock ("zero");
trial u (Bh); test v (Bh);
form m = integrate (gamma_h, u*v);
form a = integrate (gamma_h, dot(grad_s(u),grad_s(v)));
field lh = integrate (gamma_h, f(d)*v);
vector<vec<Float> > b (gamma_h.n_connected_component());
vector<Float> z (gamma_h.n_connected_component(), 0);
for (size_t i = 0; i < b.size(); i++) {
const domain& cci = gamma_h.band() ["cc"+itos(i)];
field phi_h_cci (Bh, 0);
phi_h_cci [cci] = phi_h_band [cci];
b[i] = phi_h_cci.u();
}
field c = m*field(Bh,1);
csr<Float> A = { { a.uu(), trans(b), c.u()},
{ b, 0, 0 },
{ trans(c.u()), 0, 0 } };
vec<Float> F = { lh.u(), z, 0};
A.set_symmetry(true);
solver sa = ldlt(A);
vec<Float> U = sa.solve (F);
field uh(Bh,0);
uh.set_u() = U [range(0,uh.u().size())];
dout << catchmark("phi") << phi_h
<< catchmark("u") << uh;
}
|