← Back to branch summary

~ubuntu-branches/ubuntu/trusty/rheolef/trusty

« back to all changes in this revision

Viewing changes to doc/pexamples/neumann-laplace.cc

  • Committer: Package Import Robot
  • Author(s): Pierre Saramito
  • Date: 2012-04-06 09:12:21 UTC
  • mfrom: (1.1.5)
  • Revision ID: package-import@ubuntu.com-20120406091221-m58me99p1nxqui49
Tags: 6.0-1
http://bugs.debian.org/667356
http://bugs.debian.org/661689
* New upstream release 6.0 (major changes):
  - massively distributed and parallel support
  - full FEM characteristic method (Lagrange-Gakerkin method) support
  - enhanced users documentation 
  - source code supports g++-4.7 (closes: #667356)
* debian/control: dependencies for MPI distributed solvers added
* debian/rules: build commands simplified
* debian/librheolef-dev.install: man1/* to man9/* added
* debian/changelog: package description rewritted (closes: #661689)

Show diffs side-by-side

added added

removed removed

Lines of Context:
doc/pexamples/neumann-laplace.cc
1
1
#include "rheolef.h"
2
2
using namespace rheolef;
3
3
using namespace std;
4
 
#include "neumann-laplace-assembly.h"
5
 
size_t N;
 
4
size_t d;
6
5
Float f (const point& x) { return 1; }
7
 
Float g (const point& x) { return -0.5/N; }
 
6
Float g (const point& x) { return -0.5/d; }
8
7
int main(int argc, char**argv) {
9
 
  environment distributed (argc, argv);
 
8
  environment rheolef (argc, argv);
10
9
  geo omega (argv[1]);
11
 
  N = omega.dimension();
 
10
  d = omega.dimension();
12
11
  space Xh (omega, argv[2]);
13
12
  form m (Xh, Xh, "mass");
14
13
  form a (Xh, Xh, "grad_grad");
15
 
  field b = m*field(Xh,1.0);
16
 
  field fh = interpolate(Xh, f);
17
 
  space Wh (omega["boundary"], argv[2]);
18
 
  field gh = interpolate(Wh, g);
19
 
  form mb (Wh, Xh, "mass");
20
 
  field mfh = m*fh + mb*gh;
21
 
  distributor ext_ownership = mfh.u.ownership();
22
 
  communicator comm = ext_ownership.comm();
23
 
  size_t   my_proc = comm.rank();
24
 
  size_t last_proc = comm.size() - 1;
25
 
  ext_ownership [last_proc+1]++;
26
 
  csr<Float> A = neumann_laplace_assembly (ext_ownership, a.uu, b.u);
 
14
  field b = m*field(Xh,1);
 
15
  field lh = riesz(Xh, f) + riesz(Xh, "boundary", g);
 
16
  csr<Float> A = {{  a.uu(),     b.u()}, 
 
17
                  {trans(b.u()),  0 }};
 
18
  vec<Float> B =  {  lh.u(),      0 };
27
19
  solver sa (A);
28
 
  vec<Float> B (ext_ownership);
29
 
  for (size_t i = 0; i < mfh.u.size(); i++) B[i] = mfh.u[i];
30
 
  if (my_proc == last_proc) B [B.size()-1] = 0;
31
 
  vec<Float> U (ext_ownership);
32
 
  U = sa.solve (B);
33
 
  warning_macro ("retour de solve");
 
20
  vec<Float> U = sa.solve (B);
34
21
  field uh(Xh);
35
 
  warning_macro ("coucou {1}");
36
 
  for (size_t i = 0; i < uh.u.size(); i++) uh.u[i] = U[i];
37
 
  warning_macro ("coucou {2}");
38
 
  Float lambda = (my_proc == last_proc) ? U [U.size()-1] : 0;
39
 
  warning_macro ("coucou {3}");
 
22
  uh.set_u() = U [range(0,uh.u().size())];
 
23
  Float lambda = (U.size() == uh.u().size()+1) ? U [uh.u().size()] : 0;
40
24
#ifdef _RHEOLEF_HAVE_MPI
41
 
  mpi::broadcast (comm, lambda, last_proc); 
 
25
  mpi::broadcast (U.comm(), lambda, U.comm().size() - 1); 
42
26
#endif // _RHEOLEF_HAVE_MPI
43
 
  warning_macro ("coucou {4}");
44
 
  dcout << uh 
45
 
        << "lambda" << lambda << endl;
46
 
  warning_macro ("coucou {5}");
 
27
  dout << uh 
 
28
       << "lambda" << lambda << endl;
47
29
}