diff --git a/bfps/NavierStokes.py b/bfps/NavierStokes.py
index 5f3182d3056e829333f202197c6e11a4cbd290e2..9d0a770dfeb520a25e2e971c938df0a8c181fd6d 100644
--- a/bfps/NavierStokes.py
+++ b/bfps/NavierStokes.py
@@ -41,15 +41,17 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base):
fluid_precision = 'single',
fftw_plan_rigor = 'FFTW_MEASURE',
frozen_fields = False,
- use_fftw_wisdom = True):
+ use_fftw_wisdom = True,
+ QR_stats_on = False):
+ self.QR_stats_on = QR_stats_on
+ self.frozen_fields = frozen_fields
+ self.fftw_plan_rigor = fftw_plan_rigor
super(NavierStokes, self).__init__(
name = name,
work_dir = work_dir,
simname = simname,
dtype = fluid_precision,
use_fftw_wisdom = use_fftw_wisdom)
- self.frozen_fields = frozen_fields
- self.fftw_plan_rigor = fftw_plan_rigor
self.file_datasets_grow = """
//begincpp
std::string temp_string;
@@ -119,17 +121,22 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base):
//begincpp
double *velocity_moments = new double[10*4];
double *vorticity_moments = new double[10*4];
- double *trS2_Q_R_moments = new double[10*3];
ptrdiff_t *hist_velocity = new ptrdiff_t[histogram_bins*4];
ptrdiff_t *hist_vorticity = new ptrdiff_t[histogram_bins*4];
- ptrdiff_t *hist_trS2_Q_R = new ptrdiff_t[histogram_bins*3];
- ptrdiff_t *hist_QR2D = new ptrdiff_t[QR2D_histogram_bins*QR2D_histogram_bins];
double max_estimates[4];
fs->compute_velocity(fs->cvorticity);
double *spec_velocity = new double[fs->nshells*9];
double *spec_vorticity = new double[fs->nshells*9];
fs->cospectrum(fs->cvelocity, fs->cvelocity, spec_velocity);
fs->cospectrum(fs->cvorticity, fs->cvorticity, spec_vorticity);
+ //endcpp
+ """
+ if self.QR_stats_on:
+ self.stat_src += """
+ //begincpp
+ double *trS2_Q_R_moments = new double[10*3];
+ ptrdiff_t *hist_trS2_Q_R = new ptrdiff_t[histogram_bins*3];
+ ptrdiff_t *hist_QR2D = new ptrdiff_t[QR2D_histogram_bins*QR2D_histogram_bins];
max_estimates[0] = max_trS2_estimate;
max_estimates[1] = max_Q_estimate;
max_estimates[2] = max_R_estimate;
@@ -141,6 +148,10 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base):
max_estimates,
histogram_bins,
QR2D_histogram_bins);
+ //endcpp
+ """
+ self.stat_src += """
+ //begincpp
fs->ift_velocity();
max_estimates[0] = max_velocity_estimate/sqrt(3);
max_estimates[1] = max_estimates[0];
@@ -202,14 +213,6 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base):
self.stat_src += self.create_stat_output(
'/statistics/moments/vorticity',
'vorticity_moments')
- self.stat_src += self.create_stat_output(
- '/statistics/moments/trS2_Q_R',
- 'trS2_Q_R_moments',
- size_setup ="""
- count[0] = 1;
- count[1] = 10;
- count[2] = 3;
- """)
self.stat_src += self.create_stat_output(
'/statistics/spectra/velocity_velocity',
'spec_velocity',
@@ -237,24 +240,33 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base):
'/statistics/histograms/vorticity',
'hist_vorticity',
data_type = 'H5T_NATIVE_INT64')
- self.stat_src += self.create_stat_output(
- '/statistics/histograms/trS2_Q_R',
- 'hist_trS2_Q_R',
- data_type = 'H5T_NATIVE_INT64',
- size_setup = """
- count[0] = 1;
- count[1] = histogram_bins;
- count[2] = 3;
- """)
- self.stat_src += self.create_stat_output(
- '/statistics/histograms/QR2D',
- 'hist_QR2D',
- data_type = 'H5T_NATIVE_INT64',
- size_setup = """
- count[0] = 1;
- count[1] = QR2D_histogram_bins;
- count[2] = QR2D_histogram_bins;
- """)
+ if self.QR_stats_on:
+ self.stat_src += self.create_stat_output(
+ '/statistics/moments/trS2_Q_R',
+ 'trS2_Q_R_moments',
+ size_setup ="""
+ count[0] = 1;
+ count[1] = 10;
+ count[2] = 3;
+ """)
+ self.stat_src += self.create_stat_output(
+ '/statistics/histograms/trS2_Q_R',
+ 'hist_trS2_Q_R',
+ data_type = 'H5T_NATIVE_INT64',
+ size_setup = """
+ count[0] = 1;
+ count[1] = histogram_bins;
+ count[2] = 3;
+ """)
+ self.stat_src += self.create_stat_output(
+ '/statistics/histograms/QR2D',
+ 'hist_QR2D',
+ data_type = 'H5T_NATIVE_INT64',
+ size_setup = """
+ count[0] = 1;
+ count[1] = QR2D_histogram_bins;
+ count[2] = QR2D_histogram_bins;
+ """)
self.stat_src += """
//begincpp
}
@@ -262,9 +274,14 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base):
delete[] spec_vorticity;
delete[] velocity_moments;
delete[] vorticity_moments;
- delete[] trS2_Q_R_moments;
delete[] hist_velocity;
delete[] hist_vorticity;
+ //endcpp
+ """
+ if self.QR_stats_on:
+ self.stat_src += """
+ //begincpp
+ delete[] trS2_Q_R_moments;
delete[] hist_trS2_Q_R;
delete[] hist_QR2D;
//endcpp
diff --git a/bfps/fluid_base.py b/bfps/fluid_base.py
index 54b5737df792015de8d0cddd0565c56a759ac978..f29e99e5f846188085c4a9e5c84cc7d056c4a9b2 100644
--- a/bfps/fluid_base.py
+++ b/bfps/fluid_base.py
@@ -60,6 +60,7 @@ class fluid_particle_base(bfps.code):
elif self.rtype == np.float64:
self.ctype = np.dtype(np.complex128)
self.C_dtype = 'double'
+ self.parameters['dealias_type'] = 1
self.parameters['dkx'] = 1.0
self.parameters['dky'] = 1.0
self.parameters['dkz'] = 1.0
@@ -76,13 +77,12 @@ class fluid_particle_base(bfps.code):
self.parameters['fk1'] = 3.0
self.parameters['forcing_type'] = 'linear'
self.parameters['histogram_bins'] = 256
- self.parameters['QR2D_histogram_bins'] = 64
self.parameters['max_velocity_estimate'] = 1.0
self.parameters['max_vorticity_estimate'] = 1.0
+ self.parameters['QR2D_histogram_bins'] = 64
self.parameters['max_trS2_estimate'] = 1.0
self.parameters['max_Q_estimate'] = 1.0
self.parameters['max_R_estimate'] = 1.0
- self.parameters['dealias_type'] = 1
self.fluid_includes = '#include "fluid_solver.hpp"\n'
self.fluid_variables = ''
self.fluid_definitions = ''
@@ -398,42 +398,43 @@ class fluid_particle_base(bfps.code):
for k in kspace.keys():
ofile['kspace/' + k] = kspace[k]
nshells = kspace['nshell'].shape[0]
- time_chunk = 2**20//(8*3*self.parameters['histogram_bins'])
- time_chunk = max(time_chunk, 1)
- ofile.create_dataset('statistics/histograms/trS2_Q_R',
- (1,
- self.parameters['histogram_bins'],
- 3),
- chunks = (time_chunk,
- self.parameters['histogram_bins'],
- 3),
- maxshape = (None,
- self.parameters['histogram_bins'],
- 3),
- dtype = np.int64,
- compression = 'gzip')
- time_chunk = 2**20//(8*3*10)
- time_chunk = max(time_chunk, 1)
- a = ofile.create_dataset('statistics/moments/trS2_Q_R',
- (1, 10, 3),
- chunks = (time_chunk, 10, 3),
- maxshape = (None, 10, 3),
- dtype = np.float64,
- compression = 'gzip')
- time_chunk = 2**20//(8*self.parameters['QR2D_histogram_bins']**2)
- time_chunk = max(time_chunk, 1)
- ofile.create_dataset('statistics/histograms/QR2D',
- (1,
- self.parameters['QR2D_histogram_bins'],
- self.parameters['QR2D_histogram_bins']),
- chunks = (time_chunk,
- self.parameters['QR2D_histogram_bins'],
- self.parameters['QR2D_histogram_bins']),
- maxshape = (None,
- self.parameters['QR2D_histogram_bins'],
- self.parameters['QR2D_histogram_bins']),
- dtype = np.int64,
- compression = 'gzip')
+ if self.QR_stats_on:
+ time_chunk = 2**20//(8*3*self.parameters['histogram_bins'])
+ time_chunk = max(time_chunk, 1)
+ ofile.create_dataset('statistics/histograms/trS2_Q_R',
+ (1,
+ self.parameters['histogram_bins'],
+ 3),
+ chunks = (time_chunk,
+ self.parameters['histogram_bins'],
+ 3),
+ maxshape = (None,
+ self.parameters['histogram_bins'],
+ 3),
+ dtype = np.int64,
+ compression = 'gzip')
+ time_chunk = 2**20//(8*3*10)
+ time_chunk = max(time_chunk, 1)
+ a = ofile.create_dataset('statistics/moments/trS2_Q_R',
+ (1, 10, 3),
+ chunks = (time_chunk, 10, 3),
+ maxshape = (None, 10, 3),
+ dtype = np.float64,
+ compression = 'gzip')
+ time_chunk = 2**20//(8*self.parameters['QR2D_histogram_bins']**2)
+ time_chunk = max(time_chunk, 1)
+ ofile.create_dataset('statistics/histograms/QR2D',
+ (1,
+ self.parameters['QR2D_histogram_bins'],
+ self.parameters['QR2D_histogram_bins']),
+ chunks = (time_chunk,
+ self.parameters['QR2D_histogram_bins'],
+ self.parameters['QR2D_histogram_bins']),
+ maxshape = (None,
+ self.parameters['QR2D_histogram_bins'],
+ self.parameters['QR2D_histogram_bins']),
+ dtype = np.int64,
+ compression = 'gzip')
for k in ['velocity', 'vorticity']:
time_chunk = 2**20//(8*3*self.parameters['nx'])
time_chunk = max(time_chunk, 1)