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

add simple numeric tests

parent d82f680e
No related branches found
No related tags found
No related merge requests found
......@@ -29,6 +29,9 @@ import sys
import subprocess
import pickle
import numpy as np
import matplotlib.pyplot as plt
import bfps
from bfps import fluid_resize
......@@ -212,79 +215,91 @@ def acceleration_test(c, m = 3, species = 0):
plt.close(fig)
return pid
def compare_stats(opt, c0, c1):
# plot energy and enstrophy
fig = plt.figure(figsize = (12, 12))
a = fig.add_subplot(221)
c0.set_plt_style({'label' : '1',
'dashes' : (None, None),
'color' : (1, 0, 0)})
c1.set_plt_style({'label' : '2',
'dashes' : (2, 2),
'color' : (0, 0, 1)})
for c in [c0, c1]:
a.plot(c.statistics['t'],
c.statistics['energy(t)'],
label = c.style['label'],
dashes = c.style['dashes'],
color = c.style['color'])
a.set_title('energy')
a.legend(loc = 'best')
a = fig.add_subplot(222)
for c in [c0, c1]:
a.plot(c.statistics['t'],
c.statistics['enstrophy(t)'],
dashes = c.style['dashes'],
color = c.style['color'])
a.set_title('enstrophy')
a = fig.add_subplot(223)
for c in [c0, c1]:
a.plot(c.statistics['t'],
c.statistics['kM']*c.statistics['etaK(t)'],
dashes = c.style['dashes'],
color = c.style['color'])
a.set_title('$k_M \\eta_K$')
a = fig.add_subplot(224)
for c in [c0, c1]:
a.plot(c.statistics['t'],
c.statistics['vel_max(t)'] * (c.parameters['dt'] * c.parameters['dkx'] /
(2*np.pi / c.parameters['nx'])),
dashes = c.style['dashes'],
color = c.style['color'])
a.set_title('$\\frac{\\Delta t \\| u \\|_\infty}{\\Delta x}$')
fig.savefig('plain_stats_{0}.pdf'.format(opt.precision), format = 'pdf')
def compare_stats(
opt,
c0, c1,
plots_on = False):
for key in ['energy', 'enstrophy', 'vel_max']:
print('maximum {0} difference is {1}'.format(
key,
np.max(np.abs(c0.statistics[key + '(t)'] - c0.statistics[key + '(t)']))))
for i in range(c0.particle_species):
print('maximum traj difference species {0} is {1}'.format(
i,
np.max(np.abs(c0.trajectories[i] - c1.trajectories[i]))))
if plots_on:
# plot energy and enstrophy
fig = plt.figure(figsize = (12, 12))
a = fig.add_subplot(221)
c0.set_plt_style({'label' : '1',
'dashes' : (None, None),
'color' : (1, 0, 0)})
c1.set_plt_style({'label' : '2',
'dashes' : (2, 2),
'color' : (0, 0, 1)})
for c in [c0, c1]:
a.plot(c.statistics['t'],
c.statistics['energy(t)'],
label = c.style['label'],
dashes = c.style['dashes'],
color = c.style['color'])
a.set_title('energy')
a.legend(loc = 'best')
a = fig.add_subplot(222)
for c in [c0, c1]:
a.plot(c.statistics['t'],
c.statistics['enstrophy(t)'],
dashes = c.style['dashes'],
color = c.style['color'])
a.set_title('enstrophy')
a = fig.add_subplot(223)
for c in [c0, c1]:
a.plot(c.statistics['t'],
c.statistics['kM']*c.statistics['etaK(t)'],
dashes = c.style['dashes'],
color = c.style['color'])
a.set_title('$k_M \\eta_K$')
a = fig.add_subplot(224)
for c in [c0, c1]:
a.plot(c.statistics['t'],
c.statistics['vel_max(t)'] * (c.parameters['dt'] * c.parameters['dkx'] /
(2*np.pi / c.parameters['nx'])),
dashes = c.style['dashes'],
color = c.style['color'])
a.set_title('$\\frac{\\Delta t \\| u \\|_\infty}{\\Delta x}$')
fig.savefig('plain_stats_{0}.pdf'.format(opt.precision), format = 'pdf')
fig = plt.figure(figsize = (12, 12))
a = fig.add_subplot(221)
a.plot(c0.statistics['t'],
c0.statistics['energy(t)'] - c1.statistics['energy(t)'])
a.set_title('energy')
a = fig.add_subplot(222)
a.plot(c0.statistics['t'],
c0.statistics['enstrophy(t)'] - c1.statistics['enstrophy(t)'])
a.set_title('enstrophy')
a = fig.add_subplot(223)
a.plot(c0.statistics['t'],
c0.statistics['kM']*c0.statistics['etaK(t)'] - c1.statistics['kM']*c1.statistics['etaK(t)'])
a.set_title('$k_M \\eta_K$')
a = fig.add_subplot(224)
data0 = c0.statistics['vel_max(t)'] * (c0.parameters['dt'] * c0.parameters['dkx'] /
(2*np.pi / c0.parameters['nx']))
data1 = c1.statistics['vel_max(t)'] * (c1.parameters['dt'] * c1.parameters['dkx'] /
(2*np.pi / c1.parameters['nx']))
a.plot(c0.statistics['t'],
data0 - data1)
a.set_title('$\\frac{\\Delta t \\| u \\|_\infty}{\\Delta x}$')
fig.savefig('plain_stat_diffs_{0}.pdf'.format(opt.precision), format = 'pdf')
fig = plt.figure(figsize = (12, 12))
a = fig.add_subplot(221)
a.plot(c0.statistics['t'],
c0.statistics['energy(t)'] - c1.statistics['energy(t)'])
a.set_title('energy')
a = fig.add_subplot(222)
a.plot(c0.statistics['t'],
c0.statistics['enstrophy(t)'] - c1.statistics['enstrophy(t)'])
a.set_title('enstrophy')
a = fig.add_subplot(223)
a.plot(c0.statistics['t'],
c0.statistics['kM']*c0.statistics['etaK(t)'] - c1.statistics['kM']*c1.statistics['etaK(t)'])
a.set_title('$k_M \\eta_K$')
a = fig.add_subplot(224)
data0 = c0.statistics['vel_max(t)'] * (c0.parameters['dt'] * c0.parameters['dkx'] /
(2*np.pi / c0.parameters['nx']))
data1 = c1.statistics['vel_max(t)'] * (c1.parameters['dt'] * c1.parameters['dkx'] /
(2*np.pi / c1.parameters['nx']))
a.plot(c0.statistics['t'],
data0 - data1)
a.set_title('$\\frac{\\Delta t \\| u \\|_\infty}{\\Delta x}$')
fig.savefig('plain_stat_diffs_{0}.pdf'.format(opt.precision), format = 'pdf')
# plot trajectory differences
for i in range(c0.particle_species):
fig = plt.figure(figsize=(12, 4))
for j in range(3):
a = fig.add_subplot(131 + j)
for t in range(c0.parameters['nparticles']):
a.plot(c0.trajectories[i][:, t, j] - c1.trajectories[i][:, t, j])
fig.savefig('traj_s{0}_{1}.pdf'.format(i, opt.precision), format = 'pdf')
# plot trajectory differences
for i in range(c0.particle_species):
fig = plt.figure(figsize=(12, 4))
for j in range(3):
a = fig.add_subplot(131 + j)
for t in range(c0.parameters['nparticles']):
a.plot(c0.trajectories[i][:, t, j] - c1.trajectories[i][:, t, j])
fig.savefig('traj_s{0}_{1}.pdf'.format(i, opt.precision), format = 'pdf')
return None
if __name__ == '__main__':
......
......@@ -53,6 +53,7 @@ if __name__ == '__main__':
opt = parser.parse_args(
['-n', '32',
'--run',
'--initialize',
'--ncpu', '2',
'--nparticles', '1000',
'--niter_todo', '32',
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment