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
| #include "poisson.h"
#include "utils.h"
int main()
{
size (1.[0]);
init_grid (N);
scalar a[], b[], e[];
#if 0
// this does not work with multigrid MPI or trees
foreach_dimension() {
a[right] = periodic();
}
#else
foreach_dimension()
periodic (right);
#endif
for (int n = 8; n <= 64; n *= 2) {
init_grid (n);
foreach() {
a[] = 0.;
b[] = 4.*dimension*sq(pi)*sin(2.*pi*x)*cos(2.*pi*y)*cos(2.*pi*z);
}
TOLERANCE = 1e-10;
poisson (a, b);
foreach()
e[] = a[] + sin(2.*pi*x)*cos(2.*pi*y)*cos(2.*pi*z);
stats s = statsf(e);
foreach()
e[] -= s.sum/s.volume;
fprintf (stderr, "%d %g\n", n, statsf(e).max);
}
foreach()
printf ("%g %g %g %g %g\n", x, y, z, a[], e[]);
}
|