diff --git a/src/oblique_stencil.jl b/src/oblique_stencil.jl index 076b41a..741c8b7 100644 --- a/src/oblique_stencil.jl +++ b/src/oblique_stencil.jl @@ -10,12 +10,12 @@ function dipole(x, y, z, r_treshold, direction=(0, 0, 1)) return kappa end -function stencil(; st=27, res=(1.0, 1.0, 1.0), res0=1.0, singularity_cutout=4, direction=(0, 0, 1), gridsize=(64, 64, 64)) +function stencil(; st=27, res=(1.0, 1.0, 1.0), singularity_cutout=4, direction=(0, 0, 1), gridsize=(64, 64, 64)) middle = floor.(Int, gridsize ./ 2) .+ 1 - coord = [((1:gridsize[i]) .- middle[i]) * res0 for i in 1:3] + coord = [((1:gridsize[i]) .- middle[i]) for i in 1:3] - d = [dipole(x, y, z, singularity_cutout * res0, direction) for x in coord[1], y in coord[2], z in coord[3]] + d = [dipole(x, y, z, singularity_cutout, direction) for x in coord[1], y in coord[2], z in coord[3]] d_mask = isfinite.(d) # stencil mask @@ -63,14 +63,13 @@ function stencil(; st=27, res=(1.0, 1.0, 1.0), res0=1.0, singularity_cutout=4, d y = F.Vt * d[d_mask] x = F.U * (y .* s_inv) - x_corr = x * res0^3 stencil = zeros(Float32, 3, 3, 3) ind = 1 for i in eachindex(stencil) if mask[i] - stencil[i] = 2x_corr[ind] + stencil[i] = 2x[ind] ind += 1 end end