Private types and functions

Documentation for FourierFlows.jl's internal interface.

FourierFlows

Problem

FourierFlows.ClockType
mutable struct Clock{T<:AbstractFloat}

Represents the clock of a problem.

  • dt::AbstractFloat

    the time-step

  • t::AbstractFloat

    the time

  • step::Int64

    the step number

source
FourierFlows.EquationType
struct Equation{T, TL, G<:AbstractFloat}

The equation to be solved ∂u/∂t = L*u + N(u). Array L includes the coefficients of the linear term L*u and calcN! is a function which computes the nonlinear term N(u). The struct also includes the problem's grid and the float type of the state vector (and consequently of N(u)).

  • L::Any

    array with coefficient for the linear part of the equation

  • calcN!::Function

    function that computes the nonlinear part of the equation

  • grid::AbstractGrid{G, A, Alias} where {G<:AbstractFloat, A, Alias}

    the grid

  • dims::Tuple

    the dimensions of L

  • T::Any

    the float type for the state vector

source
FourierFlows.EquationMethod
Equation(L, calcN!, grid; dims=supersize(L), T=nothing)

The equation constructor from the array L of the coefficients of the linear term, the function calcN! that computes the nonlinear term and the grid for the problem.

source
FourierFlows.ProblemType
struct Problem{T, A<:AbstractArray, Tg<:AbstractFloat, TL}

A problem that represents a partial differential equation.

  • sol::AbstractArray

    the state vector

  • clock::FourierFlows.Clock{Tg} where Tg<:AbstractFloat

    the problem's

  • eqn::FourierFlows.Equation{T, TL, Tg} where {T, Tg<:AbstractFloat, TL}

    the equation

  • grid::AbstractGrid{Tg, A, Alias} where {Tg<:AbstractFloat, A, Alias}

    the grid

  • vars::AbstractVars

    the variables

  • params::AbstractParams

    the parameters

  • timestepper::AbstractTimeStepper{A} where A<:AbstractArray

    the timestepper

source
FourierFlows.ProblemMethod
Problem(eqn::Equation, stepper, dt, grid::AbstractGrid{T}, 
        vars=EmptyVars, params=EmptyParams, dev::Device=CPU(); stepperkwargs...) where T

Construct a Problem for equation eqn using the timestepper with timestep dt, on grid and on device. Optionally, use the keyword arguments to provide variables with vars and parameters with params. The stepperkwargs are passed to the time-stepper constructor.

source

Domain

FourierFlows.getaliasedwavenumbersMethod
getaliasedwavenumbers(nk, nkr, aliased_fraction)

Returns the top aliased_fraction highest wavenumbers, both for and real FFTs, kalias and kralias respectively. For example, aliased_fraction=1/3 should return the indices of the top-most 1/6-th (in absolute value) for both positive and negative wavenumbers (i.e., 1/3 total) that should be set to zero after performing an FFT.

source
FourierFlows.makefilterMethod
makefilter(K; order=4, innerK=0.65, outerK=1)

Returns a filter acting on the non-dimensional wavenumber K that decays exponentially for K>innerK, thus removing high-wavenumber content from a spectrum it is multiplied with. The decay rate is determined by order and outerK determines the outer wavenumber at which the filter is smaller than Float64 machine precision.

source

Utilities

FourierFlows.on_gridMethod
on_grid(func, grid)

Returns an array, of the ArrayType of the device grid lives on, that contains the values of function func evaluated on the grid.

source
FourierFlows.parsevalsumMethod
parsevalsum(uh, grid)

Returns real(Σ uh) on the grid, i.e.

\[ℜ [ \sum_{𝐤} û_{𝐤} L_x L_y ] \,,\]

where $û_{𝐤} =$ uh $/($ grid.nx $e^{- i 𝐤 ⋅ 𝐱₀})$, with $𝐱₀$ the vector with components the left-most position in each direction.

source
FourierFlows.parsevalsum2Method
parsevalsum2(uh, grid)

Returns Σ |uh|² on the grid, which is equal to the domain integral of u. More specifically, it returns

\[\sum_{𝐤} |û_{𝐤}|² L_x L_y = \int u(𝐱)² \, 𝖽x 𝖽y \,,\]

where $û_{𝐤} =$ uh $/($ grid.nx $e^{- i 𝐤 ⋅ 𝐱₀})$, with $𝐱₀$ the vector with components the left-most position in each direction.

source
FourierFlows.radialspectrumMethod
radialspectrum(fh, grid; n=nothing, m=nothing, refinement=2)

Returns the radial spectrum of fh. fh lives on Cartesian wavenumber grid $(k, l)$. To compute the radial spectrum, we first interpolate $f̂(k, l)$ onto a radial wavenumber grid $(ρ, θ)$, where $ρ² = k²+l²$ and $θ = \arctan(l/k)$. Note here that $f̂ =$ fh $/($ grid.nx $e^{- i 𝐤 ⋅ 𝐱₀})$, with $𝐱₀$ the vector with components the left-most position in each direction. After interpolation, we integrate $f̂$over angles $θ$ to get ,

\[ f̂_ρ = \int f̂(ρ, θ) ρ 𝖽ρ 𝖽θ \, .\]

The resolution (n, m) for the polar wavenumber grid is n = refinement * maximum(nk, nl), m = refinement * maximum(nk, nl), where refinement = 2 by default. If fh is in conjugate symmetric form then only the upper-half plane in $θ$ is represented on the polar grid.

source

Diagnostics

FourierFlows.extend!Method
extend!(diag::Diagnostic{T,N})

Double the extend of the data, time, and steps vectors of the diagnostic diag.

source

Output

FourierFlows.uniquepathMethod
uniquepath(path)

Return path with a number appended if isfile(path) == true. The number is incremented until path does not exist.

source

Timesteppers

FourierFlows.getetdcoeffsMethod
getetdcoeffs(dt, L; ncirc=32, rcirc=1)

Calculate ETDRK4 coefficients associated with the (diagonal) linear coefficient L by integrating over a small circle in complex space.

source

Diffusion Testbed Module

FourierFlows.Diffusion.VarsType
struct Vars{Aphys, Atrans} <: AbstractVars

The variables for diffusion problem:

  • c

    tracer concentration

  • cx

    tracer concentration derivative ∂ₓc

  • ch

    Fourier transform of tracer concentration

  • cxh

    Fourier transform of tracer concentration derivative ∂ₓc

source