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
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
|
#include <assert.h>
#include "cctk.h"
#include "cctk_Parameters.h"
#include "gh.hh"
#include "carpet.hh"
#include "regrid.hh"
extern "C" {
static const char* rcsid = "$Header: /home/eschnett/C/carpet/Carpet/Carpet/CarpetRegrid/src/moving.cc,v 1.3 2004/04/18 13:29:43 schnetter Exp $";
CCTK_FILEVERSION(Carpet_CarpetRegrid_moving_cc);
}
namespace CarpetRegrid {
using namespace std;
using namespace Carpet;
int Moving (cGH const * const cctkGH,
gh<dim> const & hh,
gh<dim>::rexts & bbsss,
gh<dim>::rbnds & obss,
gh<dim>::rprocs & pss)
{
DECLARE_CCTK_ARGUMENTS;
DECLARE_CCTK_PARAMETERS;
assert (refinement_levels >= 1);
// do nothing if the levels already exist
if (reflevel == refinement_levels) return 0;
assert (bbsss.size() >= 1);
bbsss.resize (refinement_levels);
obss.resize (refinement_levels);
pss.resize (refinement_levels);
bvect const symmetric (symmetry_x, symmetry_y, symmetry_z);
ivect const zero(0);
ivect rstr = hh.baseextent.stride();
ivect rlb = hh.baseextent.lower();
ivect rub = hh.baseextent.upper();
for (size_t rl=1; rl<bbsss.size(); ++rl) {
// save old values
ivect const oldrlb = rlb;
ivect const oldrub = rub;
// refined boxes have smaller stride
assert (all(rstr%hh.reffact == 0));
rstr /= hh.reffact;
// calculate new extent
CCTK_REAL const argument = 2*M_PI * moving_circle_frequency * cctk_time;
rvect const pos
(moving_centre_x + moving_circle_radius * cos(argument),
moving_centre_y + moving_circle_radius * sin(argument),
moving_centre_z);
ivect const centre = floor(rvect(rub - rlb) * pos / rstr + 0.5) * rstr;
ivect const radius = floor(rvect(rub - rlb) * moving_region_radius / rstr + 0.5) * rstr;
rlb = oldrlb + symmetric.ifthen(zero , centre - radius);
rub = oldrlb + symmetric.ifthen(radius, centre + radius);
assert (all(rlb >= oldrlb && rub <= oldrub));
ibbox const bb (rlb, rub, rstr);
vector<ibbox> bbs (1);
bbs.at(0) = bb;
bbvect const ob (false);
gh<dim>::cbnds obs (1);
obs.at(0) = ob;
// make multiprocessor aware
gh<dim>::cprocs ps;
SplitRegions (cctkGH, bbs, obs, ps);
// make multigrid aware
vector<vector<ibbox> > bbss;
MakeMultigridBoxes (cctkGH, bbs, obs, bbss);
bbsss.at(rl) = bbss;
obss.at(rl) = obs;
pss.at(rl) = ps;
} // for rl
return 1;
}
} // namespace CarpetRegrid
|
