added spherical setup

Author: glesur

test/HD/SedovBlastWave/definitions_spherical.hpp

@@ -0,0 +1,4 @@
+#define     COMPONENTS      3
+#define     DIMENSIONS      3
+
+#define     GEOMETRY        SPHERICAL

test/HD/SedovBlastWave/idefix-spherical.ini

@@ -0,0 +1,28 @@
+[Grid]
+X1-grid    1  0.5  32   u  1.5
+X2-grid    1  0.0  32   u  1.0
+X3-grid    1  0.0  256  u  6.283185307179586
+
+[TimeIntegrator]
+CFL         0.9
+tstop       0.1
+first_dt    1.e-6
+nstages     2
+
+[Hydro]
+solver    hll
+gamma     1.666666666666666666
+
+[Setup]
+Rstart    0.03
+
+[Boundary]
+X1-beg    outflow
+X1-end    outflow
+X2-beg    axis
+X2-end    outflow
+X3-beg    periodic
+X3-end    periodic
+
+[Output]
+vtk    0.01

test/HD/SedovBlastWave/setup.cpp

@@ -20,14 +20,31 @@ void Setup::InitFlow(DataBlock &data) {
     DataBlockHost d(data);
 
     real V = 0;
+    #if GEOMETRY == SPHERICAL
+      const real x0 = 0.1;
+      const real y0 = 0.0;
+      const real z0 = 1.0;
+    #endif
     for(int k = 0; k < d.np_int[KDIR] ; k++) {
         for(int j = 0; j < d.np_int[JDIR] ; j++) {
             for(int i = 0; i < d.np_int[IDIR] ; i++) {
-              real x = d.x[IDIR](i);
-              real y = d.x[JDIR](j);
-              real z = d.x[KDIR](k);
-              real r=sqrt(x*x+y*y+z*z);
-
+              #if GEOMETRY == CARTESIAN
+                real x = d.x[IDIR](i);
+                real y = d.x[JDIR](j);
+                real z = d.x[KDIR](k);
+                real r=sqrt(x*x+y*y+z*z);
+
+              #elif GEOMETRY == SPHERICAL
+                real r = d.x[IDIR](i);
+                real th = d.x[JDIR](j);
+                real phi = d.x[KDIR](k);
+
+                real x = r*sin(th)*cos(phi) - x0;
+                real y = r*sin(th)*sin(phi) - y0;
+                real z = r*cos(th) - z0;
+
+                r=sqrt(x*x+y*y+z*z);
+              #endif
               if(r<Rstart) {
                 V += d.dV(k,j,i);
               }
@@ -41,9 +58,23 @@ void Setup::InitFlow(DataBlock &data) {
     for(int k = 0; k < d.np_tot[KDIR] ; k++) {
         for(int j = 0; j < d.np_tot[JDIR] ; j++) {
             for(int i = 0; i < d.np_tot[IDIR] ; i++) {
-              real x = d.x[IDIR](i);
-              real y = d.x[JDIR](j);
-              real z = d.x[KDIR](k);
+              #if GEOMETRY == CARTESIAN
+                real x = d.x[IDIR](i);
+                real y = d.x[JDIR](j);
+                real z = d.x[KDIR](k);
+                real r=sqrt(x*x+y*y+z*z);
+
+              #elif GEOMETRY == SPHERICAL
+                real r = d.x[IDIR](i);
+                real th = d.x[JDIR](j);
+                real phi = d.x[KDIR](k);
+
+                real x = r*sin(th)*cos(phi) - x0;
+                real y = r*sin(th)*sin(phi) - y0;
+                real z = r*cos(th) - z0;
+
+                r=sqrt(x*x+y*y+z*z);
+              #endif
 
               // Sedov Blast Wave Following Stone+2018, 3.4.4
               d.Vc(RHO,k,j,i) = 1.0;
@@ -53,7 +84,6 @@ void Setup::InitFlow(DataBlock &data) {
 
               d.Vc(PRS,k,j,i) = 0.01;
 
-              real r=sqrt(x*x+y*y+z*z);
               if(r<Rstart) {
                 d.Vc(PRS,k,j,i) = (gamma-1)/V;
               }