Best approach to smooth array? None of my approaches are very useful

A “low-tech” but relatively fast and simple approach may be to consider explicit diffusion as smoothing mechanism. Let’s say you want to smooth a field A, then following would do the job in 2D (and can be easily extended to n-D, GPU parallel processing, …):

using Plots
@views function do_smoothing()
    nsmooth = 50
    A       = rand(100,100)
    CFL     = 0.5/2.1/length(size(A))
    for ismooth = 1:nsmooth
        A[2:end-1,2:end-1] .= A[2:end-1,2:end-1] + CFL*(diff(diff(A[:,2:end-1],dims=1),dims=1) + diff(diff(A[2:end-1,:],dims=2),dims=2))
        A[1,:] = A[2,:]; A[end,:] = A[end-1,:]; A[:,1] = A[:,2]; A[:,end] = A[:,end-1]
    end
    display(heatmap(A'))
end

@time do_smoothing()

Then you could as well decouple the diffusion in different dimensions and apply various number of smoothing steps individually for different dimensions.


EDIT: To put this in context with the moving mean filter, it would e.g. take about 10 diffusion smoothing steps nsmooth=10 to achieve a ~10-cell window smoothing:


modifying the above code as such

using Plots
@views function do_smoothing()
    nsmooth = 10
    A0      = 0.0*rand(100,100)
    A0[40:60,40:60] .= 1.0
    A       = copy(A0)
    CFL     = 0.5/2.1/length(size(A))
    for ismooth = 1:nsmooth
        A[2:end-1,2:end-1] .= A[2:end-1,2:end-1] + CFL*(diff(diff(A[:,2:end-1],dims=1),dims=1) + diff(diff(A[2:end-1,:],dims=2),dims=2))
        A[1,:] = A[2,:]; A[end,:] = A[end-1,:]; A[:,1] = A[:,2]; A[:,end] = A[:,end-1]
    end
    p=plot(A0[:,50],label="initial", linewidth=3); plot!(A[:,50], linewidth=3, markersize=5, markershape=:circle, label="smooth", framestyle=:box)
    display(p)
end

@time do_smoothing()