test_Parseval.py

#! /usr/bin/env python

import os
import numpy as np
import h5py
import sys

import bfps
from bfps import DNS

import matplotlib.pyplot as plt

def main():
    niterations = 10
    c = DNS()
    c.launch(
            ['NSVE',
             '--nx', '32',
             '--ny', '32',
             '--nz', '32',
             '--forcing_type', 'linear',
             '--np', '4',
             '--ntpp', '1',
             '--niter_todo', '{0}'.format(niterations),
             '--niter_out', '{0}'.format(niterations),
             '--niter_stat', '1',
             '--wd', './'] +
             sys.argv[1:])
    c.compute_statistics()
    print((c.statistics['energy(t)'] - c.statistics['renergy(t)']) / c.statistics['renergy(t)'])
    energyk = c.statistics['energy(k)']
    nshell = c.get_data_file()['kspace/nshell'].value
    renergy = np.mean(c.statistics['renergy(t)'])
    print(renergy, np.trapz(energyk[:-2], c.statistics['kshell'][:-2]))

    f = plt.figure()
    a = f.add_subplot(111)
    a.plot(c.statistics['kshell'], energyk)
    a.plot(c.statistics['kshell'], (energyk / nshell)*(4*np.pi*c.statistics['kshell']**2))
    a.set_yscale('log')
    a.set_xscale('log')
    f.savefig('spectrum.pdf')
    return None

if __name__ == '__main__':
    main()