FourierFlows/GeophysicalFlows.jl

construct a random field with prescribed spectral slope

Opened this issue · 2 comments

I would like to construct a random initial field with a prescribed spectral slope x. If Energy = \int E(k) k dk then E(k) ~ k^{-x}.
As suggested by Navid, maybe we can iterate on peakedisotropicspectrum to implement this.

I think this should do

function slopedisotropicspectrum(grid::TwoDGrid{T, A}, slope::Real; mask=ones(size(grid.Krsq)), allones=false) where {T, A}
 if grid.Lx !== grid.Ly
   error("the domain is not square")
 else
   modk = sqrt.(grid.Krsq)
   modψ = A(zeros(T, (grid.nk, grid.nl)))
   modψ = @. modk^((slope-1)/2)
   #modψ[1, 1] = 0.0
   CUDA.@allowscalar modψ[1, 1] = 0.0

   phases = randn(Complex{T}, size(grid.Krsq))
   phases_real, phases_imag = real.(phases), imag.(phases)
   phases = A(phases_real) + im * A(phases_imag)
   ψh = @. phases * modψ
   if allones; ψh = modψ; end
   ψh = ψh .* A(mask)
   q = A(irfft(ψh, grid.nx))
 end

 return q
end

So what you want is a function that returns a field q for which if energy is $E = \int q(x, y)^2 \mathrm{d}x \mathrm{d}y = \int \mathcal{E}(k) k \mathrm{d}k$ you want $\mathcal{E}(k) \propto k^{-\chi}$. Am I correct?