3D Navier-Stokes Equations - Rising Thermal Bubble
Location: hypar/Examples/3D/NavierStokes3D/RisingThermalBubble_Config1 (This directory contains all the input files needed to run this case. If there is a Run.m, run it in MATLAB to quickly set up, run, and visualize the example).
Governing equations: 3D Navier-Stokes Equations (navierstokes3d.h)
Domain: \(0 \le x,y,z < 1000\,{\rm m}\), "slip-wall" (_SLIP_WALL_) boundaries everywhere, with zero wall velocity.
Reference:
- Kelly, J. F., Giraldo, F. X., "Continuous and discontinuous Galerkin methods for a scalable three-dimensional nonhydrostatic atmospheric model: Limited-area mode", J. Comput. Phys., 231, 2012, pp. 7988-8008 (see section 5.1.2).
- Giraldo, F. X., Kelly, J. F., Constantinescu, E. M., "Implicit-Explicit Formulations of a Three-Dimensional Nonhydrostatic Unified Model of the Atmosphere (NUMA)", SIAM J. Sci. Comput., 35 (5), 2013, pp. B1162-B1194 (see section 4.1).
Initial solution: A warm bubble in cool ambient atmosphere. Note that in this example, the gravitational forces and rising of the bubble is along the y-axis.
Other parameters (all dimensional quantities are in SI units):
- Specific heat ratio \(\gamma = 1.4\) (NavierStokes3D::gamma)
- Universal gas constant \(R = 287.058\) (NavierStokes3D::R)
- Gravitational force per unit mass \(g = 9.8\) along y-axis (NavierStokes3D::grav_y)
- Reference density (at zero altitude) \(\rho_{ref} = 1.1612055171196529\) (NavierStokes3D::rho0)
- Reference pressure (at zero altitude) \(P_{ref} = 100000\) (NavierStokes3D::p0)
- Hydrostatic balance type 2 (NavierStokes3D::HB)
Numerical Method:
- Spatial discretization (hyperbolic): 5th order WENO (Interp1PrimFifthOrderWENO())
- Time integration: SSP RK3 (TimeRK(), _RK_SSP3_)
Input files required:
solver.inp
begin
ndims 3
nvars 5
size 201 201 201
iproc 4 4 4
ghost 3
n_iter 20000
restart_iter 0
time_scheme rk
time_scheme_type ssprk3
hyp_space_scheme weno5
hyp_flux_split no
hyp_interp_type components
dt 0.01
conservation_check no
screen_op_iter 200
file_op_iter 2000
input_mode serial
ip_file_type binary
output_mode serial
op_file_format binary
op_overwrite no
model navierstokes3d
end
boundary.inp
6
slip-wall 0 1 0 0 0 1000.0 0 1000.0
0.0 0.0 0.0
slip-wall 0 -1 0 0 0 1000.0 0 1000.0
0.0 0.0 0.0
slip-wall 1 1 0 1000.0 0 0 0 1000.0
0.0 0.0 0.0
slip-wall 1 -1 0 1000.0 0 0 0 1000.0
0.0 0.0 0.0
slip-wall 2 1 0 1000.0 0 1000.0 0 0
0.0 0.0 0.0
slip-wall 2 -1 0 1000.0 0 1000.0 0 0
0.0 0.0 0.0
physics.inp
begin
gamma 1.4
upwinding rusanov
gravity 0.0 9.8 0.0
rho_ref 1.1612055171196529
p_ref 100000.0
R 287.058
HB 2
end
weno.inp (optional)
begin
mapped 0
borges 0
yc 1
no_limiting 0
epsilon 0.000001
p 2.0
rc 0.3
xi 0.001
end
To generate initial.inp (initial solution), compile and run the following code in the run directory.
/*
This code generates the initial solution for the rising
thermal bubble case for the 3D Navier-Stokes equations.
Note: this code allocates the global domain, so be careful
when using it for a large number of grid points.
*/
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <string.h>
{
return(exp(a*log(x)));
}
int main()
{
double gamma = 1.4;
double R = 287.058;
double rho_ref = 1.1612055171196529;
double p_ref = 100000.0;
double grav_x = 0.0;
double grav_y = 9.8;
double grav_z = 0.0;
int HB = 0;
int NI,NJ,NK,ndims;
char ip_file_type[50]; strcpy(ip_file_type,"ascii");
FILE *in;
printf("Reading file \"solver.inp\"...\n");
in = fopen("solver.inp","r");
if (!in) {
printf("Error: Input file \"solver.inp\" not found.\n");
return(0);
} else {
char word[500];
fscanf(in,"%s",word);
if (!strcmp(word, "begin")) {
while (strcmp(word, "end")) {
fscanf(in,"%s",word);
if (!strcmp(word, "ndims")) fscanf(in,"%d",&ndims);
else if (!strcmp(word, "size")) {
fscanf(in,"%d",&NI);
fscanf(in,"%d",&NJ);
fscanf(in,"%d",&NK);
} else if (!strcmp(word, "ip_file_type")) fscanf(in,"%s",ip_file_type);
}
} else printf("Error: Illegal format in solver.inp. Crash and burn!\n");
}
fclose(in);
printf("Reading file \"physics.inp\"...\n");
in = fopen("physics.inp","r");
if (!in) {
printf("Error: Input file \"physics.inp\" not found.\n");
return(0);
} else {
char word[500];
fscanf(in,"%s",word);
if (!strcmp(word, "begin")) {
while (strcmp(word, "end")) {
fscanf(in,"%s",word);
if (!strcmp(word, "rho_ref")) fscanf(in,"%lf",&rho_ref);
else if (!strcmp(word, "p_ref" )) fscanf(in,"%lf",&p_ref );
else if (!strcmp(word, "gamma" )) fscanf(in,"%lf",&gamma );
else if (!strcmp(word, "R" )) fscanf(in,"%lf",&R );
else if (!strcmp(word, "HB" )) fscanf(in,"%d" ,&HB );
else if (!strcmp(word, "gravity")) {
fscanf(in,"%lf",&grav_x );
fscanf(in,"%lf",&grav_y );
fscanf(in,"%lf",&grav_z );
}
}
} else printf("Error: Illegal format in physics.inp. Crash and burn!\n");
}
fclose(in);
if (ndims != 3) {
printf("ndims is not 3 in solver.inp. this code is to generate 3D initial conditions\n");
return(0);
}
if (HB != 2) {
printf("Error: Specify \"HB\" as 2 in physics.inp.\n");
}
if ((grav_x != 0.0) || (grav_z != 0.0)) {
printf("Error: gravity must be zero along x and z for this example.\n");
return(0);
}
printf("Grid:\t\t\t%d X %d X %d\n",NI,NJ,NK);
printf("Reference density and pressure: %lf, %lf.\n",rho_ref,p_ref);
int i,j,k;
double dx = 1000.0 / ((double)(NI-1));
double dy = 1000.0 / ((double)(NJ-1));
double dz = 1000.0 / ((double)(NK-1));
double *x, *y, *z, *U;
x = (double*) calloc (NI , sizeof(double));
y = (double*) calloc (NJ , sizeof(double));
z = (double*) calloc (NK , sizeof(double));
U = (double*) calloc (5*NI*NJ*NK, sizeof(double));
/* Initial perturbation center */
double xc = 500;
double yc = 260;
double zc = 500;
double Cp = gamma * R / (gamma-1.0);
/* initial perturbation parameters */
double pi = 4.0*atan(1.0);
double rc = 250.0;
double T_ref = p_ref / (R * rho_ref);
for (i = 0; i < NI; i++){
for (j = 0; j < NJ; j++){
for (k = 0; k < NK; k++){
x[i] = i*dx;
y[j] = j*dy;
z[k] = k*dz;
int p = i + NI*j + NI*NJ*k;
/* temperature peturbation */
double r = sqrt((x[i]-xc)*(x[i]-xc)+(y[j]-yc)*(y[j]-yc)+(z[k]-zc)*(z[k]-zc));
double dtheta = (r>rc ? 0.0 : (0.5*(1.0+cos(pi*r/rc))) );
double theta = T_ref + dtheta;
double Pexner = 1.0 - (grav_y*y[j])/(Cp*T_ref);
double rho = (p_ref/(R*theta)) * raiseto(Pexner, (1.0/(gamma-1.0)));
double E = rho * (R/(gamma-1.0)) * theta*Pexner;
U[5*p+0] = rho;
U[5*p+1] = 0.0;
U[5*p+2] = 0.0;
U[5*p+3] = 0.0;
U[5*p+4] = E;
}
}
}
FILE *out;
if (!strcmp(ip_file_type,"ascii")) {
printf("Error: sorry, ASCII format initial solution file not implemented. ");
printf("Please choose \"ip_file_format\" in solver.inp as \"binary\".\n");
return(0);
} else if ((!strcmp(ip_file_type,"binary")) || (!strcmp(ip_file_type,"bin"))) {
printf("Writing binary initial solution file initial.inp\n");
out = fopen("initial.inp","wb");
fwrite(x,sizeof(double),NI,out);
fwrite(y,sizeof(double),NJ,out);
fwrite(z,sizeof(double),NK,out);
fwrite(U,sizeof(double),5*NI*NJ*NK,out);
fclose(out);
}
free(x);
free(y);
free(z);
free(U);
return(0);
}
Output:
Note that iproc is set to
4 4 4
in solver.inp (i.e., 4 processors along x, 4 processors along y, and 4 processor along z). Thus, this example should be run with 64 MPI ranks (or change iproc).
After running the code, there should be 11 output files op_00000.bin, op_00001.bin, ... op_00010.bin; the first one is the solution at \(t=0\) and the final one is the solution at \(t=200\,{\rm s}\). Since HyPar::op_overwrite is set to no in solver.inp, separate files are written for solutions at each output time. All the files are binary (HyPar::op_file_format is set to binary in solver.inp).
The following code (Examples/3D/NavierStokes3D/RisingThermalBubble_Config1/aux/PostProcess.c) can be used to convert the binary solution file (with conserved variables \(\rho,\rho u,\rho v,\rho w,e\)) to Tecplot or plain text files with the primitive and reference variables \(\rho,u,v,w,P,\theta,\rho_0,P_0,\pi,\theta_0\) where the subscript \(0\) indicates the hydrostatic mean value.
/*
Rising Thermal Bubble:-
The code takes a binary solution file (that contains the
conserved variable (rho, rho*u, rho*v, rho*w, e) as its input
and calculates the primitive atmospheric flow variables:
rho, u, v, w, P, theta, pi, rho0, P0, theta0, pi0
and writes them to a tecplot or text file.
*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
typedef struct _parameters_{
double grav_x, grav_y, grav_z, R, gamma, P_ref, rho_ref;
int HB;
} Parameters;
#define _ArraySetValue_(x,size,value) \
{ \
int arraycounter; \
for (arraycounter = 0; arraycounter < (size); arraycounter++) x[arraycounter] = (value); \
}
#define _ArrayIncrementIndex_(N,imax,i,done) \
{ \
int arraycounter = 0; \
while (arraycounter < (N)) { \
if (i[arraycounter] == imax[arraycounter]-1) { \
i[arraycounter] = 0; \
arraycounter++; \
} else { \
i[arraycounter]++; \
break; \
} \
} \
if (arraycounter == (N)) done = 1; \
else done = 0; \
}
#define _ArrayIndex1D_(N,imax,i,ghost,index) \
{ \
index = i[N-1]+(ghost); \
int arraycounter; \
for (arraycounter = (N)-2; arraycounter > -1; arraycounter--) { \
index = ((index*(imax[arraycounter]+2*(ghost))) + (i[arraycounter]+(ghost))); \
} \
}
void IncrementFilename(char *f)
{
if (f[7] == '9') {
f[7] = '0';
if (f[6] == '9') {
f[6] = '0';
if (f[5] == '9') {
f[5] = '0';
if (f[4] == '9') {
f[4] = '0';
if (f[3] == '9') {
f[3] = '0';
fprintf(stderr,"Warning: file increment hit max limit. Resetting to zero.\n");
} else {
f[3]++;
}
} else {
f[4]++;
}
} else {
f[5]++;
}
} else {
f[6]++;
}
} else {
f[7]++;
}
}
{
return(exp(a*log(x)));
}
{
if (ndims !=3) {
fprintf(stderr,"Error in WriteTecplot3D(): This functions is hardcoded for 3-dimensional ");
fprintf(stderr,"problems only. Instead, ndims=%d.\n",ndims);
return(1);
}
int i;
int imax = dim[0];
int jmax = dim[1];
int kmax = dim[2];
printf("Writing tecplot solution file %s.\n",f);
FILE *out;
out = fopen(f,"w");
if (!out) {
fprintf(stderr,"Error: could not open %s for writing.\n",f);
return(1);
}
/* writing tecplot data file headers */
fprintf(out,"VARIABLES=\"I\",\"J\",\"K\",\"X\",\"Y\",\"Z\",");
char varname[3] = "00";
for (i = 0; i < nvars; i++) {
fprintf(out,"\"%s\",",varname);
if (varname[1] == '9') { varname[0]++; varname[1] = '0'; }
else varname[1]++;
}
fprintf(out,"\n");
fprintf(out,"ZONE I=%d,J=%d,K=%d,F=POINT\n",imax,jmax,kmax);
/* writing the data */
int done = 0; _ArraySetValue_(index,ndims,0);
while (!done) {
int i, p;
_ArrayIndex1D_(ndims,dim,index,0,p);
for (i=0; i<ndims; i++) fprintf(out,"%4d ",index[i]);
for (i=0; i<ndims; i++) {
int j,offset = 0; for (j=0; j<i; j++) offset += dim[j];
fprintf(out,"%+E ",x[offset+index[i]]);
}
for (i=0; i<nvars; i++) fprintf(out,"%+E ",u[nvars*p+i]);
fprintf(out,"\n");
_ArrayIncrementIndex_(ndims,dim,index,done);
}
fclose(out);
return(0);
}
{
FILE *out;
out = fopen(f,"w");
if (!out) {
fprintf(stderr,"Error: could not open %s for writing.\n",f);
return(1);
}
int done = 0; _ArraySetValue_(index,ndims,0);
while (!done) {
int i, p;
_ArrayIndex1D_(ndims,dim,index,0,p);
for (i=0; i<ndims; i++) fprintf(out,"%4d ",index[i]);
for (i=0; i<ndims; i++) {
int j,offset = 0; for (j=0; j<i; j++) offset += dim[j];
fprintf(out,"%+E ",x[offset+index[i]]);
}
for (i=0; i<nvars; i++) fprintf(out,"%+E ",u[nvars*p+i]);
fprintf(out,"\n");
_ArrayIncrementIndex_(ndims,dim,index,done);
}
fclose(out);
return(0);
}
int PostProcess(char *fname, char *oname, void *p, int flag)
{
Parameters *params = (Parameters*) p;
FILE *in; in = fopen(fname,"rb");
if (!in) return(-1);
printf("Reading file %s.\n",fname);
int ndims, nvars;
double *U,*x;
/* read the file headers */
fread(&ndims,sizeof(int),1,in);
fread(&nvars,sizeof(int),1,in);
/* some checks */
if (ndims != 3) {
printf("Error: ndims in %s not equal to 3!\n",fname);
return(1);
}
if (nvars != 5) {
printf("Error: nvars in %s not equal to 5!\n",fname);
return(1);
}
/* read dimensions */
int dims[ndims];
fread(dims,sizeof(int),ndims,in);
printf("Dimensions: %d x %d x %d\n",dims[0],dims[1],dims[2]);
printf("Nvars : %d\n",nvars);
/* allocate grid and solution arrays */
x = (double*) calloc (dims[0]+dims[1]+dims[2] ,sizeof(double));
U = (double*) calloc (dims[0]*dims[1]*dims[2]*nvars ,sizeof(double));
/* read grid and solution */
fread(x,sizeof(double),dims[0]+dims[1]+dims[2] ,in);
fread(U,sizeof(double),dims[0]*dims[1]*dims[2]*nvars,in);
/* done reading */
fclose(in);
int imax = dims[0];
int jmax = dims[1];
int kmax = dims[2];
/* allocate primitive variable array (rho, u, v, w, P, theta, rho0, P0, pi0, theta0) */
int evars = 5;
double *Q = (double*) calloc ((nvars+evars)*imax*jmax*kmax,sizeof(double));
/* calculate primitive variables */
int i, j, k;
double *X = x;
double *Y = x+imax;
double *Z = x+imax+jmax;
double grav_y = params->grav_y;
double R = params->R;
double gamma = params->gamma;
double P_ref = params->P_ref;
double rho_ref = params->rho_ref;
double T_ref = P_ref / (R*rho_ref);
double inv_gamma_m1 = 1.0 / (gamma-1.0);
double Cp = gamma * inv_gamma_m1 * R;
for (i=0; i<imax; i++) {
for (j=0; j<jmax; j++) {
for (k=0; k<kmax; k++) {
int p = i + imax*j + imax*jmax*k;
double rho0, theta0, Pexner, P0;
theta0 = T_ref;
Pexner = 1.0 - (grav_y*Y[j])/(Cp*T_ref);
rho0 = (P_ref/(R*theta0)) * raiseto(Pexner, inv_gamma_m1);
P0 = P_ref * raiseto(Pexner, gamma*inv_gamma_m1);
double rho, uvel, vvel, wvel, E, P, theta;
rho = U[nvars*p+0];
uvel = U[nvars*p+1] / rho;
vvel = U[nvars*p+2] / rho;
wvel = U[nvars*p+3] / rho;
E = U[nvars*p+4];
P = (gamma-1.0) * (E - 0.5*rho*(uvel*uvel+vvel*vvel+wvel*wvel));
theta = (E-0.5*rho*(uvel*uvel+vvel*vvel+wvel*wvel))/(Pexner*rho) * ((gamma-1.0)/R);
Q[(nvars+evars)*p+0] = rho;
Q[(nvars+evars)*p+1] = uvel;
Q[(nvars+evars)*p+2] = vvel;
Q[(nvars+evars)*p+3] = wvel;
Q[(nvars+evars)*p+4] = P;
Q[(nvars+evars)*p+5] = theta;
Q[(nvars+evars)*p+6] = rho0;
Q[(nvars+evars)*p+7] = P0;
Q[(nvars+evars)*p+8] = Pexner;
Q[(nvars+evars)*p+9] = theta0;
}
}
}
int dim[ndims], index[ndims];
dim[0] = imax;
dim[1] = jmax;
dim[2] = kmax;
/* write Tecplot/Text file */
if (flag) WriteTecplot3D(ndims,nvars+evars,&dim[0],x,Q,oname,&index[0]);
else WriteText (ndims,nvars+evars,&dim[0],x,Q,oname,&index[0]);
/* clean up */
free(U);
free(Q);
free(x);
}
int main()
{
FILE *out1, *out2, *in, *inputs;
char filename[50], op_file_format[50], tecfile[50], overwrite[50];
int flag;
printf("Write tecplot file (1) or plain text file (0): ");
scanf("%d",&flag);
if ((flag != 1) && (flag != 0)) {
printf("Error: Invalid input. Should be 1 or 0.\n");
return(0);
}
printf("Reading solver.inp.\n");
inputs = fopen("solver.inp","r");
if (!inputs) {
fprintf(stderr,"Error: File \"solver.inp\" not found.\n");
return(1);
} else {
char word[100];
fscanf(inputs,"%s",word);
if (!strcmp(word, "begin")){
while (strcmp(word, "end")){
fscanf(inputs,"%s",word);
if (!strcmp(word, "op_file_format" )) fscanf(inputs,"%s" ,op_file_format);
else if (!strcmp(word, "op_overwrite" )) fscanf(inputs,"%s" ,overwrite );
}
}
fclose(inputs);
}
if (strcmp(op_file_format,"binary") && strcmp(op_file_format,"bin")) {
printf("Error: solution output needs to be in binary files.\n");
return(0);
}
Parameters params;
/* default values */
params.grav_x = 0.0;
params.grav_y = 9.8;
params.grav_z = 0.0;
params.R = 287.058;
params.gamma = 1.4;
params.P_ref = 100000.0;
params.rho_ref = 100000.0 / (params.R * 300.0);
params.HB = 0;
/* read these parameters from file */
printf("Reading physics.inp.\n");
inputs = fopen("physics.inp","r");
if (!inputs) {
fprintf(stderr,"Error: File \"physics.inp\" not found.\n");
return(1);
} else {
char word[100];
fscanf(inputs,"%s",word);
if (!strcmp(word, "begin")){
while (strcmp(word, "end")){
fscanf(inputs,"%s",word);
if (!strcmp(word, "gamma")) fscanf(inputs,"%lf",¶ms.gamma);
else if (!strcmp(word, "gravity")) {
fscanf(inputs,"%lf",¶ms.grav_x);
fscanf(inputs,"%lf",¶ms.grav_y);
fscanf(inputs,"%lf",¶ms.grav_z);
} else if (!strcmp(word,"p_ref")) fscanf(inputs,"%lf",¶ms.P_ref);
else if (!strcmp(word,"rho_ref")) fscanf(inputs,"%lf",¶ms.rho_ref);
else if (!strcmp(word,"HB")) fscanf(inputs,"%d",¶ms.HB);
}
}
fclose(inputs);
}
if (!strcmp(overwrite,"no")) {
strcpy(filename,"op_00000.bin");
while(1) {
/* set filename */
strcpy(tecfile,filename);
tecfile[9] = 'd';
tecfile[10] = 'a';
tecfile[11] = 't';
int err = PostProcess(filename, tecfile, ¶ms, flag);
if (err == -1) {
printf("No more files found. Exiting.\n");
break;
}
IncrementFilename(filename);
}
} else if (!strcmp(overwrite,"yes")) {
strcpy(filename,"op.bin");
/* set filename */
strcpy(tecfile,filename);
tecfile[3] = 'd';
tecfile[4] = 'a';
tecfile[5] = 't';
int err = PostProcess(filename, tecfile, ¶ms, flag);
if (err == -1) {
printf("Error: op.bin not found.\n");
return(0);
}
}
return(0);
}
The following figure shows the potential temperature iso-surface for the initial and final solutions (plotted in VisIt):
The file Extras/ExtractSlice.c can be used to extract a slice perpendicular to any dimension at a specified location along that dimension. The extract slice is written out in the same binary format as the original solutions files (with the same names op_xxxxx.bin) in a subdirectory called slices (Note: make the subdirectory called slices before running this code). The following code (Examples/3D/NavierStokes3D/RisingThermalBubble_Config1/aux/PostProcessSlice.c) can then be used (in the slices subdirectory) to convert the binary slice solution file (with conserved variables \(\rho,\rho u,\rho v,\rho w,e\)) to Tecplot or plain text files with the primitive and reference variables \(\rho,u,v,w,P,\theta,\rho_0,P_0,\pi,\theta_0\). Note that it needs the relevant solver.inp and physics.inp that can be created as follows:
- Copy the original solver.inp and physics.inp to the slices subdirectory.
- In solver.inp, set ndims as 2, and remove the component of size and iproc corresponding to the dimension being eliminated while extracting the slices (in this case, it is z or the 3rd component).
- In physics.inp, remove the component of gravity corresponding to the dimension perpendicular to the slice.
/*
Rising Thermal Bubble:-
The code takes a binary solution file (that contains the
conserved variable (rho, rho*u, rho*v, e) as its input
and calculates the primitive atmospheric flow variables:
rho, u, v, P, theta, pi, rho0, P0, theta0, pi0
and writes them to a tecplot file.
*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
typedef struct _parameters_{
double grav_x, grav_y, R, gamma, P_ref, rho_ref;
int HB;
} Parameters;
void IncrementFilename(char *f)
{
if (f[7] == '9') {
f[7] = '0';
if (f[6] == '9') {
f[6] = '0';
if (f[5] == '9') {
f[5] = '0';
if (f[4] == '9') {
f[4] = '0';
if (f[3] == '9') {
f[3] = '0';
fprintf(stderr,"Warning: file increment hit max limit. Resetting to zero.\n");
} else {
f[3]++;
}
} else {
f[4]++;
}
} else {
f[5]++;
}
} else {
f[6]++;
}
} else {
f[7]++;
}
}
{
return(exp(a*log(x)));
}
{
printf("\tWriting tecplot solution file %s.\n",f);
FILE *out;
out = fopen(f,"w");
if (!out) {
fprintf(stderr,"Error: could not open %s for writing.\n",f);
return;
}
double *X = x;
double *Y = x+imax;
/* writing tecplot data file headers */
fprintf(out,"VARIABLES=\"I\",\"J\",\"X\",\"Y\",");
fprintf(out,"\"RHO\",\"U\",\"V\",\"W\",\"P\",");
fprintf(out,"\"THETA\",\"RHO0\",\"P0\",");
fprintf(out,"\"PI0\",\"THETA0\",\n");
fprintf(out,"ZONE I=%d,J=%d,F=POINT\n",imax,jmax);
/* writing the data */
int i,j;
for (j=0; j<jmax; j++) {
for (i=0; i<imax; i++) {
int v, p = i + imax*j;
fprintf(out,"%4d %4d ",i,j);
fprintf(out,"%1.16E %1.16E ",X[i],Y[j]);
for (v=0; v<nvars; v++) fprintf(out,"%1.16E ",u[nvars*p+v]);
fprintf(out,"\n");
}
}
fclose(out);
return;
}
void WriteText2D(int nvars,int imax, int jmax,double *x,double *u,char *f)
{
printf("\tWriting text solution file %s.\n",f);
FILE *out;
out = fopen(f,"w");
if (!out) {
fprintf(stderr,"Error: could not open %s for writing.\n",f);
return;
}
double *X = x;
double *Y = x+imax;
/* writing the data */
int i,j;
for (j=0; j<jmax; j++) {
for (i=0; i<imax; i++) {
int v, p = i + imax*j;
fprintf(out,"%4d %4d ",i,j);
fprintf(out,"%1.16E %1.16E ",X[i],Y[j]);
for (v=0; v<nvars; v++) fprintf(out,"%1.16E ",u[nvars*p+v]);
fprintf(out,"\n");
}
}
fclose(out);
return;
}
int PostProcess(char *fname, char *oname, void *p, int flag)
{
Parameters *params = (Parameters*) p;
FILE *in; in = fopen(fname,"rb");
if (!in) return(-1);
printf("Reading file %s.\n",fname);
int ndims, nvars;
double *U,*x;
/* read the file headers */
fread(&ndims,sizeof(int),1,in);
fread(&nvars,sizeof(int),1,in);
/* some checks */
if (ndims != 2) {
printf("Error: ndims in %s not equal to 2!\n",fname);
return(1);
}
if (nvars != 5) {
printf("Error: nvars in %s not equal to 5!\n",fname);
return(1);
}
/* read dimensions */
int dims[ndims];
fread(dims,sizeof(int),ndims,in);
printf("Dimensions: %d x %d\n",dims[0],dims[1]);
printf("Nvars : %d\n",nvars);
/* allocate grid and solution arrays */
x = (double*) calloc (dims[0]+dims[1] ,sizeof(double));
U = (double*) calloc (dims[0]*dims[1]*nvars ,sizeof(double));
/* read grid and solution */
fread(x,sizeof(double),dims[0]+dims[1] ,in);
fread(U,sizeof(double),dims[0]*dims[1]*nvars,in);
/* done reading */
fclose(in);
int imax = dims[0];
int jmax = dims[1];
/* allocate primitive variable array (rho, u, v, w, P, theta, rho0, P0, pi0, theta0) */
int evars = 5;
double *Q = (double*) calloc ((nvars+evars)*imax*jmax,sizeof(double));
/* calculate primitive variables */
int i, j;
double *X = x;
double *Y = x+imax;
double grav_x = params->grav_x;
double grav_y = params->grav_y;
double R = params->R;
double gamma = params->gamma;
double P_ref = params->P_ref;
double rho_ref = params->rho_ref;
double T_ref = P_ref / (R*rho_ref);
double inv_gamma_m1 = 1.0 / (gamma-1.0);
double Cp = gamma * inv_gamma_m1 * R;
for (i=0; i<imax; i++) {
for (j=0; j<jmax; j++) {
int p = i + imax*j;
double rho0, theta0, Pexner, P0;
theta0 = T_ref;
Pexner = 1.0 - (grav_x*X[i]+grav_y*Y[j])/(Cp*T_ref);
rho0 = (P_ref/(R*theta0)) * raiseto(Pexner, inv_gamma_m1);
P0 = P_ref * raiseto(Pexner, gamma*inv_gamma_m1);
double rho, uvel, vvel, wvel, E, P, theta;
rho = U[nvars*p+0];
uvel = U[nvars*p+1] / rho;
vvel = U[nvars*p+2] / rho;
wvel = U[nvars*p+3] / rho;
E = U[nvars*p+4];
P = (gamma-1.0) * (E - 0.5*rho*(uvel*uvel+vvel*vvel+wvel*wvel));
theta = (E-0.5*rho*(uvel*uvel+vvel*vvel+wvel*wvel))/(Pexner*rho) * ((gamma-1.0)/R);
Q[(nvars+evars)*p+0] = rho;
Q[(nvars+evars)*p+1] = uvel;
Q[(nvars+evars)*p+2] = vvel;
Q[(nvars+evars)*p+3] = wvel;
Q[(nvars+evars)*p+4] = P;
Q[(nvars+evars)*p+5] = theta;
Q[(nvars+evars)*p+6] = rho0;
Q[(nvars+evars)*p+7] = P0;
Q[(nvars+evars)*p+8] = Pexner;
Q[(nvars+evars)*p+9] = theta0;
}
}
/* write Tecplot/Text file */
if (flag) WriteTecplot2D(nvars+evars,imax,jmax,x,Q,oname);
else WriteText2D (nvars+evars,imax,jmax,x,Q,oname);
/* clean up */
free(U);
free(Q);
free(x);
}
int main()
{
FILE *out1, *out2, *in, *inputs;
char filename[50], op_file_format[50], tecfile[50], overwrite[50];
int flag;
printf("Write tecplot file (1) or plain text file (0): ");
scanf("%d",&flag);
if ((flag != 1) && (flag != 0)) {
printf("Error: Invalid input. Should be 1 or 0.\n");
return(0);
}
printf("Reading solver.inp.\n");
inputs = fopen("solver.inp","r");
if (!inputs) {
fprintf(stderr,"Error: File \"solver.inp\" not found.\n");
return(1);
} else {
char word[100];
fscanf(inputs,"%s",word);
if (!strcmp(word, "begin")){
while (strcmp(word, "end")){
fscanf(inputs,"%s",word);
if (!strcmp(word, "op_file_format" )) fscanf(inputs,"%s" ,op_file_format);
else if (!strcmp(word, "op_overwrite" )) fscanf(inputs,"%s" ,overwrite );
}
}
fclose(inputs);
}
if (strcmp(op_file_format,"binary") && strcmp(op_file_format,"bin")) {
printf("Error: solution output needs to be in binary files.\n");
return(0);
}
Parameters params;
/* default values */
params.grav_x = 0.0;
params.grav_y = 9.8;
params.R = 287.058;
params.gamma = 1.4;
params.P_ref = 100000.0;
params.rho_ref = 100000.0 / (params.R * 300.0);
params.HB = 0;
/* read these parameters from file */
printf("Reading physics.inp.\n");
inputs = fopen("physics.inp","r");
if (!inputs) {
fprintf(stderr,"Error: File \"physics.inp\" not found.\n");
return(1);
} else {
char word[100];
fscanf(inputs,"%s",word);
if (!strcmp(word, "begin")){
while (strcmp(word, "end")){
fscanf(inputs,"%s",word);
if (!strcmp(word, "gamma")) fscanf(inputs,"%lf",¶ms.gamma);
else if (!strcmp(word, "gravity")) {
fscanf(inputs,"%lf",¶ms.grav_x);
fscanf(inputs,"%lf",¶ms.grav_y);
} else if (!strcmp(word,"p_ref")) fscanf(inputs,"%lf",¶ms.P_ref);
else if (!strcmp(word,"rho_ref")) fscanf(inputs,"%lf",¶ms.rho_ref);
else if (!strcmp(word,"HB")) fscanf(inputs,"%d",¶ms.HB);
}
}
fclose(inputs);
}
if (params.HB != 2) {
printf("Error: \"HB\" must be specified as 2 in physics.inp.\n");
return(0);
}
if (!strcmp(overwrite,"no")) {
strcpy(filename,"op_00000.bin");
while(1) {
/* set filename */
strcpy(tecfile,filename);
tecfile[9] = 'd';
tecfile[10] = 'a';
tecfile[11] = 't';
int err = PostProcess(filename, tecfile, ¶ms, flag);
if (err == -1) {
printf("No more files found. Exiting.\n");
break;
}
IncrementFilename(filename);
}
} else if (!strcmp(overwrite,"yes")) {
strcpy(filename,"op.bin");
/* set filename */
strcpy(tecfile,filename);
tecfile[3] = 'd';
tecfile[4] = 'a';
tecfile[5] = 't';
int err = PostProcess(filename, tecfile, ¶ms, flag);
if (err == -1) {
printf("Error: op.bin not found.\n");
return(0);
}
}
return(0);
}
The following figure shows the potential temperature \(\theta\) along a slice at \(z=500\,{\rm m}\) (plotted in VisIt):
Expected screen output:
HyPar - Parallel (MPI) version with 64 processes
Reading solver inputs from file "solver.inp".
No. of dimensions : 3
No. of variables : 5
Domain size : 201 201 201
Processes along each dimension : 4 4 4
No. of ghosts pts : 3
No. of iter. : 20000
Restart iteration : 0
Time integration scheme : rk (ssprk3)
Spatial discretization scheme (hyperbolic) : weno5
Split hyperbolic flux term? : no
Interpolation type for hyperbolic term : components
Spatial discretization type (parabolic ) : nonconservative-1stage
Spatial discretization scheme (parabolic ) : 2
Time Step : 1.000000E-02
Check for conservation : no
Screen output iterations : 200
File output iterations : 2000
Initial solution file type : binary
Initial solution read mode : serial
Solution file write mode : serial
Solution file format : binary
Overwrite solution file : no
Physical model : navierstokes3d
Partitioning domain.
Allocating data arrays.
Reading array from binary file initial.inp (Serial mode).
Volume integral of the initial solution:
0: 1.1315344432084937E+09
1: 0.0000000000000000E+00
2: 0.0000000000000000E+00
3: 0.0000000000000000E+00
4: 2.3971700067997931E+14
Reading boundary conditions from "boundary.inp".
Boundary slip-wall: Along dimension 0 and face +1
Boundary slip-wall: Along dimension 0 and face -1
Boundary slip-wall: Along dimension 1 and face +1
Boundary slip-wall: Along dimension 1 and face -1
Boundary slip-wall: Along dimension 2 and face +1
Boundary slip-wall: Along dimension 2 and face -1
6 boundary condition(s) read.
Initializing solvers.
Reading WENO parameters from weno.inp.
Initializing physics. Model = "navierstokes3d"
Reading physical model inputs from file "physics.inp".
Setting up time integration.
Solving in time (from 0 to 20000 iterations)
Writing solution file op_00000.bin.
Iteration: 200 Time: 2.000E+00 Max CFL: 6.945E-01 Max Diff. No.: -1.000E+00 Norm: 1.0902E-02
Iteration: 400 Time: 4.000E+00 Max CFL: 6.945E-01 Max Diff. No.: -1.000E+00 Norm: 9.2576E-03
Iteration: 600 Time: 6.000E+00 Max CFL: 6.945E-01 Max Diff. No.: -1.000E+00 Norm: 8.7706E-03
Iteration: 800 Time: 8.000E+00 Max CFL: 6.945E-01 Max Diff. No.: -1.000E+00 Norm: 1.1474E-02
Iteration: 1000 Time: 1.000E+01 Max CFL: 6.945E-01 Max Diff. No.: -1.000E+00 Norm: 1.0926E-02
Iteration: 1200 Time: 1.200E+01 Max CFL: 6.945E-01 Max Diff. No.: -1.000E+00 Norm: 1.0338E-02
Iteration: 1400 Time: 1.400E+01 Max CFL: 6.945E-01 Max Diff. No.: -1.000E+00 Norm: 1.2228E-02
Iteration: 1600 Time: 1.600E+01 Max CFL: 6.945E-01 Max Diff. No.: -1.000E+00 Norm: 1.1628E-02
Iteration: 1800 Time: 1.800E+01 Max CFL: 6.946E-01 Max Diff. No.: -1.000E+00 Norm: 1.1837E-02
Iteration: 2000 Time: 2.000E+01 Max CFL: 6.946E-01 Max Diff. No.: -1.000E+00 Norm: 1.1553E-02
Writing solution file op_00001.bin.
Iteration: 2200 Time: 2.200E+01 Max CFL: 6.946E-01 Max Diff. No.: -1.000E+00 Norm: 1.0333E-02
Iteration: 2400 Time: 2.400E+01 Max CFL: 6.946E-01 Max Diff. No.: -1.000E+00 Norm: 1.1710E-02
Iteration: 2600 Time: 2.600E+01 Max CFL: 6.946E-01 Max Diff. No.: -1.000E+00 Norm: 8.4824E-03
Iteration: 2800 Time: 2.800E+01 Max CFL: 6.946E-01 Max Diff. No.: -1.000E+00 Norm: 1.0401E-02
Iteration: 3000 Time: 3.000E+01 Max CFL: 6.946E-01 Max Diff. No.: -1.000E+00 Norm: 9.8131E-03
Iteration: 3200 Time: 3.200E+01 Max CFL: 6.946E-01 Max Diff. No.: -1.000E+00 Norm: 7.4465E-03
Iteration: 3400 Time: 3.400E+01 Max CFL: 6.946E-01 Max Diff. No.: -1.000E+00 Norm: 8.5672E-03
Iteration: 3600 Time: 3.600E+01 Max CFL: 6.947E-01 Max Diff. No.: -1.000E+00 Norm: 1.0509E-02
Iteration: 3800 Time: 3.800E+01 Max CFL: 6.947E-01 Max Diff. No.: -1.000E+00 Norm: 7.2593E-03
Iteration: 4000 Time: 4.000E+01 Max CFL: 6.947E-01 Max Diff. No.: -1.000E+00 Norm: 1.1689E-02
Writing solution file op_00002.bin.
Iteration: 4200 Time: 4.200E+01 Max CFL: 6.947E-01 Max Diff. No.: -1.000E+00 Norm: 1.0779E-02
Iteration: 4400 Time: 4.400E+01 Max CFL: 6.947E-01 Max Diff. No.: -1.000E+00 Norm: 1.0578E-02
Iteration: 4600 Time: 4.600E+01 Max CFL: 6.948E-01 Max Diff. No.: -1.000E+00 Norm: 1.1789E-02
Iteration: 4800 Time: 4.800E+01 Max CFL: 6.948E-01 Max Diff. No.: -1.000E+00 Norm: 1.1421E-02
Iteration: 5000 Time: 5.000E+01 Max CFL: 6.949E-01 Max Diff. No.: -1.000E+00 Norm: 1.1930E-02
Iteration: 5200 Time: 5.200E+01 Max CFL: 6.949E-01 Max Diff. No.: -1.000E+00 Norm: 1.1638E-02
Iteration: 5400 Time: 5.400E+01 Max CFL: 6.950E-01 Max Diff. No.: -1.000E+00 Norm: 1.1239E-02
Iteration: 5600 Time: 5.600E+01 Max CFL: 6.950E-01 Max Diff. No.: -1.000E+00 Norm: 1.0715E-02
Iteration: 5800 Time: 5.800E+01 Max CFL: 6.951E-01 Max Diff. No.: -1.000E+00 Norm: 9.7812E-03
Iteration: 6000 Time: 6.000E+01 Max CFL: 6.952E-01 Max Diff. No.: -1.000E+00 Norm: 8.9185E-03
Writing solution file op_00003.bin.
Iteration: 6200 Time: 6.200E+01 Max CFL: 6.952E-01 Max Diff. No.: -1.000E+00 Norm: 1.1160E-02
Iteration: 6400 Time: 6.400E+01 Max CFL: 6.953E-01 Max Diff. No.: -1.000E+00 Norm: 7.4901E-03
Iteration: 6600 Time: 6.600E+01 Max CFL: 6.953E-01 Max Diff. No.: -1.000E+00 Norm: 9.5477E-03
Iteration: 6800 Time: 6.800E+01 Max CFL: 6.954E-01 Max Diff. No.: -1.000E+00 Norm: 9.6800E-03
Iteration: 7000 Time: 7.000E+01 Max CFL: 6.954E-01 Max Diff. No.: -1.000E+00 Norm: 8.0774E-03
Iteration: 7200 Time: 7.200E+01 Max CFL: 6.955E-01 Max Diff. No.: -1.000E+00 Norm: 1.0148E-02
Iteration: 7400 Time: 7.400E+01 Max CFL: 6.955E-01 Max Diff. No.: -1.000E+00 Norm: 1.1495E-02
Iteration: 7600 Time: 7.600E+01 Max CFL: 6.956E-01 Max Diff. No.: -1.000E+00 Norm: 1.0026E-02
Iteration: 7800 Time: 7.800E+01 Max CFL: 6.956E-01 Max Diff. No.: -1.000E+00 Norm: 1.1869E-02
Iteration: 8000 Time: 8.000E+01 Max CFL: 6.957E-01 Max Diff. No.: -1.000E+00 Norm: 1.1001E-02
Writing solution file op_00004.bin.
Iteration: 8200 Time: 8.200E+01 Max CFL: 6.957E-01 Max Diff. No.: -1.000E+00 Norm: 1.1812E-02
Iteration: 8400 Time: 8.400E+01 Max CFL: 6.958E-01 Max Diff. No.: -1.000E+00 Norm: 1.1184E-02
Iteration: 8600 Time: 8.600E+01 Max CFL: 6.958E-01 Max Diff. No.: -1.000E+00 Norm: 1.1299E-02
Iteration: 8800 Time: 8.800E+01 Max CFL: 6.959E-01 Max Diff. No.: -1.000E+00 Norm: 1.1730E-02
Iteration: 9000 Time: 9.000E+01 Max CFL: 6.959E-01 Max Diff. No.: -1.000E+00 Norm: 8.9314E-03
Iteration: 9200 Time: 9.200E+01 Max CFL: 6.960E-01 Max Diff. No.: -1.000E+00 Norm: 1.0053E-02
Iteration: 9400 Time: 9.400E+01 Max CFL: 6.960E-01 Max Diff. No.: -1.000E+00 Norm: 9.9710E-03
Iteration: 9600 Time: 9.600E+01 Max CFL: 6.960E-01 Max Diff. No.: -1.000E+00 Norm: 8.7663E-03
Iteration: 9800 Time: 9.800E+01 Max CFL: 6.961E-01 Max Diff. No.: -1.000E+00 Norm: 9.2397E-03
Iteration: 10000 Time: 1.000E+02 Max CFL: 6.961E-01 Max Diff. No.: -1.000E+00 Norm: 1.1179E-02
Writing solution file op_00005.bin.
Iteration: 10200 Time: 1.020E+02 Max CFL: 6.962E-01 Max Diff. No.: -1.000E+00 Norm: 6.9634E-03
Iteration: 10400 Time: 1.040E+02 Max CFL: 6.962E-01 Max Diff. No.: -1.000E+00 Norm: 1.0934E-02
Iteration: 10600 Time: 1.060E+02 Max CFL: 6.962E-01 Max Diff. No.: -1.000E+00 Norm: 1.0509E-02
Iteration: 10800 Time: 1.080E+02 Max CFL: 6.963E-01 Max Diff. No.: -1.000E+00 Norm: 9.9973E-03
Iteration: 11000 Time: 1.100E+02 Max CFL: 6.963E-01 Max Diff. No.: -1.000E+00 Norm: 1.1854E-02
Iteration: 11200 Time: 1.120E+02 Max CFL: 6.963E-01 Max Diff. No.: -1.000E+00 Norm: 1.1314E-02
Iteration: 11400 Time: 1.140E+02 Max CFL: 6.964E-01 Max Diff. No.: -1.000E+00 Norm: 1.0630E-02
Iteration: 11600 Time: 1.160E+02 Max CFL: 6.964E-01 Max Diff. No.: -1.000E+00 Norm: 1.1933E-02
Iteration: 11800 Time: 1.180E+02 Max CFL: 6.964E-01 Max Diff. No.: -1.000E+00 Norm: 1.0585E-02
Iteration: 12000 Time: 1.200E+02 Max CFL: 6.965E-01 Max Diff. No.: -1.000E+00 Norm: 1.1356E-02
Writing solution file op_00006.bin.
Iteration: 12200 Time: 1.220E+02 Max CFL: 6.965E-01 Max Diff. No.: -1.000E+00 Norm: 1.0765E-02
Iteration: 12400 Time: 1.240E+02 Max CFL: 6.965E-01 Max Diff. No.: -1.000E+00 Norm: 8.7541E-03
Iteration: 12600 Time: 1.260E+02 Max CFL: 6.965E-01 Max Diff. No.: -1.000E+00 Norm: 1.0741E-02
Iteration: 12800 Time: 1.280E+02 Max CFL: 6.966E-01 Max Diff. No.: -1.000E+00 Norm: 7.8668E-03
Iteration: 13000 Time: 1.300E+02 Max CFL: 6.966E-01 Max Diff. No.: -1.000E+00 Norm: 9.9230E-03
Iteration: 13200 Time: 1.320E+02 Max CFL: 6.966E-01 Max Diff. No.: -1.000E+00 Norm: 1.0421E-02
Iteration: 13400 Time: 1.340E+02 Max CFL: 6.966E-01 Max Diff. No.: -1.000E+00 Norm: 9.6055E-03
Iteration: 13600 Time: 1.360E+02 Max CFL: 6.966E-01 Max Diff. No.: -1.000E+00 Norm: 9.7125E-03
Iteration: 13800 Time: 1.380E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 1.1042E-02
Iteration: 14000 Time: 1.400E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 8.9203E-03
Writing solution file op_00007.bin.
Iteration: 14200 Time: 1.420E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 1.2056E-02
Iteration: 14400 Time: 1.440E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 1.0913E-02
Iteration: 14600 Time: 1.460E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 1.0782E-02
Iteration: 14800 Time: 1.480E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 1.1767E-02
Iteration: 15000 Time: 1.500E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 9.7494E-03
Iteration: 15200 Time: 1.520E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 1.1396E-02
Iteration: 15400 Time: 1.540E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 1.0354E-02
Iteration: 15600 Time: 1.560E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 9.7353E-03
Iteration: 15800 Time: 1.580E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 1.0257E-02
Iteration: 16000 Time: 1.600E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 9.3197E-03
Writing solution file op_00008.bin.
Iteration: 16200 Time: 1.620E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 8.5178E-03
Iteration: 16400 Time: 1.640E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 1.0796E-02
Iteration: 16600 Time: 1.660E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 9.3926E-03
Iteration: 16800 Time: 1.680E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 1.0289E-02
Iteration: 17000 Time: 1.700E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 1.0896E-02
Iteration: 17200 Time: 1.720E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 1.0072E-02
Iteration: 17400 Time: 1.740E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 1.0509E-02
Iteration: 17600 Time: 1.760E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 1.1708E-02
Iteration: 17800 Time: 1.780E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 1.0457E-02
Iteration: 18000 Time: 1.800E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 1.1534E-02
Writing solution file op_00009.bin.
Iteration: 18200 Time: 1.820E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 1.0566E-02
Iteration: 18400 Time: 1.840E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 9.7589E-03
Iteration: 18600 Time: 1.860E+02 Max CFL: 6.968E-01 Max Diff. No.: -1.000E+00 Norm: 1.0785E-02
Iteration: 18800 Time: 1.880E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 9.8248E-03
Iteration: 19000 Time: 1.900E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 1.0260E-02
Iteration: 19200 Time: 1.920E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 9.7764E-03
Iteration: 19400 Time: 1.940E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 9.4451E-03
Iteration: 19600 Time: 1.960E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 9.2906E-03
Iteration: 19800 Time: 1.980E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 9.9346E-03
Iteration: 20000 Time: 2.000E+02 Max CFL: 6.967E-01 Max Diff. No.: -1.000E+00 Norm: 1.0322E-02
Writing solution file op_00010.bin.
Completed time integration (Final time: 200.000000).
Computed errors:
L1 Error : 0.0000000000000000E+00
L2 Error : 0.0000000000000000E+00
Linfinity Error : 0.0000000000000000E+00
Conservation Errors:
0.0000000000000000E+00
0.0000000000000000E+00
0.0000000000000000E+00
0.0000000000000000E+00
0.0000000000000000E+00
Solver runtime (in seconds): 3.6566252696000003E+04
Total runtime (in seconds): 3.6567538389000001E+04
Deallocating arrays.
Finished.