Program Listing for File SimpleGridTurbulence.h

Return to documentation for file (include/crpropa/magneticField/turbulentField/SimpleGridTurbulence.h)

#ifndef CRPROPA_SIMPLEGRIDTURBULENCE_H
#define CRPROPA_SIMPLEGRIDTURBULENCE_H

#ifdef CRPROPA_HAVE_FFTW3F

#include "crpropa/magneticField/turbulentField/GridTurbulence.h"

#include "kiss/logger.h"
#include "kiss/string.h"

namespace crpropa {
class SimpleTurbulenceSpectrum : public TurbulenceSpectrum {
        const int constScaleBendover = 1000;    // define the bandover scale as 1000 * lMax to ensure k * lBendover >> 1. The bendover scale is necessary for the implementation of PlaneWaveTurbulence.
  public:
        SimpleTurbulenceSpectrum(double Brms, double lMin, double lMax,
                                 double sIndex = 5. / 3)
            : TurbulenceSpectrum(Brms, lMin, lMax, constScaleBendover * lMax, sIndex, 0) {}
        ~SimpleTurbulenceSpectrum() {}

        double energySpectrum(double k) const {
                return std::pow(k, -getSindex() - 2);
        }

        double getCorrelationLength() const {
                return turbulentCorrelationLength(getLmin(), getLmax(),
                                                  getSindex());
        }
        static double turbulentCorrelationLength(double lMin, double lMax,
                                                 double s) {
                double r = lMin / lMax;
                return lMax / 2 * (s - 1) / s * (1 - pow(r, s)) / (1 - pow(r, s - 1));
        }
};

class SimpleGridTurbulence : public GridTurbulence {
  public:
        SimpleGridTurbulence(const SimpleTurbulenceSpectrum &spectrum,
                             const GridProperties &gridProp, unsigned int seed = 0);

        static void initTurbulence(ref_ptr<Grid3f> grid, double Brms, double lMin,
                                   double lMax, double alpha, int seed);
};

// Compatibility with old functions from GridTurbulence:

inline double turbulentCorrelationLength(double lMin, double lMax,
                                         double alpha = -11 / 3.) {
        KISS_LOG_WARNING
            << "turbulentCorrelationLength is deprecated and will be "
               "removed in the future. Replace it with a more appropriate "
               "turbulent field model and call getCorrelationLength().";
        return SimpleTurbulenceSpectrum::turbulentCorrelationLength(lMin, lMax,
                                                                    -alpha - 2);
}

inline void initTurbulence(ref_ptr<Grid3f> grid, double Brms, double lMin,
                           double lMax, double alpha = -11 / 3., int seed = 0) {
        KISS_LOG_WARNING
            << "initTurbulence is deprecated and will be removed in the future. "
               "Replace it with a more appropriate turbulent field model instance.";
        SimpleGridTurbulence::initTurbulence(grid, Brms, lMin, lMax, alpha, seed);
}

} // namespace crpropa

#endif // CRPROPA_HAVE_FFTW3F

#endif // CRPROPA_SIMPLEGRIDTURBULENCE_H