Skip to content
Snippets Groups Projects
Commit 0a02b749 authored by Cristian Lalescu's avatar Cristian Lalescu
Browse files

make QR stats optional

parent 8342a247
No related branches found
No related tags found
No related merge requests found
...@@ -41,15 +41,17 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base): ...@@ -41,15 +41,17 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base):
fluid_precision = 'single', fluid_precision = 'single',
fftw_plan_rigor = 'FFTW_MEASURE', fftw_plan_rigor = 'FFTW_MEASURE',
frozen_fields = False, 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__( super(NavierStokes, self).__init__(
name = name, name = name,
work_dir = work_dir, work_dir = work_dir,
simname = simname, simname = simname,
dtype = fluid_precision, dtype = fluid_precision,
use_fftw_wisdom = use_fftw_wisdom) use_fftw_wisdom = use_fftw_wisdom)
self.frozen_fields = frozen_fields
self.fftw_plan_rigor = fftw_plan_rigor
self.file_datasets_grow = """ self.file_datasets_grow = """
//begincpp //begincpp
std::string temp_string; std::string temp_string;
...@@ -119,17 +121,22 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base): ...@@ -119,17 +121,22 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base):
//begincpp //begincpp
double *velocity_moments = new double[10*4]; double *velocity_moments = new double[10*4];
double *vorticity_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_velocity = new ptrdiff_t[histogram_bins*4];
ptrdiff_t *hist_vorticity = 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]; double max_estimates[4];
fs->compute_velocity(fs->cvorticity); fs->compute_velocity(fs->cvorticity);
double *spec_velocity = new double[fs->nshells*9]; double *spec_velocity = new double[fs->nshells*9];
double *spec_vorticity = new double[fs->nshells*9]; double *spec_vorticity = new double[fs->nshells*9];
fs->cospectrum(fs->cvelocity, fs->cvelocity, spec_velocity); fs->cospectrum(fs->cvelocity, fs->cvelocity, spec_velocity);
fs->cospectrum(fs->cvorticity, fs->cvorticity, spec_vorticity); 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[0] = max_trS2_estimate;
max_estimates[1] = max_Q_estimate; max_estimates[1] = max_Q_estimate;
max_estimates[2] = max_R_estimate; max_estimates[2] = max_R_estimate;
...@@ -141,6 +148,10 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base): ...@@ -141,6 +148,10 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base):
max_estimates, max_estimates,
histogram_bins, histogram_bins,
QR2D_histogram_bins); QR2D_histogram_bins);
//endcpp
"""
self.stat_src += """
//begincpp
fs->ift_velocity(); fs->ift_velocity();
max_estimates[0] = max_velocity_estimate/sqrt(3); max_estimates[0] = max_velocity_estimate/sqrt(3);
max_estimates[1] = max_estimates[0]; max_estimates[1] = max_estimates[0];
...@@ -202,14 +213,6 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base): ...@@ -202,14 +213,6 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base):
self.stat_src += self.create_stat_output( self.stat_src += self.create_stat_output(
'/statistics/moments/vorticity', '/statistics/moments/vorticity',
'vorticity_moments') '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( self.stat_src += self.create_stat_output(
'/statistics/spectra/velocity_velocity', '/statistics/spectra/velocity_velocity',
'spec_velocity', 'spec_velocity',
...@@ -237,24 +240,33 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base): ...@@ -237,24 +240,33 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base):
'/statistics/histograms/vorticity', '/statistics/histograms/vorticity',
'hist_vorticity', 'hist_vorticity',
data_type = 'H5T_NATIVE_INT64') data_type = 'H5T_NATIVE_INT64')
self.stat_src += self.create_stat_output( if self.QR_stats_on:
'/statistics/histograms/trS2_Q_R', self.stat_src += self.create_stat_output(
'hist_trS2_Q_R', '/statistics/moments/trS2_Q_R',
data_type = 'H5T_NATIVE_INT64', 'trS2_Q_R_moments',
size_setup = """ size_setup ="""
count[0] = 1; count[0] = 1;
count[1] = histogram_bins; count[1] = 10;
count[2] = 3; count[2] = 3;
""") """)
self.stat_src += self.create_stat_output( self.stat_src += self.create_stat_output(
'/statistics/histograms/QR2D', '/statistics/histograms/trS2_Q_R',
'hist_QR2D', 'hist_trS2_Q_R',
data_type = 'H5T_NATIVE_INT64', data_type = 'H5T_NATIVE_INT64',
size_setup = """ size_setup = """
count[0] = 1; count[0] = 1;
count[1] = QR2D_histogram_bins; count[1] = histogram_bins;
count[2] = QR2D_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 += """ self.stat_src += """
//begincpp //begincpp
} }
...@@ -262,9 +274,14 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base): ...@@ -262,9 +274,14 @@ class NavierStokes(bfps.fluid_base.fluid_particle_base):
delete[] spec_vorticity; delete[] spec_vorticity;
delete[] velocity_moments; delete[] velocity_moments;
delete[] vorticity_moments; delete[] vorticity_moments;
delete[] trS2_Q_R_moments;
delete[] hist_velocity; delete[] hist_velocity;
delete[] hist_vorticity; 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_trS2_Q_R;
delete[] hist_QR2D; delete[] hist_QR2D;
//endcpp //endcpp
......
...@@ -60,6 +60,7 @@ class fluid_particle_base(bfps.code): ...@@ -60,6 +60,7 @@ class fluid_particle_base(bfps.code):
elif self.rtype == np.float64: elif self.rtype == np.float64:
self.ctype = np.dtype(np.complex128) self.ctype = np.dtype(np.complex128)
self.C_dtype = 'double' self.C_dtype = 'double'
self.parameters['dealias_type'] = 1
self.parameters['dkx'] = 1.0 self.parameters['dkx'] = 1.0
self.parameters['dky'] = 1.0 self.parameters['dky'] = 1.0
self.parameters['dkz'] = 1.0 self.parameters['dkz'] = 1.0
...@@ -76,13 +77,12 @@ class fluid_particle_base(bfps.code): ...@@ -76,13 +77,12 @@ class fluid_particle_base(bfps.code):
self.parameters['fk1'] = 3.0 self.parameters['fk1'] = 3.0
self.parameters['forcing_type'] = 'linear' self.parameters['forcing_type'] = 'linear'
self.parameters['histogram_bins'] = 256 self.parameters['histogram_bins'] = 256
self.parameters['QR2D_histogram_bins'] = 64
self.parameters['max_velocity_estimate'] = 1.0 self.parameters['max_velocity_estimate'] = 1.0
self.parameters['max_vorticity_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_trS2_estimate'] = 1.0
self.parameters['max_Q_estimate'] = 1.0 self.parameters['max_Q_estimate'] = 1.0
self.parameters['max_R_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_includes = '#include "fluid_solver.hpp"\n'
self.fluid_variables = '' self.fluid_variables = ''
self.fluid_definitions = '' self.fluid_definitions = ''
...@@ -398,42 +398,43 @@ class fluid_particle_base(bfps.code): ...@@ -398,42 +398,43 @@ class fluid_particle_base(bfps.code):
for k in kspace.keys(): for k in kspace.keys():
ofile['kspace/' + k] = kspace[k] ofile['kspace/' + k] = kspace[k]
nshells = kspace['nshell'].shape[0] nshells = kspace['nshell'].shape[0]
time_chunk = 2**20//(8*3*self.parameters['histogram_bins']) if self.QR_stats_on:
time_chunk = max(time_chunk, 1) time_chunk = 2**20//(8*3*self.parameters['histogram_bins'])
ofile.create_dataset('statistics/histograms/trS2_Q_R', time_chunk = max(time_chunk, 1)
(1, ofile.create_dataset('statistics/histograms/trS2_Q_R',
self.parameters['histogram_bins'], (1,
3), self.parameters['histogram_bins'],
chunks = (time_chunk, 3),
self.parameters['histogram_bins'], chunks = (time_chunk,
3), self.parameters['histogram_bins'],
maxshape = (None, 3),
self.parameters['histogram_bins'], maxshape = (None,
3), self.parameters['histogram_bins'],
dtype = np.int64, 3),
compression = 'gzip') dtype = np.int64,
time_chunk = 2**20//(8*3*10) compression = 'gzip')
time_chunk = max(time_chunk, 1) time_chunk = 2**20//(8*3*10)
a = ofile.create_dataset('statistics/moments/trS2_Q_R', time_chunk = max(time_chunk, 1)
(1, 10, 3), a = ofile.create_dataset('statistics/moments/trS2_Q_R',
chunks = (time_chunk, 10, 3), (1, 10, 3),
maxshape = (None, 10, 3), chunks = (time_chunk, 10, 3),
dtype = np.float64, maxshape = (None, 10, 3),
compression = 'gzip') dtype = np.float64,
time_chunk = 2**20//(8*self.parameters['QR2D_histogram_bins']**2) compression = 'gzip')
time_chunk = max(time_chunk, 1) time_chunk = 2**20//(8*self.parameters['QR2D_histogram_bins']**2)
ofile.create_dataset('statistics/histograms/QR2D', time_chunk = max(time_chunk, 1)
(1, ofile.create_dataset('statistics/histograms/QR2D',
self.parameters['QR2D_histogram_bins'], (1,
self.parameters['QR2D_histogram_bins']), self.parameters['QR2D_histogram_bins'],
chunks = (time_chunk, 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'], maxshape = (None,
self.parameters['QR2D_histogram_bins']), self.parameters['QR2D_histogram_bins'],
dtype = np.int64, self.parameters['QR2D_histogram_bins']),
compression = 'gzip') dtype = np.int64,
compression = 'gzip')
for k in ['velocity', 'vorticity']: for k in ['velocity', 'vorticity']:
time_chunk = 2**20//(8*3*self.parameters['nx']) time_chunk = 2**20//(8*3*self.parameters['nx'])
time_chunk = max(time_chunk, 1) time_chunk = max(time_chunk, 1)
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment