Commit 3c3a23a0 authored by Cristian Lalescu's avatar Cristian Lalescu
Browse files

tweak statistics cache

basic stats can now be read even in the absence of the main simname.h5
file, using just the cache file.
parent 39a161c2
Pipeline #30110 passed with stage
in 15 minutes and 7 seconds
...@@ -233,76 +233,102 @@ class DNS(_code): ...@@ -233,76 +233,102 @@ class DNS(_code):
""" """
if len(list(self.statistics.keys())) > 0: if len(list(self.statistics.keys())) > 0:
return None return None
self.read_parameters() if not os.path.exists(self.get_data_file_name()):
with self.get_data_file() as data_file: if os.path.exists(self.get_cache_file_name()):
if 'moments' not in data_file['statistics'].keys(): self.read_parameters(fname = self.get_cache_file_name())
return None with self.get_cache_file() as pp_file:
iter0 = min((data_file['statistics/moments/velocity'].shape[0] * for k in ['t',
self.parameters['niter_stat']-1), 'energy(t)',
iter0) 'energy(k)',
if type(iter1) == type(None): 'enstrophy(t)',
iter1 = data_file['iteration'].value 'enstrophy(k)',
else:
iter1 = min(data_file['iteration'].value, iter1)
ii0 = iter0 // self.parameters['niter_stat']
ii1 = iter1 // self.parameters['niter_stat']
self.statistics['kshell'] = data_file['kspace/kshell'].value
self.statistics['nshell'] = data_file['kspace/nshell'].value
for kk in [-1, -2]:
if (self.statistics['kshell'][kk] == 0):
self.statistics['kshell'][kk] = np.nan
self.statistics['kM'] = data_file['kspace/kM'].value
self.statistics['dk'] = data_file['kspace/dk'].value
computation_needed = True
pp_file = h5py.File(self.get_postprocess_file_name(), 'a')
if 'ii0' in pp_file.keys():
computation_needed = not (ii0 == pp_file['ii0'].value and
ii1 == pp_file['ii1'].value)
if computation_needed:
for k in ['t', 'vel_max(t)', 'renergy(t)',
'energy(t)', 'enstrophy(t)',
'energy(k)', 'enstrophy(k)',
'energy(t, k)',
'enstrophy(t, k)',
'R_ij(t)', 'R_ij(t)',
'ii0', 'ii1', 'iter0', 'iter1']: 'vel_max(t)',
'renergy(t)']:
if k in pp_file.keys(): if k in pp_file.keys():
del pp_file[k] self.statistics[k] = pp_file[k].value
if computation_needed: self.statistics['kM'] = pp_file['kspace/kM'].value
pp_file['iter0'] = iter0 self.statistics['dk'] = pp_file['kspace/dk'].value
pp_file['iter1'] = iter1 self.statistics['kshell'] = pp_file['kspace/kshell'].value
pp_file['ii0'] = ii0 self.statistics['nshell'] = pp_file['kspace/nshell'].value
pp_file['ii1'] = ii1 else:
pp_file['t'] = (self.parameters['dt']* self.read_parameters()
self.parameters['niter_stat']* with self.get_data_file() as data_file:
(np.arange(ii0, ii1+1).astype(np.float))) if 'moments' not in data_file['statistics'].keys():
phi_ij = data_file['statistics/spectra/velocity_velocity'][ii0:ii1+1] return None
pp_file['R_ij(t)'] = np.sum(phi_ij, axis = 1) iter0 = min((data_file['statistics/moments/velocity'].shape[0] *
energy_tk = ( self.parameters['niter_stat']-1),
phi_ij[:, :, 0, 0] + iter0)
phi_ij[:, :, 1, 1] + if type(iter1) == type(None):
phi_ij[:, :, 2, 2])/2 iter1 = data_file['iteration'].value
pp_file['energy(t)'] = np.sum(energy_tk, axis = 1) else:
pp_file['energy(k)'] = np.mean(energy_tk, axis = 0)*(4*np.pi*self.statistics['kshell']**2) / (self.statistics['dk']*self.statistics['nshell']) iter1 = min(data_file['iteration'].value, iter1)
enstrophy_tk = ( ii0 = iter0 // self.parameters['niter_stat']
data_file['statistics/spectra/vorticity_vorticity'][ii0:ii1+1, :, 0, 0] + ii1 = iter1 // self.parameters['niter_stat']
data_file['statistics/spectra/vorticity_vorticity'][ii0:ii1+1, :, 1, 1] + self.statistics['kshell'] = data_file['kspace/kshell'].value
data_file['statistics/spectra/vorticity_vorticity'][ii0:ii1+1, :, 2, 2])/2 self.statistics['nshell'] = data_file['kspace/nshell'].value
pp_file['enstrophy(t)'] = np.sum(enstrophy_tk, axis = 1) for kk in [-1, -2]:
pp_file['enstrophy(k)'] = np.mean(enstrophy_tk, axis = 0)*(4*np.pi*self.statistics['kshell']**2) / self.statistics['nshell'] if (self.statistics['kshell'][kk] == 0):
pp_file['vel_max(t)'] = data_file['statistics/moments/velocity'][ii0:ii1+1, 9, 3] self.statistics['kshell'][kk] = np.nan
pp_file['renergy(t)'] = data_file['statistics/moments/velocity'][ii0:ii1+1, 2, 3]/2 self.statistics['kM'] = data_file['kspace/kM'].value
for k in ['t', self.statistics['dk'] = data_file['kspace/dk'].value
'energy(t)', computation_needed = True
'energy(k)', pp_file = h5py.File(self.get_postprocess_file_name(), 'a')
'enstrophy(t)', if not ('parameters' in pp_file.keys()):
'enstrophy(k)', data_file.copy('parameters', pp_file)
'R_ij(t)', data_file.copy('kspace', pp_file)
'vel_max(t)', if 'ii0' in pp_file.keys():
'renergy(t)']: computation_needed = not (ii0 == pp_file['ii0'].value and
if k in pp_file.keys(): ii1 == pp_file['ii1'].value)
self.statistics[k] = pp_file[k].value if computation_needed:
self.compute_time_averages() for k in ['t', 'vel_max(t)', 'renergy(t)',
'energy(t)', 'enstrophy(t)',
'energy(k)', 'enstrophy(k)',
'energy(t, k)',
'enstrophy(t, k)',
'R_ij(t)',
'ii0', 'ii1', 'iter0', 'iter1']:
if k in pp_file.keys():
del pp_file[k]
if computation_needed:
pp_file['iter0'] = iter0
pp_file['iter1'] = iter1
pp_file['ii0'] = ii0
pp_file['ii1'] = ii1
pp_file['t'] = (self.parameters['dt']*
self.parameters['niter_stat']*
(np.arange(ii0, ii1+1).astype(np.float)))
phi_ij = data_file['statistics/spectra/velocity_velocity'][ii0:ii1+1]
pp_file['R_ij(t)'] = np.sum(phi_ij, axis = 1)
energy_tk = (
phi_ij[:, :, 0, 0] +
phi_ij[:, :, 1, 1] +
phi_ij[:, :, 2, 2])/2
pp_file['energy(t)'] = np.sum(energy_tk, axis = 1)
pp_file['energy(k)'] = np.mean(energy_tk, axis = 0)*(4*np.pi*self.statistics['kshell']**2) / (self.statistics['dk']*self.statistics['nshell'])
enstrophy_tk = (
data_file['statistics/spectra/vorticity_vorticity'][ii0:ii1+1, :, 0, 0] +
data_file['statistics/spectra/vorticity_vorticity'][ii0:ii1+1, :, 1, 1] +
data_file['statistics/spectra/vorticity_vorticity'][ii0:ii1+1, :, 2, 2])/2
pp_file['enstrophy(t)'] = np.sum(enstrophy_tk, axis = 1)
pp_file['enstrophy(k)'] = np.mean(enstrophy_tk, axis = 0)*(4*np.pi*self.statistics['kshell']**2) / (self.statistics['dk']*self.statistics['nshell'])
pp_file['vel_max(t)'] = data_file['statistics/moments/velocity'][ii0:ii1+1, 9, 3]
pp_file['renergy(t)'] = data_file['statistics/moments/velocity'][ii0:ii1+1, 2, 3]/2
for k in ['t',
'energy(t)',
'energy(k)',
'enstrophy(t)',
'enstrophy(k)',
'R_ij(t)',
'vel_max(t)',
'renergy(t)']:
if k in pp_file.keys():
self.statistics[k] = pp_file[k].value
# sanity check --- Parseval theorem check
assert(np.max(np.abs(
self.statistics['renergy(t)'] -
self.statistics['energy(t)']) / self.statistics['energy(t)']) < 1e-5)
self.compute_time_averages()
return None return None
def compute_Reynolds_stress_invariants( def compute_Reynolds_stress_invariants(
self): self):
......
...@@ -584,79 +584,105 @@ class NavierStokes(_fluid_particle_base): ...@@ -584,79 +584,105 @@ class NavierStokes(_fluid_particle_base):
""" """
if len(list(self.statistics.keys())) > 0: if len(list(self.statistics.keys())) > 0:
return None return None
self.read_parameters() if not os.path.exists(self.get_data_file_name()):
with self.get_data_file() as data_file: if os.path.exists(self.get_cache_file_name()):
if 'moments' not in data_file['statistics'].keys(): self.read_parameters(fname = self.get_cache_file_name())
return None with self.get_cache_file() as pp_file:
iter0 = min((data_file['statistics/moments/velocity'].shape[0] * for k in ['t',
self.parameters['niter_stat']-1), 'energy(t)',
iter0) 'energy(k)',
if type(iter1) == type(None): 'enstrophy(t)',
iter1 = data_file['iteration'].value 'enstrophy(k)',
else:
iter1 = min(data_file['iteration'].value, iter1)
ii0 = iter0 // self.parameters['niter_stat']
ii1 = iter1 // self.parameters['niter_stat']
self.statistics['kshell'] = data_file['kspace/kshell'].value
self.statistics['nshell'] = data_file['kspace/nshell'].value
for kk in [-1, -2]:
if (self.statistics['kshell'][kk] == 0):
self.statistics['kshell'][kk] = np.nan
self.statistics['kM'] = data_file['kspace/kM'].value
self.statistics['dk'] = data_file['kspace/dk'].value
computation_needed = True
pp_file = h5py.File(self.get_postprocess_file_name(), 'a')
if 'ii0' in pp_file.keys():
computation_needed = not (ii0 == pp_file['ii0'].value and
ii1 == pp_file['ii1'].value)
if computation_needed:
for k in ['t', 'vel_max(t)', 'renergy(t)',
'energy(t)', 'enstrophy(t)',
'energy(k)', 'enstrophy(k)',
'energy(t, k)',
'enstrophy(t, k)',
'R_ij(t)', 'R_ij(t)',
'ii0', 'ii1', 'iter0', 'iter1']: 'vel_max(t)',
'renergy(t)']:
if k in pp_file.keys(): if k in pp_file.keys():
del pp_file[k] self.statistics[k] = pp_file[k].value
if computation_needed: self.statistics['kM'] = pp_file['kspace/kM'].value
pp_file['iter0'] = iter0 self.statistics['dk'] = pp_file['kspace/dk'].value
pp_file['iter1'] = iter1 self.statistics['kshell'] = pp_file['kspace/kshell'].value
pp_file['ii0'] = ii0 self.statistics['nshell'] = pp_file['kspace/nshell'].value
pp_file['ii1'] = ii1 else:
pp_file['t'] = (self.parameters['dt']* self.read_parameters()
self.parameters['niter_stat']* with self.get_data_file() as data_file:
(np.arange(ii0, ii1+1).astype(np.float))) if 'moments' not in data_file['statistics'].keys():
phi_ij = data_file['statistics/spectra/velocity_velocity'][ii0:ii1+1] return None
pp_file['R_ij(t)'] = np.sum(phi_ij, axis = 1) iter0 = min((data_file['statistics/moments/velocity'].shape[0] *
energy_tk = ( self.parameters['niter_stat']-1),
phi_ij[:, :, 0, 0] + iter0)
phi_ij[:, :, 1, 1] + if type(iter1) == type(None):
phi_ij[:, :, 2, 2])/2 iter1 = data_file['iteration'].value
pp_file['energy(t)'] = np.sum(energy_tk, axis = 1) else:
pp_file['energy(k)'] = np.mean(energy_tk, axis = 0)*(4*np.pi*self.statistics['kshell']**2) / (self.statistics['dk']*self.statistics['nshell']) iter1 = min(data_file['iteration'].value, iter1)
enstrophy_tk = ( ii0 = iter0 // self.parameters['niter_stat']
data_file['statistics/spectra/vorticity_vorticity'][ii0:ii1+1, :, 0, 0] + ii1 = iter1 // self.parameters['niter_stat']
data_file['statistics/spectra/vorticity_vorticity'][ii0:ii1+1, :, 1, 1] + self.statistics['kshell'] = data_file['kspace/kshell'].value
data_file['statistics/spectra/vorticity_vorticity'][ii0:ii1+1, :, 2, 2])/2 self.statistics['nshell'] = data_file['kspace/nshell'].value
pp_file['enstrophy(t)'] = np.sum(enstrophy_tk, axis = 1) for kk in [-1, -2]:
pp_file['enstrophy(k)'] = np.mean(enstrophy_tk, axis = 0)*(4*np.pi*self.statistics['kshell']**2) / (self.statistics['dk']*self.statistics['nshell']) if (self.statistics['kshell'][kk] == 0):
pp_file['vel_max(t)'] = data_file['statistics/moments/velocity'][ii0:ii1+1, 9, 3] self.statistics['kshell'][kk] = np.nan
pp_file['renergy(t)'] = data_file['statistics/moments/velocity'][ii0:ii1+1, 2, 3]/2 self.statistics['kM'] = data_file['kspace/kM'].value
if 'trS2_Q_R' in data_file['statistics/moments'].keys(): self.statistics['dk'] = data_file['kspace/dk'].value
pp_file['mean_trS2(t)'] = data_file['statistics/moments/trS2_Q_R'][:, 1, 0] computation_needed = True
for k in ['t', pp_file = h5py.File(self.get_postprocess_file_name(), 'a')
'energy(t)', if not ('parameters' in pp_file.keys()):
'energy(k)', data_file.copy('parameters', pp_file)
'enstrophy(t)', data_file.copy('kspace', pp_file)
'enstrophy(k)', if 'ii0' in pp_file.keys():
'R_ij(t)', computation_needed = not (ii0 == pp_file['ii0'].value and
'vel_max(t)', ii1 == pp_file['ii1'].value)
'renergy(t)', if computation_needed:
'mean_trS2(t)']: for k in ['t', 'vel_max(t)', 'renergy(t)',
if k in pp_file.keys(): 'energy(t)', 'enstrophy(t)',
self.statistics[k] = pp_file[k].value 'energy(k)', 'enstrophy(k)',
self.compute_time_averages() 'energy(t, k)',
'enstrophy(t, k)',
'R_ij(t)',
'ii0', 'ii1', 'iter0', 'iter1']:
if k in pp_file.keys():
del pp_file[k]
if computation_needed:
pp_file['iter0'] = iter0
pp_file['iter1'] = iter1
pp_file['ii0'] = ii0
pp_file['ii1'] = ii1
pp_file['t'] = (self.parameters['dt']*
self.parameters['niter_stat']*
(np.arange(ii0, ii1+1).astype(np.float)))
phi_ij = data_file['statistics/spectra/velocity_velocity'][ii0:ii1+1]
pp_file['R_ij(t)'] = np.sum(phi_ij, axis = 1)
energy_tk = (
phi_ij[:, :, 0, 0] +
phi_ij[:, :, 1, 1] +
phi_ij[:, :, 2, 2])/2
pp_file['energy(t)'] = np.sum(energy_tk, axis = 1)
pp_file['energy(k)'] = np.mean(energy_tk, axis = 0)*(4*np.pi*self.statistics['kshell']**2) / (self.statistics['dk']*self.statistics['nshell'])
enstrophy_tk = (
data_file['statistics/spectra/vorticity_vorticity'][ii0:ii1+1, :, 0, 0] +
data_file['statistics/spectra/vorticity_vorticity'][ii0:ii1+1, :, 1, 1] +
data_file['statistics/spectra/vorticity_vorticity'][ii0:ii1+1, :, 2, 2])/2
pp_file['enstrophy(t)'] = np.sum(enstrophy_tk, axis = 1)
pp_file['enstrophy(k)'] = np.mean(enstrophy_tk, axis = 0)*(4*np.pi*self.statistics['kshell']**2) / (self.statistics['dk']*self.statistics['nshell'])
pp_file['vel_max(t)'] = data_file['statistics/moments/velocity'][ii0:ii1+1, 9, 3]
pp_file['renergy(t)'] = data_file['statistics/moments/velocity'][ii0:ii1+1, 2, 3]/2
if 'trS2_Q_R' in data_file['statistics/moments'].keys():
pp_file['mean_trS2(t)'] = data_file['statistics/moments/trS2_Q_R'][:, 1, 0]
for k in ['t',
'energy(t)',
'energy(k)',
'enstrophy(t)',
'enstrophy(k)',
'R_ij(t)',
'vel_max(t)',
'renergy(t)',
'mean_trS2(t)']:
if k in pp_file.keys():
self.statistics[k] = pp_file[k].value
# sanity check --- Parseval theorem check
assert(np.max(np.abs(
self.statistics['renergy(t)'] -
self.statistics['energy(t)']) / self.statistics['energy(t)']) < 1e-5)
self.compute_time_averages()
return None return None
def compute_Reynolds_stress_invariants( def compute_Reynolds_stress_invariants(
self): self):
......
...@@ -224,8 +224,10 @@ class _base(object): ...@@ -224,8 +224,10 @@ class _base(object):
ofile[group + '/' + k][...] = parameters[k] ofile[group + '/' + k][...] = parameters[k]
ofile.close() ofile.close()
return None return None
def read_parameters(self): def read_parameters(self, fname = None):
with h5py.File(os.path.join(self.work_dir, self.simname + '.h5'), 'r') as data_file: if type(fname) == type(None):
fname = os.path.join(self.work_dir, self.simname + '.h5')
with h5py.File(fname, 'r') as data_file:
for k in data_file['parameters'].keys(): for k in data_file['parameters'].keys():
if k in self.parameters.keys(): if k in self.parameters.keys():
if type(self.parameters[k]) in [int, str, float]: if type(self.parameters[k]) in [int, str, float]:
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment