-
Notifications
You must be signed in to change notification settings - Fork 38
/
GiRaFFE_NRPy_Metric_Face_Values.py
115 lines (102 loc) · 5.23 KB
/
GiRaFFE_NRPy_Metric_Face_Values.py
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
# Step 0: Add NRPy's directory to the path
# https://stackoverflow.com/questions/16780014/import-file-from-parent-directory
import os,sys
nrpy_dir_path = os.path.join("..")
if nrpy_dir_path not in sys.path:
sys.path.append(nrpy_dir_path)
from outputC import outCfunction, add_to_Cfunction_dict # NRPy+: Core C code output module
import cmdline_helper as cmd # NRPy+: Multi-platform Python command-line interface
prefunc = """// Side note: the following values could be used for cell averaged gfs:
// am2=-1.0/12.0, am1=7.0/12.0, a0=7.0/12.0, a1=-1.0/12.0
// However, since the metric gfs store the grid point values instead of the cell average,
// the following coefficients should be used:
// am2 = -1/16, am1 = 9/16, a0 = 9/16, a1 = -1/16
// This will yield the third-order-accurate face values at m-1/2,
// using values specified at {m-2,m-1,m,m+1}
#define AM2 -0.0625
#define AM1 0.5625
#define A0 0.5625
#define A1 -0.0625
#define COMPUTE_FCVAL(METRICm2,METRICm1,METRIC,METRICp1) (AM2*(METRICm2) + AM1*(METRICm1) + A0*(METRIC) + A1*(METRICp1))
const int metric_gfs_list[10] = {GAMMADD00GF,
GAMMADD01GF,
GAMMADD02GF,
GAMMADD11GF,
GAMMADD12GF,
GAMMADD22GF,
BETAU0GF,
BETAU1GF,
BETAU2GF,
ALPHAGF};
const int metric_gfs_face_list[10] = {GAMMA_FACEDD00GF,
GAMMA_FACEDD01GF,
GAMMA_FACEDD02GF,
GAMMA_FACEDD11GF,
GAMMA_FACEDD12GF,
GAMMA_FACEDD22GF,
BETA_FACEU0GF,
BETA_FACEU1GF,
BETA_FACEU2GF,
ALPHA_FACEGF};
const int num_metric_gfs = 10;
"""
def GiRaFFE_NRPy_FCVAL(Ccodesdir):
cmd.mkdir(Ccodesdir)
# Write out the code to a file.
with open(os.path.join(Ccodesdir,"interpolate_metric_gfs_to_cell_faces.h"),"w") as file:
file.write(prefunc)
desc = "Interpolate metric gridfunctions to cell faces"
name = "interpolate_metric_gfs_to_cell_faces"
interp_Cfunc = outCfunction(
outfile = "returnstring", desc=desc, name=name,
params ="const paramstruct *params,REAL *auxevol_gfs,const int flux_dirn",
preloop =""" int in_gf,out_gf;
REAL Qm2,Qm1,Qp0,Qp1;
""" ,
body =""" for(int gf = 0;gf < num_metric_gfs;gf++) {
in_gf = metric_gfs_list[gf];
out_gf = metric_gfs_face_list[gf];
for (int i2 = 2;i2 < Nxx_plus_2NGHOSTS2-1;i2++) {
for (int i1 = 2;i1 < Nxx_plus_2NGHOSTS1-1;i1++) {
for (int i0 = 2;i0 < Nxx_plus_2NGHOSTS0-1;i0++) {
Qm2 = auxevol_gfs[IDX4S(in_gf,i0-2*kronecker_delta[flux_dirn][0],i1-2*kronecker_delta[flux_dirn][1],i2-2*kronecker_delta[flux_dirn][2])];
Qm1 = auxevol_gfs[IDX4S(in_gf,i0-kronecker_delta[flux_dirn][0],i1-kronecker_delta[flux_dirn][1],i2-kronecker_delta[flux_dirn][2])];
Qp0 = auxevol_gfs[IDX4S(in_gf,i0,i1,i2)];
Qp1 = auxevol_gfs[IDX4S(in_gf,i0+kronecker_delta[flux_dirn][0],i1+kronecker_delta[flux_dirn][1],i2+kronecker_delta[flux_dirn][2])];
auxevol_gfs[IDX4S(out_gf,i0,i1,i2)] = COMPUTE_FCVAL(Qm2,Qm1,Qp0,Qp1);
}
}
}
}
""",
rel_path_to_Cparams=os.path.join("../"))
with open(os.path.join(Ccodesdir,"interpolate_metric_gfs_to_cell_faces.h"),"a") as file:
file.write(interp_Cfunc)
def add_to_Cfunction_dict__GiRaFFE_NRPy_FCVAL(includes=None):
desc = "Interpolate metric gridfunctions to cell faces"
name = "interpolate_metric_gfs_to_cell_faces"
params ="const paramstruct *params,REAL *auxevol_gfs,const int flux_dirn"
preloop =""" int in_gf,out_gf;
REAL Qm2,Qm1,Qp0,Qp1;
"""
body =""" for(int gf = 0;gf < num_metric_gfs;gf++) {
in_gf = metric_gfs_list[gf];
out_gf = metric_gfs_face_list[gf];
for (int i2 = 2;i2 < Nxx_plus_2NGHOSTS2-1;i2++) {
for (int i1 = 2;i1 < Nxx_plus_2NGHOSTS1-1;i1++) {
for (int i0 = 2;i0 < Nxx_plus_2NGHOSTS0-1;i0++) {
Qm2 = auxevol_gfs[IDX4S(in_gf,i0-2*kronecker_delta[flux_dirn][0],i1-2*kronecker_delta[flux_dirn][1],i2-2*kronecker_delta[flux_dirn][2])];
Qm1 = auxevol_gfs[IDX4S(in_gf,i0-kronecker_delta[flux_dirn][0],i1-kronecker_delta[flux_dirn][1],i2-kronecker_delta[flux_dirn][2])];
Qp0 = auxevol_gfs[IDX4S(in_gf,i0,i1,i2)];
Qp1 = auxevol_gfs[IDX4S(in_gf,i0+kronecker_delta[flux_dirn][0],i1+kronecker_delta[flux_dirn][1],i2+kronecker_delta[flux_dirn][2])];
auxevol_gfs[IDX4S(out_gf,i0,i1,i2)] = COMPUTE_FCVAL(Qm2,Qm1,Qp0,Qp1);
}
}
}
}
"""
add_to_Cfunction_dict(
includes=includes,
desc=desc,
name=name, params=params,
prefunc = prefunc, preloop = preloop, body=body)