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
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
|
/*@@
@file RadiationBoundaryWrappers.c
@date Mon Mar 15 15:09:00 1999
@author Gerd Lanfermann
@desc
This file contains the C part of the BC routines. They are
called by C directly or by F via the wrappers below.
@enddesc
@@*/
/*#define DEBUG_BOUND*/
#include <stdio.h>
#include <assert.h>
#include <stdlib.h>
#include <ctype.h>
#include <stdarg.h>
/*#include <math.h>*/
#include <string.h>
#include "cctk.h"
#include "cctk_Flesh.h"
#include "cctk_parameters.h"
#include "cctk_Comm.h"
#include "cctk_Coord.h"
#include "cctk_Groups.h"
#include "cctk_GHExtensions.h"
#include "CactusBase/CartGrid3D/src/Symmetry.h"
#include "cctk_Misc.h"
#include "cctk_WarnLevel.h"
#include "cctk_FortranString.h"
/*@@
@routine ApplyRadiativeBC
@date Mon Mar 15 15:21:45 1999
@author Gerd Lanfermann
@desc
Radiative BC require GFs at past timesteps:
Call with implementation::group and implementation_group_previous
eg: ApplyFlatBC(admgerd::metric,admgerd::metric_previous)
It is not checked, that these two groups make sense!
Actual assigment is carried out in Fortran by FortranRadiativeBC
@enddesc
@calls
@calledby
@history
@endhistory
@@*/
int ApplyRadiativeBC(cGH *GH, CCTK_REAL var0, int *sw,
char *name, char *name_p) {
DECLARE_CCTK_PARAMETERS
SymmetryGHex *sGHex; /* the symmetry boundary GHextension, needed for ref*/
int xi,yi,zi,ri; /* index of the cartesian coordinate fields */
int first,last,index,first_p,last_p;
int num,num_p; /* index number of the group (not the gf)*/
int doBC[6]; /* flags if lower/upper BCs are applied (1) or not (0) */
/* indexing as in bbox: 0 xlow, 1 xup, 2 ylow, ...(C!) */
/* this has 3 dimension hardwired! will break for GH->dim>3 */
int type; /* 0 groups; 1 variables; -1 neither */
int retval; /* negative if boundary condition is not applied */
/* Get the pointer to the SymmetryGHextension */
sGHex = (SymmetryGHex*)malloc(sizeof(SymmetryGHex));
sGHex = (SymmetryGHex*)GH->extensions[CCTK_GHExtensionHandle("Symmetry")];
/* Decide if we have a group or a variable, and get the index */
/* type = 1 (group), type = 0 (var), type = -1 (neither) */
num = CCTK_GroupIndex(name);
if (num < 0)
{
num = CCTK_VarIndex(name);
if (num >= 0)
{
type = 1; /* Variable */
}
else
{
type = -1;
{
char *message;
message = (char *)malloc( 1024*sizeof(char) );
sprintf(message,"Name (%s) in ApplyRadiativeBC is not a group or variable",
name);
CCTK_WARN(1,message);
free(message);
retval = -1;
return retval;
}
}
}
else
{
type = 0; /* Group */
}
num_p = CCTK_GroupIndex(name_p);
if (num_p < 0)
{
num_p = CCTK_VarIndex(name_p);
if (num_p >= 0)
{
if (type != 1) /* Variable */
{
CCTK_WARN(1,"Group and Variable given in ApplyRadiativeBC");
retval = -1;
return retval;
}
}
else
{
char *message;
message = (char *)malloc( 1024*sizeof(char) );
sprintf(message,
"Name (%s) in ApplyRadiativeBC is not a group or variable", name_p);
CCTK_WARN(1,message);
free(message);
retval = -1;
return retval;
}
}
else
{
if (type != 0) /* Group */
{
CCTK_WARN(1,"Group and Variable given in ApplyRadiativeBC");
retval = -1;
return retval;
}
}
/*get the index (offset) of the first GF in the group and how many Vars there are*/
if (type == 0)
{
first = CCTK_FirstVarIndexI(num);
last = first+CCTK_NumVarsInGroupI(num)-1;
first_p = CCTK_FirstVarIndexI(num_p);
last_p = first_p+CCTK_NumVarsInGroupI(num_p)-1;
if (first == -1 || last == -1 || first_p == -1)
{
CCTK_WARN(1,"Invalid group number used");
retval = -1;
return retval;
}
if (last-first != last_p-first_p)
{
CCTK_WARN(1,"cctk_Groups.have different number of variables in ApplyRadiativeBC");
retval = -1;
return retval;
}
}
else
{
first = num;
last = num;
first_p = num_p;
last_p = num_p;
}
#ifdef DEBUG_BOUND
CCTK_PRINTSEPARATOR
printf("In RadiationBoundaryWrappers\n----------------------------\n");
if (type == 0)
{
printf("Group: %s, has indices %d-%d\n",CCTK_GroupName(num),first,last);
printf("Group: %s, has indices %d-%d\n",CCTK_GroupName(num_p),first_p,last_p);
}
else
{
printf("Variable: %s\n",CCTK_VarName(num));
printf("Variable: %s\n",CCTK_VarName(num_p));
}
#endif
/* Radiative boundaries need to know about the underlying Cartesian coordinates */
xi = CCTK_CoordIndex("x");
yi = CCTK_CoordIndex("y");
zi = CCTK_CoordIndex("z");
ri = CCTK_CoordIndex("r");
/* loop over the vars in the group (index + offset) */
for (index=0;index<last-first+1;index++) {
/* and check that we actually have a grid function (and not a scalar)*/
if (CCTK_GroupTypeFromVarI(first+index) == GROUP_GF) {
/* The rule is: IF we have no symmetries for the lower grid boundary,
(GFSym==GFSYM_NOSYM or ==GFSYM_UNSET) AND we have a lower(upper) bound
AND we have enough gridpoints THEN apply BC */
doBC[0]=(((sGHex->GFSym[first+index][0]==GFSYM_NOSYM)||
(sGHex->GFSym[first+index][0]==GFSYM_UNSET)) &&
GH->cctk_lsh[1]>1 && GH->cctk_bbox[0]);
doBC[2]=(((sGHex->GFSym[first+index][2]==GFSYM_NOSYM)||
(sGHex->GFSym[first+index][2]==GFSYM_UNSET)) &&
GH->cctk_lsh[1]>1 && GH->cctk_bbox[2]);
doBC[4]=(((sGHex->GFSym[first+index][4]==GFSYM_NOSYM)||
(sGHex->GFSym[first+index][4]==GFSYM_UNSET)) &&
GH->cctk_lsh[2]>1 && GH->cctk_bbox[4]);
doBC[1] = GH->cctk_lsh[0]>1 && GH->cctk_bbox[1];
doBC[3] = GH->cctk_lsh[1]>1 && GH->cctk_bbox[3];
doBC[5] = GH->cctk_lsh[2]>1 && GH->cctk_bbox[5];
#ifdef DEBUG_BOUND
printf("Applying radiative boundary condition to %s (%d)\n",CCTK_FullName(first+index),first+index);
printf("var0 is %f, sw is %d %d %d\n",var0,sw[0],sw[1],sw[2]);
#endif
/* Now make the call to F-boundary routine, passing all relevant information,
and let F do the physics */
FORTRAN_NAME(FortranRadiativeBC)(
(GH->cctk_lsh), /* yzx-size of PE local grid */
(GH->cctk_delta_space), /* grid spacing */
&(GH->cctk_delta_time), /* dt */
(GH->data[index+first] [0]), /* pointer to start of data array GF[] */
(GH->data[index+first_p][0]), /* pointer to start of prev data array */
GH->data[xi][0],GH->data[yi][0], /* pointer to the XYZR fields */
GH->data[zi][0],GH->data[ri][0],
sw, /* stencil_width */
(doBC), /* 6 dim. array to flag applying low/up BCs */
&(var0)); /* radiative BC limit ? */
}
}
#ifdef DEBUG_BOUND
CCTK_PRINTSEPARATOR
#endif
}
void FMODIFIER FORTRAN_NAME(ApplyRadiativeBC)(int *retval, cGH *GH,CCTK_REAL *var0,int *sw, TWO_FORTSTRINGS_ARGS) {
TWO_FORTSTRINGS_CREATE(name,name_p)
*retval = ApplyRadiativeBC(GH,*var0,sw,name,name_p);
free(name);
free(name_p);
}
|
