refactor

DualTVDD/src/fd.jl

 module FD
 
+# local modules
+"Finite Differences"
+module FDiff
+
+using OffsetArrays: OffsetVector
+
+const forward = OffsetVector([-1, 1], -1)
+const backward = OffsetVector([-1, 1], -2)
+const central = OffsetVector([-1, 0, 1], -2)
+
+end # module FDiff
+
+
+# imports
 using StaticArrays: SVector
 using OffsetArrays
 import Colors: Gray
 import Plots: plot
+# local imports
+using ..Filtering
+using .FDiff
 
-
+# export
 export Grid, MetaGrid, GridView, GridFunction
 export grid
 
 
+# abstract types
 """
     AbstractGrid{d}
 
@@ -22,21 +40,6 @@ Base.IndexStyle(::Type{<:AbstractGrid}) = IndexCartesian()
 Base.parent(grid::AbstractGrid) = grid
 
 
-"Finite Differences"
-module FDiff
-
-using OffsetArrays
-
-const forward = OffsetVector([-1, 1], -1)
-const backward = OffsetVector([-1, 1], -2)
-const central = OffsetVector([-1, 0, 1], -2)
-
-end
-
-
-using ..Filtering
-using .FDiff
-
 """
     _fdkern(fd, ndims, dir)
 
@@ -65,6 +68,8 @@ struct Grid{d} <: AbstractGrid{d}
     domain::NTuple{d, FDStepRange}
 end
 Grid(domain::FDStepRange...) = Grid(domain)
+
+# TODO: somehow display(grid) is incorrect for d = 1 ?
 Base.show(io::IO, grid::Grid) =
     print("Grid from $(first.(grid.domain)) to $(last.(grid.domain)) ",
           "with spacing $(step.(grid.domain))")
@@ -97,7 +102,7 @@ Base.getindex(grid::GridView, idx...) = getindex(grid.grid, (first.(grid.indices
 Base.parent(grid::GridView) = grid.grid
 
 
-# MetaGrid
+# MetaGrid (partitioned Grid)
 
 
 struct MetaGrid{d} <: AbstractGrid{d}
@@ -106,14 +111,16 @@ struct MetaGrid{d} <: AbstractGrid{d}
     overlap::Int
 end
 
+MetaGrid(domain::NTuple{d, FDStepRange}, pnum, overlap) where d = MetaGrid(Grid(domain), pnum, overlap)
+
 """
     MetaGrid(domain::NTuple{d, FDStepRange}, pnum::NTuple{d, Int}, overlap::Int)
+    MetaGrid(grid::Grid{d}, pnum::NTuple{d, Int}, overlap::Int)
 
 Create a grid on `domain`, partitioned along dimensions according to `pnum`
 respecting `overlap`.
 """
-function MetaGrid(domain::NTuple{d, FDStepRange}, pnum::NTuple{d, Int}, overlap::Int) where d
-    grid = Grid(domain)
+function MetaGrid(grid::Grid{d}, pnum::NTuple{d, Int}, overlap::Int) where d
     tsize = size(grid) .+ (pnum .- 1) .* overlap
 
     psize = tsize .÷ pnum
@@ -143,7 +150,8 @@ Base.view(grid::MetaGrid, idx...) = GridView(grid, indices(grid, idx...))
 parts(grid::MetaGrid) = CartesianIndices(grid.indices)
 
 
-## GridFunction
+# GridFunction
+
 
 abstract type AbstractGridFunction{G,T,d,n} <: AbstractArray{T,1} end
 
@@ -376,8 +384,14 @@ function my_norm(f::GridFunction)
     GridFunction(f.grid, sqrt.(sum(f.data.^2, dims=d+1)), Val(rangedim(f)))
 end
 
+"""
+    f ⊙ g
+
+Pointwise multiplication of two gridfunctions.
+"""
+⊙(f, g) = pwmul(f, g)
 function pwmul(f::GridFunction{G,T,d,1}, g::GridFunction{G,T,d}) where {G,T,d}
-    # TODO: maybe comprehension?
+    # TODO: maybe comprehension? Or make broadcast work
     data = similar(g.data)
     for I in CartesianIndices(f.grid)
         data[Tuple(I)..., :] = f.data[Tuple(I)..., 1] * g.data[Tuple(I)..., :]
@@ -385,9 +399,6 @@ function pwmul(f::GridFunction{G,T,d,1}, g::GridFunction{G,T,d}) where {G,T,d}
     GridFunction(f.grid, data, Val(rangedim(g)))
 end
 
-
-
-
 """
     neuman_kernel(kernel, bidx)