Skip to content
Snippets Groups Projects
ornstein_uhlenbeck_test.cpp 2.82 KiB
#include <string>
#include <cmath>
#include "ornstein_uhlenbeck_test.hpp"
#include "scope_timer.hpp"


template <typename rnumber>
int ornstein_uhlenbeck_test<rnumber>::initialize(void)
{
    TIMEZONE("ornstein_uhlenbeck_test::initialize");
    this->read_parameters();
    this->ou = new ornstein_uhlenbeck_process<rnumber,FFTW>(
        "test",
        this->nx,this->ny,this->nz,
        this->ou_kmin, this->ou_kmax, this->ou_energy_amplitude,
        this->dkx,this->dky,this->dkz,
        FFTW_ESTIMATE);
    if (this->myrank == 0)
    {
        hid_t stat_file = H5Fopen(
                (this->simname + std::string(".h5")).c_str(),
                H5F_ACC_RDWR,
                H5P_DEFAULT);
        this->ou->kk->store(stat_file);
        H5Fclose(stat_file);
    }
    return EXIT_SUCCESS;
}

template <typename rnumber>
int ornstein_uhlenbeck_test<rnumber>::finalize(void)
{
    TIMEZONE("ornstein_uhlenbeck_test::finalize");
    delete this->ou;
    return EXIT_SUCCESS;
}

template <typename rnumber>
int ornstein_uhlenbeck_test<rnumber>::read_parameters()
{
    TIMEZONE("ornstein_uhlenbeck_test::read_parameters");
    this->test::read_parameters();
    // hid_t parameter_file = H5Fopen(
    //         (this->simname + std::string(".h5")).c_str(),
    //         H5F_ACC_RDONLY,
    //         H5P_DEFAULT);
    // this->ou_kmin = hdf5_tools::read_value<double>(parameter_file, "/parameters/ou_kmin");
    // this->ou_kmax = hdf5_tools::read_value<double>(parameter_file, "/parameters/ou_kmax");
    // this->ou_energy_amplitude = hdf5_tools::read_value<double>(parameter_file, "/parameters/ou_energy_amplitude");
    // H5Fclose(parameter_file);
    //TODO
    this->ou_kmin = 3;
    this->ou_kmax = 14;
    this->ou_energy_amplitude = 1;


    return EXIT_SUCCESS;
}

template <typename rnumber>
int ornstein_uhlenbeck_test<rnumber>::do_work(void)
{
    TIMEZONE("ornstein_uhlenbeck_test::do_work");
    std::string filename = this->simname + std::string("_fields.h5");
    for (int step = 0; step < 10000; step++)
    {
        this->ou->step(1e-3);
    }
    hid_t out_file;
    if (this->myrank == 0) {
      out_file = H5Fopen("test_ou_output.h5", H5F_ACC_RDWR, H5P_DEFAULT);
    }

    this->ou->kk->template cospectrum<rnumber,THREE>(
      this->ou->ou_field->get_cdata(),
      this->ou->ou_field->get_cdata(),
      out_file,"ou_spectra",0);

    this->ou->ou_field->symmetrize();
    this->ou->ou_field->ift();
    std::vector<double> me;
    me.resize(3);
    std::fill(me.begin(),me.end(),30);
    // this->ou->ou_field->template compute_rspace_stats<rnumber, FFTW, THREE>(out_file, "ou", 0, me);
    this->ou->ou_field->compute_rspace_stats(out_file, "ou", 0, me);

    this->ou->ou_field->io("ou_field.h5", "ou_field", 0, false);

    return EXIT_SUCCESS;
}

template class ornstein_uhlenbeck_test<float>;
template class ornstein_uhlenbeck_test<double>;