See 2D (1D-1V) Vlasov Equation - Two-Stream Instability to familiarize yourself with this case. This example uses a DMD object that has already been trained (see 2D (1D-1V) Vlasov Equation - Two-Stream Instability (Time-Windowed DMD)).
Location: hypar/Examples/2D/Vlasov1D1V/TwoStreamInstability_libROM_DMD_Predict (This directory contains all the input files needed to run this case.)
Governing equations: 2D (1D-1V) Vlasov Equation (vlasov.h)
Reduced Order Modeling: This example predicts the solution from a trained DMD object. The code does not solve the PDE by discretizing in space and integrating in time.
Domain:
Initial solution: \(f\left(x,v\right) = \frac{4}{\pi T}\left(1+\frac{1}{10}\cos\left(2k\pi\frac{x}{L}\right)\right)\left[\exp\left(-\frac{\left(v-2\right)^2}{2T}\right) + \exp\left(-\frac{\left(v+2\right)^2}{2T}\right)\right]\), \(k=1,T=1,L=2\pi\).
Self-consistent electric field is computed by solving the Poisson equation in a periodic domain using Fourier transforms. This examples requires HyPar to be compiled with FFTW (http://www.fftw.org/).
Reduced Order Modeling:
Note:
In this mode, HyPar will run just like an usual PDE simulation, except that it will swap out the numerical spatial discretization and time integration with the ROM-based prediction. The input files and output files will be the same as a regular simulation with the following comments:
- Numerical method inputs are ignored (eg. those that specify spatial discretization and time integration methods).
- In solver.inp, the values for dt, n_iter, and file_op_iter is used only to compute the simulation times at which to compute and write the solution. The time step size, dt, need not respect any CFL criterion.
- HyPar::ConservationCheck is set to "no" since it is not possible to compute conservation loss for a general domain (because boundary fluxes are not being computed).
Input files required:
librom.inp
begin
mode predict
dmd_dirname DMD
end
DMD Object(s) :
The trained DMD object(s) must be located in the directory specified in librom.inp as dmd_dirname (DMDROMObject::m_dirname). For this example, they were generated using 2D (1D-1V) Vlasov Equation - Two-Stream Instability (Time-Windowed DMD).
solver.inp
begin
ndims 2
nvars 1
size 128 128
iproc 4 4
ghost 3
n_iter 25
dt 0.2
screen_op_iter 1
file_op_iter 1
op_file_format binary
ip_file_type binary
op_overwrite no
model vlasov
end
boundary.inp
4
periodic 0 1 0.0000000000000000e+00 0.0000000000000000e+00 -7.0000000000000000e+00 7.0000000000000000e+00
periodic 0 -1 0.0000000000000000e+00 0.0000000000000000e+00 -7.0000000000000000e+00 7.0000000000000000e+00
dirichlet 1 1 0.0000000000000000e+00 6.28318530717959e+000 0.0000000000000000e+00 0.0000000000000000e+00
0
dirichlet 1 -1 0.0000000000000000e+00 6.28318530717959e+000 0.0000000000000000e+00 0.0000000000000000e+00
0
physics.inp
begin
self_consistent_electric_field 1
end
To generate initial.inp (initial solution), compile and run the following code in the run directory.
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <string.h>
{
double pi = 4.0 * atan(1.0);
double two_pi = 8.0 * atan(1.0);
double L = two_pi;
double T = 1.0;
double k = 1.0;
int NI,NJ,ndims,n_iter;
FILE *in;
char ip_file_type[50];
strcpy(ip_file_type,"ascii");
printf("Reading file \"solver.inp\"...\n");
in = fopen("solver.inp","r");
if (!in) {
fprintf(stderr,"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);
} else if (!strcmp(word, "ip_file_type")) {
fscanf(in,"%s",ip_file_type);
}
}
} else {
fprintf(stderr,"Error: Illegal format in solver.inp. Crash and burn!\n");
return(0);
}
}
fclose(in);
if (ndims != 2) {
fprintf(stderr,"ndims is not 2 in solver.inp. this code is to generate 2D initial conditions\n");
return(0);
}
printf("Grid: %d, %d\n",NI,NJ);
int i,j;
double dx = two_pi / ((double)NI);
double dv = (14.0 * T) / ((double)NJ);
double start_x = 0.0;
double start_v = -7.0 * T;
double *x, *v, *f;
x = (double*) calloc (NI , sizeof(double));
v = (double*) calloc (NJ , sizeof(double));
f = (double*) calloc (NI*NJ, sizeof(double));
for (i = 0; i < NI; i++){
for (j = 0; j < NJ; j++){
x[i] = start_x + i*dx;
v[j] = start_v + j*dv;
int p = NJ*i + j;
double temp1 = cos(2.0 * k * pi * x[i] / L);
double temp2 = exp(- (v[j] -2.0) * (v[j] - 2.0) / (2.0 * T));
double temp3 = exp(- (v[j] + 2.0) * (v[j] + 2.0) / (2.0 * T));
double temp4 = sqrt(2.0 * pi * T);
f[p] = 16.0 * (1.0 + 0.1 * temp1) * 0.5 * (temp2 / temp4 + temp3 / temp4);
}
}
FILE *out;
if (!strcmp(ip_file_type,"ascii")) {
out = fopen("initial.inp","w");
for (i = 0; i < NI; i++) fprintf(out,"%1.16e ",x[i]);
fprintf(out,"\n");
for (j = 0; j < NJ; j++) fprintf(out,"%1.16e ",v[j]);
fprintf(out,"\n");
for (j = 0; j < NJ; j++) {
for (i = 0; i < NI; i++) {
int p = NJ*i + j;
fprintf(out,"%1.16e ", f[p]);
}
}
fprintf(out,"\n");
fclose(out);
} else if ((!strcmp(ip_file_type,"binary")) || (!strcmp(ip_file_type,"bin"))) {
printf("Writing binary exact solution file initial.inp\n");
out = fopen("initial.inp","wb");
fwrite(x, sizeof(double),NI,out);
fwrite(v, sizeof(double),NJ,out);
double *F = (double*) calloc (NI*NJ,sizeof(double));
for (i=0; i < NI; i++) {
for (j = 0; j < NJ; j++) {
int p = NJ*i + j;
int q = NI*j + i;
F[q+0] = f[p];
}
}
fwrite(F, sizeof(double),NI*NJ,out);
free(F);
fclose(out);
}
free(x);
free(v);
free(f);
return(0);
}
Output:
Note that iproc is set to
4 4
in solver.inp (i.e., 4 processors along x, and 4 processors along y). Thus, this example should be run with 16 MPI ranks (or change iproc).
After running the code, there should be the following output files:
- 26 output files op_00000.bin, op_00001.bin, ... op_00025.bin; these are the predicted solutions from the DMD object(s).
The first of these file sets is the solution at \(t=0\) and the final one is the solution at \(t=5\). Since HyPar::op_overwrite is set to no in solver.inp, a separate file is written for solutions at each output time. All the files are binary (HyPar::op_file_format is set to binary in solver.inp).
The provided Python script (Examples/Python/plotSolution_2DBinary.py) can be used to generate plots from the binary files that compare the HyPar and DMD solutions. Alternatively, HyPar::op_file_format can be set to tecplot2d or text, and Tecplot/VisIt or something similar can be used to plot the resulting files.
The following plot shows the evolution of the distribution function:
Expected screen output:
HyPar - Parallel (MPI) version with 16 processes
Compiled with PETSc time integration.
Allocated simulation object(s).
Reading solver inputs from file "solver.inp".
No. of dimensions : 2
No. of variables : 1
Domain size : 128 128
Processes along each dimension : 4 4
Exact solution domain size : 128 128
No. of ghosts pts : 3
No. of iter. : 25
Restart iteration : 0
Time integration scheme : euler
Spatial discretization scheme (hyperbolic) : 1
Split hyperbolic flux term? : no
Interpolation type for hyperbolic term : characteristic
Spatial discretization type (parabolic ) : nonconservative-1stage
Spatial discretization scheme (parabolic ) : 2
Time Step : 2.000000E-01
Check for conservation : no
Screen output iterations : 1
File output iterations : 1
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 : vlasov
Partitioning domain and allocating data arrays.
Reading array from binary file initial.inp (Serial mode).
Volume integral of the initial solution:
0: 1.0053093535576996E+02
Reading boundary conditions from boundary.inp.
Boundary periodic: Along dimension 0 and face +1
Boundary periodic: Along dimension 0 and face -1
Boundary dirichlet: Along dimension 1 and face +1
Boundary dirichlet: Along dimension 1 and face -1
4 boundary condition(s) read.
Initializing solvers.
Initializing physics. Model = "vlasov"
Reading physical model inputs from file "physics.inp".
Setting up libROM interface.
libROM inputs and parameters:
reduced model dimensionality: 1646581949
sampling frequency: 1076887551
mode: predict
type: DMD
save to file: true
libROM DMD inputs:
number of samples per window: 2147483647
directory name for DMD onjects: DMD
libROMInterface::loadROM() - loading ROM objects.
Loading DMD object (DMD/dmdobj_0000), time window=[0.00e+00,1.00e+00].
Loading DMD object (DMD/dmdobj_0001), time window=[1.00e+00,2.00e+00].
Loading DMD object (DMD/dmdobj_0002), time window=[2.00e+00,3.00e+00].
Loading DMD object (DMD/dmdobj_0003), time window=[3.00e+00,4.00e+00].
Loading DMD object (DMD/dmdobj_0004), time window=[4.00e+00,-1.00e+00].
libROM: Predicted solution at time 0.0000e+00 using ROM, wallclock time: 0.013092.
Writing solution file op_00000.bin.
libROM: Predicted solution at time 2.0000e-01 using ROM, wallclock time: 0.152247.
Writing solution file op_00001.bin.
libROM: Predicted solution at time 4.0000e-01 using ROM, wallclock time: 0.149397.
Writing solution file op_00002.bin.
libROM: Predicted solution at time 6.0000e-01 using ROM, wallclock time: 0.216336.
Writing solution file op_00003.bin.
libROM: Predicted solution at time 8.0000e-01 using ROM, wallclock time: 0.176951.
Writing solution file op_00004.bin.
libROM: Predicted solution at time 1.0000e+00 using ROM, wallclock time: 0.227423.
Writing solution file op_00005.bin.
libROM: Predicted solution at time 1.2000e+00 using ROM, wallclock time: 0.195707.
Writing solution file op_00006.bin.
libROM: Predicted solution at time 1.4000e+00 using ROM, wallclock time: 0.089968.
Writing solution file op_00007.bin.
libROM: Predicted solution at time 1.6000e+00 using ROM, wallclock time: 0.241131.
Writing solution file op_00008.bin.
libROM: Predicted solution at time 1.8000e+00 using ROM, wallclock time: 0.189393.
Writing solution file op_00009.bin.
libROM: Predicted solution at time 2.0000e+00 using ROM, wallclock time: 0.213828.
Writing solution file op_00010.bin.
libROM: Predicted solution at time 2.2000e+00 using ROM, wallclock time: 0.161835.
Writing solution file op_00011.bin.
libROM: Predicted solution at time 2.4000e+00 using ROM, wallclock time: 0.296023.
Writing solution file op_00012.bin.
libROM: Predicted solution at time 2.6000e+00 using ROM, wallclock time: 0.378588.
Writing solution file op_00013.bin.
libROM: Predicted solution at time 2.8000e+00 using ROM, wallclock time: 0.191621.
Writing solution file op_00014.bin.
libROM: Predicted solution at time 3.0000e+00 using ROM, wallclock time: 0.231661.
Writing solution file op_00015.bin.
libROM: Predicted solution at time 3.2000e+00 using ROM, wallclock time: 0.198773.
Writing solution file op_00016.bin.
libROM: Predicted solution at time 3.4000e+00 using ROM, wallclock time: 0.176326.
Writing solution file op_00017.bin.
libROM: Predicted solution at time 3.6000e+00 using ROM, wallclock time: 0.256226.
Writing solution file op_00018.bin.
libROM: Predicted solution at time 3.8000e+00 using ROM, wallclock time: 0.217788.
Writing solution file op_00019.bin.
libROM: Predicted solution at time 4.0000e+00 using ROM, wallclock time: 0.195570.
Writing solution file op_00020.bin.
libROM: Predicted solution at time 4.2000e+00 using ROM, wallclock time: 0.147288.
Writing solution file op_00021.bin.
libROM: Predicted solution at time 4.4000e+00 using ROM, wallclock time: 0.155655.
Writing solution file op_00022.bin.
libROM: Predicted solution at time 4.6000e+00 using ROM, wallclock time: 0.197037.
Writing solution file op_00023.bin.
libROM: Predicted solution at time 4.8000e+00 using ROM, wallclock time: 0.175624.
Writing solution file op_00024.bin.
libROM: Predicted solution at time 5.0000e+00 using ROM, wallclock time: 0.136981.
Writing solution file op_00025.bin.
libROM: total prediction/query wallclock time: 4.982469 (seconds).
Solver runtime (in seconds): 1.1819117000000000E+01
Total runtime (in seconds): 1.2562611000000000E+01
Deallocating arrays.
Finished.
