Peano 4
Loading...
Searching...
No Matches
FiniteVolumeCCZ4.cpp
Go to the documentation of this file.
1#include "FiniteVolumeCCZ4.h"
3//#include "exahype2/CellAccess.h"
4
5#include <algorithm>
6
7#include "Constants.h"
8
9#include <limits>
10
11#include <stdio.h>
12#include <string.h>
14
15#ifdef IncludeTwoPunctures
17#endif
18
19tarch::logging::Log applications::exahype2::ccz4::FiniteVolumeCCZ4::_log( "applications::exahype2::ccz4::FiniteVolumeCCZ4" );
20
21#ifdef IncludeTwoPunctures
22//pre-process, solve the puncture equations
23void applications::exahype2::ccz4::FiniteVolumeCCZ4::prepare(TP::TwoPunctures* tp){
24 //first we set the parameter. TODO:find a way to read parameter from python script
25 //int swi=0;//0--single black hole, 2--BBH hoc, 2--BBH rotation, 3--GW150914
26 //bbhtype=0 head on, 1-b=3, 2-b=2
27
28 if (CCZ4swi==0){
29 tp->par_b=1.0;
30 tp->center_offset[0]=-1.0; tp->center_offset[1]=0.0; tp->center_offset[2]=0.0;
31 tp->target_M_plus=1.0;//adm mass
32 tp->par_P_plus[0]=0.0; tp->par_P_plus[1]=0.0; tp->par_P_plus[2]=0.0;//linear momentum
33 tp->par_S_plus[0]=0.0; tp->par_S_plus[1]=0.0; tp->par_S_plus[2]=0.0;//spin
34 tp->target_M_minus=0.0;//adm mass
35 tp->par_P_minus[0]=0.0; tp->par_P_minus[1]=0.0; tp->par_P_minus[2]=0.0;//linear momentum
36 tp->par_S_minus[0]=0.0; tp->par_S_minus[1]=0.0; tp->par_S_minus[2]=0.0; //spin
37 tp->grid_setup_method="evaluation"; //evaluation or Taylor expansion
38 tp->TP_epsilon=1e-6;}
39
40 if (CCZ4swi==4){
41 tp->par_b=1.0;
42 tp->center_offset[0]=-3.0; tp->center_offset[1]=-2.0; tp->center_offset[2]=0.0;
43 tp->target_M_plus=1.0;//adm mass
44 tp->par_P_plus[0]=0.1; tp->par_P_plus[1]=0.1; tp->par_P_plus[2]=0.0;//linear momentum
45 tp->par_S_plus[0]=0.0; tp->par_S_plus[1]=0.0; tp->par_S_plus[2]=0.0;//spin
46 tp->target_M_minus=0.0;//adm mass
47 tp->par_P_minus[0]=0.0; tp->par_P_minus[1]=0.0; tp->par_P_minus[2]=0.0;//linear momentum
48 tp->par_S_minus[0]=0.0; tp->par_S_minus[1]=0.0; tp->par_S_minus[2]=0.0; //spin
49 tp->grid_setup_method="evaluation"; //evaluation or Taylor expansion
50 tp->TP_epsilon=1e-6;}
51
52 if (CCZ4swi==2 and CCZ4BBHType==0){
53 tp->par_b=2.0;
54 tp->center_offset[0]=0.0; tp->center_offset[1]=0.0; tp->center_offset[2]=0.0;
55 tp->target_M_plus=0.5;//adm mass
56 tp->par_P_plus[0]=0.0; tp->par_P_plus[1]=0.0; tp->par_P_plus[2]=0.0;//linear momentum
57 tp->par_S_plus[0]=0.0; tp->par_S_plus[1]=0.0; tp->par_S_plus[2]=0.0;//spin
58 tp->target_M_minus=0.5;//adm mass
59 tp->par_P_minus[0]=0.0; tp->par_P_minus[1]=0.0; tp->par_P_minus[2]=0.0;//linear momentum
60 tp->par_S_minus[0]=0.0; tp->par_S_minus[1]=0.0; tp->par_S_minus[2]=0.0; //spin
61 tp->grid_setup_method="evaluation"; //evaluation or Taylor expansion
62 tp->TP_epsilon=1e-6;}
63
64 //following data reference: https://journals.aps.org/prd/pdf/10.1103/PhysRevD.69.024006
65 if (CCZ4swi==2 and CCZ4BBHType==1){
66 tp->par_b=3.0;
67 tp->center_offset[0]=0.0; tp->center_offset[1]=0.0; tp->center_offset[2]=0.0;
68 tp->give_bare_mass=true;//use puncture mass instead of adm mass
69 tp->par_m_plus=0.47656; tp->par_m_minus=0.47656;
70 //tp->target_M_plus=999;//adm mass
71 tp->par_P_plus[0]=0.0; tp->par_P_plus[1]=0.13808; tp->par_P_plus[2]=0.0;//linear momentum
72 tp->par_S_plus[0]=0.0; tp->par_S_plus[1]=0.0; tp->par_S_plus[2]=0.0;//spin
73 //tp->target_M_minus=999;//adm mass
74 tp->par_P_minus[0]=0.0; tp->par_P_minus[1]=-0.13808; tp->par_P_minus[2]=0.0;//linear momentum
75 tp->par_S_minus[0]=0.0; tp->par_S_minus[1]=0.0; tp->par_S_minus[2]=0.0; //spin
76 tp->grid_setup_method="evaluation"; //evaluation or Taylor expansion
77 tp->TP_epsilon=1e-6;}
78
79 if (CCZ4swi==2 and CCZ4BBHType==2){
80 tp->par_b=2.0;
81 tp->center_offset[0]=0.0; tp->center_offset[1]=0.0; tp->center_offset[2]=0.0;
82 tp->give_bare_mass=true;//use puncture mass instead of adm mass
83 tp->par_m_plus=0.46477; tp->par_m_minus=0.46477;
84 //tp->target_M_plus=999;//adm mass
85 tp->par_P_plus[0]=0.0; tp->par_P_plus[1]=0.19243; tp->par_P_plus[2]=0.0;//linear momentum
86 tp->par_S_plus[0]=0.0; tp->par_S_plus[1]=0.0; tp->par_S_plus[2]=0.0;//spin
87 //tp->target_M_minus=999;//adm mass
88 tp->par_P_minus[0]=0.0; tp->par_P_minus[1]=-0.19243; tp->par_P_minus[2]=0.0;//linear momentum
89 tp->par_S_minus[0]=0.0; tp->par_S_minus[1]=0.0; tp->par_S_minus[2]=0.0; //spin
90 tp->grid_setup_method="evaluation"; //evaluation or Taylor expansion
91 tp->TP_epsilon=1e-6;}
92
93 if (CCZ4swi==3){
94 double D=10.0, q=36.0/29.0, chip=0.31, chim=-0.46, M=1.0;
95 double Pr=-0.00084541526517121, Pphi=0.09530152296974252;
96 double mp=M*q/(1+q), mm=M*1/(1+q);
97 tp->par_b=5.0;
98 tp->center_offset[0]=D*mm-D/2; tp->center_offset[1]=0.0; tp->center_offset[2]=0.0;
99 tp->target_M_plus=mp;//adm mass
100 tp->par_P_plus[0]=Pr; tp->par_P_plus[1]=Pphi; tp->par_P_plus[2]=0.0;//linear momentum
101 tp->par_S_plus[0]=0.0; tp->par_S_plus[1]=0.0; tp->par_S_plus[2]=chip*mp*mp;//spin
102 tp->target_M_minus=1/(1+q);//adm mass
103 tp->par_P_minus[0]=-Pr; tp->par_P_minus[1]=-Pphi; tp->par_P_minus[2]=0.0;//linear momentum
104 tp->par_S_minus[0]=0.0; tp->par_S_minus[1]=0.0; tp->par_S_minus[2]=chim*mm*mm; //spin
105 tp->grid_setup_method="evaluation"; //evaluation or Taylor expansion
106 tp->TP_epsilon=1e-6;}
107 tp->PrintParameters();
108
109 //then solve the equation
110 tp->Run();
111}
112#endif
113
115 if ( Scenario==0 || Scenario==1 ) {
116 const char* name = "GaugeWave";//nothing to do here for now
117 int length = strlen(name);
118 //initparameters_(&length, name);
119 }
120 #ifdef IncludeTwoPunctures
121 if ( Scenario==2 ) {
122 prepare(_tp);//we solve the puncture equation here.
123 //exit(0);
124 }
125 #endif
126 else {
127 //std::cerr << "initial scenario " << Scenario << " is not supported" << std::endl << std::endl << std::endl;
128 }
129}
130
132 double * __restrict__ Q,
135 bool gridIsConstructed
136) {
137 logTraceInWith3Arguments( "initialCondition(...)", volumeX, volumeH, gridIsConstructed );
138
139 if ( Scenario==0 ) {
141 }
142 else if ( Scenario==3 ) {
144 }
145 else if ( Scenario==1 ) {
147 }
148 else if ( Scenario==4 ) {
150 }
151 #ifdef IncludeTwoPunctures
152 else if ( Scenario==2 ) {
153
154 // We use the bool to trigger the hgh res interpolation once the grid is constructed
155 applications::exahype2::ccz4::ApplyTwoPunctures(Q, volumeX, 0, _tp, not gridIsConstructed); //we interpolate for real IC here.
156 }
157 #endif
158 else {
159 logError( "initialCondition(...)", "initial scenario " << Scenario << " is not supported" );
160 }
161
162 for (int i=0; i<NumberOfUnknowns; i++) {
163 assertion2( std::isfinite(Q[i]), volumeX, i );
164 }
165
166 for (int i=NumberOfUnknowns; i<NumberOfUnknowns+NumberOfAuxiliaryVariables; i++) {
167 Q[i] = 0.0;
168 }
169
170
171/*
172 else {
173 enforceCCZ4constraints(Q);
174 }
175*/
176 logTraceOut( "initialCondition(...)" );
177}
178
179
181 const double * __restrict__ Q, // Q[59+0]
184 double t,
185 double dt,
186 double * __restrict__ S // S[59
187) {
188#if !defined(GPUOffloadingOMP)
189 logTraceInWith4Arguments( "sourceTerm(...)", volumeX, volumeH, t, dt );
190 for(int i=0; i<NumberOfUnknowns; i++){
191 assertion3( std::isfinite(Q[i]), i, volumeX, t );
192 }
193#endif
194 // for (int i=NumberOfUnknowns; i<NumberOfUnknowns+NumberOfAuxiliaryVariables; i++) {
195 // S[i] = 0.0;
196 //}
197 applications::exahype2::ccz4::FiniteVolumeCCZ4::sourceTerm(Q, volumeX, volumeH, t, dt, S, Offloadable::Yes);
198 /*
199 if (CCZ4smoothing>0){// then we apply simple laplacian smoothing
200 constexpr int NumberOfRefinementLayers = 3;
201 double Radius[NumberOfRefinementLayers] = {5.0*1.8, 3.0*1.8, 1.5*1.8};
202 double radius=volumeX(0)*volumeX(0)+volumeX(1)*volumeX(1)+volumeX(2)*volumeX(2); radius=pow(radius,0.5);
203 ::exahype2::CellAccess access(
204 Q, // make it initially point to input
205 1, // halo size, which you know as user
206 NumberOfUnknowns, // defined in superclass
207 NumberOfAuxiliaryVariables // defined in superclass
208 NumberOfDoFsPerAxisInCell, // defined in superclass
209 );
210 double laplacian=0.0;
211 for (int unknown=0; unknown<NumberOfUnknowns; unknown++) {
212 laplacian = -6.0 * access.centre(unknown)
213 + 1.0 * access.left(0,unknown)+ 1.0 * access.right(0,unknown)
214 + 1.0 * access.left(1,unknown)+ 1.0 * access.right(1,unknown)
215 + 1.0 * access.left(2,unknown)+ 1.0 * access.right(2,unknown);
216 laplacian /= volumeH(0);
217 laplacian /= volumeH(0);
218 double coef=CCZ4smoothing*(std::exp(-5.0*std::abs(radius-Radius[0])) + std::exp(-5.0*std::abs(radius-Radius[1])) + std::exp(-5.0*std::abs(radius-Radius[2])) )/3.0;
219 coef=CCZ4smoothing;
220 S[unknown]+=coef*laplacian;
221 }
222 }*/
223
224#if !defined(GPUOffloadingOMP)
225 for(int i=0; i<NumberOfUnknowns; i++){
226 nonCriticalAssertion3( std::isfinite(S[i]), i, volumeX, t );
227 }
228 logTraceOut( "sourceTerm(...)" );
229#endif
230}
231
232
234 const double * __restrict__ Qinside, // Qinside[59+0]
235 double * __restrict__ Qoutside, // Qoutside[59+0]
236 const tarch::la::Vector<Dimensions,double>& faceCentre,
238 double t,
239 int normal
240) {
241 logTraceInWith4Arguments( "boundaryConditions(...)", faceCentre, volumeH, t, normal );
242 for(int i=0; i<NumberOfUnknowns+NumberOfAuxiliaryVariables; i++) {
243 assertion4( Qinside[i]==Qinside[i], faceCentre, t, normal, i );
244 Qoutside[i]=Qinside[i];
245 }
246 logTraceOut( "boundaryConditions(...)" );
247}
248
249
251 const double * __restrict__ Q, // Q[59+0],
252 const tarch::la::Vector<Dimensions,double>& faceCentre,
254 double t,
255 double dt,
256 int normal
257)
258{
259 return maxEigenvalue(Q, faceCentre, volumeH, t, dt, normal, Offloadable::Yes);
260}
261
262
264 const double * __restrict__ Q, // Q[59+0],
265 const double * __restrict__ deltaQ, // [59+0]
266 const tarch::la::Vector<Dimensions,double>& faceCentre,
268 double t,
269 double dt,
270 int normal,
271 double * __restrict__ BgradQ // BgradQ[59]
272) {
273#if !defined(GPUOffloadingOMP)
274 logTraceInWith4Arguments( "nonconservativeProduct(...)", faceCentre, volumeH, t, normal );
275 assertion( normal>=0 );
276 assertion( normal<Dimensions );
277#endif
278 nonconservativeProduct(Q, deltaQ, faceCentre, volumeH, t, dt, normal, BgradQ, Offloadable::Yes);
279
280#if !defined(GPUOffloadingOMP)
281 for (int i=0; i<NumberOfUnknowns; i++) {
282 nonCriticalAssertion4( std::isfinite(BgradQ[i]), i, faceCentre, t, normal );
283 }
284 logTraceOut( "nonconservativeProduct(...)" );
285#endif
286}
287
288#if defined(GPUOffloadingOMP)
289#pragma omp declare target
290#endif
292 const double * __restrict__ Q, // Q[59+0],
293 const double * __restrict__ deltaQ, // [59+0]
294 const tarch::la::Vector<Dimensions,double>& faceCentre,
296 double t,
297 double dt,
298 int normal,
299 double * __restrict__ BgradQ, // BgradQ[59]
300 Offloadable
301)
302{
303 double gradQSerialised[NumberOfUnknowns*3];
304 for (int i=0; i<NumberOfUnknowns; i++) {
305 gradQSerialised[i+0*NumberOfUnknowns] = 0.0;
306 gradQSerialised[i+1*NumberOfUnknowns] = 0.0;
307 gradQSerialised[i+2*NumberOfUnknowns] = 0.0;
308
309 gradQSerialised[i+normal*NumberOfUnknowns] = deltaQ[i];
310 }
311 //for (int i=NumberOfUnknowns; i<NumberOfUnknowns+NumberOfAuxiliaryVariables; i++) {
312 // BgradQ[i] = 0.0;
313 //}
314 ncp(BgradQ, Q, gradQSerialised, normal, CCZ4LapseType, CCZ4ds, CCZ4c, CCZ4e, CCZ4f, CCZ4bs, CCZ4sk, CCZ4xi, CCZ4mu);
315 //applications::exahype2::ccz4::FiniteVolumeCCZ4::nonconservativeProduct(Q, deltaQ, faceCentre, volumeH, t, normal, BgradQ);
316}
317#if defined(GPUOffloadingOMP)
318#pragma omp end declare target
319#endif
320
322 const double * __restrict__ Q,
323 const tarch::la::Vector<Dimensions,double>& volumeCentre,
325 double t
326) {
328
329 const double radius = tarch::la::norm2( volumeCentre );
330 //
331 // see documentation in header file
332 //
333 if (CCZ4ReSwi==1) { //radius based
334 if (radius<0.1) {
336 }
337 }
338 if (CCZ4ReSwi==2) { //
339 if (tarch::la::equals(t,0.0)){ //as we use a quantity calculated in postpocessing, we need to provide criterion at the first timestep
340 constexpr int NumberOfRefinementLayers = 2;
341// current "standard" refinement pattern
342 double Radius[NumberOfRefinementLayers] = {4, 2.5};
343 double MaxH[NumberOfRefinementLayers] = {0.15,0.04};
345 for (int i=0; i<NumberOfRefinementLayers; i++) {
346 if (radius<Radius[i] and tarch::la::max(volumeH)>MaxH[i]) {
348 }
349 }
350 }
351 }
352 if (CCZ4ReSwi==3) { //binary black holes
353 //double radius1=(volumeCentre(0)-4.251)*(volumeCentre(0)-4.251)+volumeCentre(1)*volumeCentre(1)+volumeCentre(2)*volumeCentre(2);
354 //double radius2=(volumeCentre(0)+4.251)*(volumeCentre(0)+4.251)+volumeCentre(1)*volumeCentre(1)+volumeCentre(2)*volumeCentre(2);
355
356 tarch::la::Vector<Dimensions,double> leftBH = {-4.241, 0.0, 0.0};
357 tarch::la::Vector<Dimensions,double> rightBH = { 4.241, 0.0, 0.0};
358
359 const double radius1 = tarch::la::norm2( volumeCentre - leftBH );
360 const double radius2 = tarch::la::norm2( volumeCentre - rightBH );
361
362 if (tarch::la::equals(t,0.0)){ //as we use a quantity calculated in postpocessing, we need to provide criterion at the first timestep
363 if ( ((radius1<5) or (radius2<5)) and (volumeH(0)>1.0)) { result=::exahype2::RefinementCommand::Refine; }
364 else if ((radius1<2.5) or (radius2<2.5)) { result=::exahype2::RefinementCommand::Refine; }
366 } else {
367 if ((Q[65]>0.1) and (volumeH(0)>1.0)) { result=::exahype2::RefinementCommand::Refine; }
368 else if (Q[65]>0.2) { result=::exahype2::RefinementCommand::Refine; }
370 }
371 }
372 if (CCZ4ReSwi==4){ //single black hole, higher level of amr
373 if (tarch::la::equals(t,0.0)){ //as we use a quantity calculated in postpocessing, we need to provide criterion at the first timestep
374 constexpr int NumberOfRefinementLayers = 4;
375 // with some alter
376 double Radius[NumberOfRefinementLayers] = {5.0, 3.0, 1.5, 0.5};
377 double MaxH[NumberOfRefinementLayers] = {0.3, 0.15, 0.04, 0.01};
378
380 for (int i=0; i<NumberOfRefinementLayers; i++) {
381 if (radius<Radius[i] and tarch::la::max(volumeH)>MaxH[i]) {
383 }
384 }
385 }
386 }
387 if (CCZ4ReSwi==8){ //single black hole, higher level of amr
388 if (tarch::la::equals(t,0.0)){ //as we use a quantity calculated in postpocessing, we need to provide criterion at the first timestep
389 constexpr int NumberOfRefinementLayers = 1;
390 double Radius[NumberOfRefinementLayers] = {2};
391 double MaxH[NumberOfRefinementLayers] = {0.04};
392
394 for (int i=0; i<NumberOfRefinementLayers; i++) {
395 if (radius<Radius[i] and tarch::la::max(volumeH)>MaxH[i]) {
397 }
398 }
399 }
400 }
401 if (CCZ4ReSwi==5){ //single black hole, decent refinement
402 if (tarch::la::equals(t,0.0)){ //as we use a quantity calculated in postpocessing, we need to provide criterion at the first timestep
403 constexpr int NumberOfRefinementLayers = 3;
404// current "standard" refinement pattern
405 double Radius[NumberOfRefinementLayers] = {3.0, 1.0, 0.5};
406 double MaxH[NumberOfRefinementLayers] = {0.2, 0.06, 0.02};
407// with some alter
408// double Radius[NumberOfRefinementLayers] = {5.0, 3.0};
409// double MaxH[NumberOfRefinementLayers] = {0.3, 0.15};
410
412 for (int i=0; i<NumberOfRefinementLayers; i++) {
413 if (radius<Radius[i] and tarch::la::max(volumeH)>MaxH[i]) {
415 }
416 }
417 }
418 }
419 if (CCZ4ReSwi==6){
420 if (tarch::la::equals(t,0.0)){ //as we use a quantity calculated in postpocessing, we need to provide criterion at the first timestep
421 constexpr int NumberOfRefinementLayers = 2;
422 double Radius[NumberOfRefinementLayers] = {4, 3};
423 double MaxH[NumberOfRefinementLayers] = {0.15,0.04};
424
426 for (int i=0; i<NumberOfRefinementLayers; i++) {
427 if (radius<Radius[i] and tarch::la::max(volumeH)>MaxH[i]) {
429 }
430 }
431 }
432 }
433 if (CCZ4ReSwi==7){
434 if (tarch::la::equals(t,0.0)){ //as we use a quantity calculated in postpocessing, we need to provide criterion at the first timestep
435 constexpr int NumberOfRefinementLayers = 2;
436 double Radius[NumberOfRefinementLayers] = {3.0, 1.5};
437 double MaxH[NumberOfRefinementLayers] = {0.15,0.04};
439 for (int i=0; i<NumberOfRefinementLayers; i++) {
440 if (radius<Radius[i] and tarch::la::max(volumeH)>MaxH[i]) {
442 }
443 }
444 }
445 }
446 return result;
447}
448
449
450
#define assertion2(expr, param0, param1)
#define assertion4(expr, param0, param1, param2, param3)
#define assertion3(expr, param0, param1, param2)
#define assertion(expr)
#define logError(methodName, logMacroMessageStream)
Wrapper macro around tarch::tarch::logging::Log to improve logging.
Definition Log.h:464
#define logTraceOut(methodName)
Definition Log.h:379
#define logTraceInWith4Arguments(methodName, argument0, argument1, argument2, argument3)
Definition Log.h:373
#define logTraceInWith3Arguments(methodName, argument0, argument1, argument2)
Definition Log.h:372
#define nonCriticalAssertion3(expr, param0, param1, param2)
#define nonCriticalAssertion4(expr, param0, param1, param2, param3)
virtual double maxEigenvalue(const double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &faceCentre, const tarch::la::Vector< Dimensions, double > &volumeH, double t, double dt, int normal) override
Determine max eigenvalue over Jacobian in a given point with solution values (states) Q.
void sourceTerm(const double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &volumeCentre, const tarch::la::Vector< Dimensions, double > &volumeH, double t, double dt, double *__restrict__ S) override
::exahype2::RefinementCommand refinementCriterion(const double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &volumeCentre, const tarch::la::Vector< Dimensions, double > &volumeH, double t) override
Refinement criterion.
virtual void nonconservativeProduct(const double *__restrict__ Q, const double *__restrict__ deltaQ, const tarch::la::Vector< Dimensions, double > &faceCentre, const tarch::la::Vector< Dimensions, double > &volumeH, double t, double dt, int normal, double *__restrict__ BgradQ) override
void initialCondition(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &volumeCentre, const tarch::la::Vector< Dimensions, double > &volumeH, bool gridIsConstructed) override
virtual void boundaryConditions(const double *__restrict__ Qinside, double *__restrict__ Qoutside, const tarch::la::Vector< Dimensions, double > &faceCentre, const tarch::la::Vector< Dimensions, double > &volumeH, double t, int normal) override
Log Device.
Definition Log.h:516
void linearWave(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &X, double t)
void diagonal_gaugeWave(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t)
void flat(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &X, double t)
KeywordToAvoidDuplicateSymbolsForInlinedFunctions double maxEigenvalue(const double *const Q, int normal, const double CCZ4e, const double CCZ4ds, const double CCZ4GLMc, const double CCZ4GLMd) InlineMethod
KeywordToAvoidDuplicateSymbolsForInlinedFunctions void ncp(double *BgradQ, const double *const Q, const double *const gradQSerialised, const int normal, const int CCZ4LapseType, const double CCZ4ds, const double CCZ4c, const double CCZ4e, const double CCZ4f, const double CCZ4bs, const double CCZ4sk, const double CCZ4xi, const double CCZ4mu, const double CCZ4SO) InlineMethod
void gaugeWave(double *__restrict__ Q, const tarch::la::Vector< Dimensions, double > &x, double t)
float M
Definition csv_plot.py:10
double max(double a, double b, double c)
I need the maximum of three values all the time, to I decided to write a function for this.
Definition Scalar.cpp:8
bool equals(const Matrix< Rows, Cols, Scalar > &lhs, const Matrix< Rows, Cols, Scalar > &rhs, const Scalar &tolerance=NUMERICAL_ZERO_DIFFERENCE)
Compares to matrices on equality by means of a numerical accuracy.
Scalar norm2(const Vector< Size, Scalar > &vector)
Computes the 2-norm of the vector, i.e.
double par_S_plus[3]
std::string grid_setup_method
double par_P_plus[3]
double par_S_minus[3]
double center_offset[3]
double target_M_minus
double par_P_minus[3]
Simple vector class.
Definition Vector.h:134