bug fix and doc update

docs/src/index.md

@@ -16,6 +16,7 @@ SeparableFunctions
 ## Generic Ways to Create Seperable and Radial Functions 
 ```@docs
 calculate_broadcasted
+calculate_separables
 calculate_separables_nokw
 separable_view
 separable_create

src/SeparableFunctions.jl

@@ -33,7 +33,7 @@ using ChainRulesCore # for adjoint definition
 # using ZygoteRules
 # using Zygote # to use rrule_via_ad
 
-export calculate_separables_nokw, separable_view, separable_create
+export calculate_separables_nokw, calculate_separables, separable_view, separable_create
 export calculate_broadcasted
 export copy_corners!
 export propagator_col, propagator_col!, phase_kz_col, phase_kz_col!

src/general.jl

@@ -92,6 +92,45 @@ function ChainRulesCore.rrule(config::RuleConfig{>:HasReverseMode}, ::typeof(in_
 end
 
 
+"""
+    calculate_separables([::Type{AT},] fct, sz::NTuple{N, Int}, args...; 
+                                    all_axes = (similar_arr_type(AT, eltype(AT), Val(1)))(undef, sum(sz)),
+                                    defaults=NamedTuple(), 
+                                    offset = sz.÷2 .+1,
+                                    scale = one(real(eltype(AT))),
+                                    kwargs...) where {AT, N}
+
+creates a list of one-dimensional vectors, which can be combined to yield a separable array. In a way this can be seen as a half-way Lazy operation.
+The (potentially heavy) work of calculating the one-dimensional functions is done now but the memory-heavy calculation of the array is done later.
+This function is used in `separable_view` and `separable_create`.
+For automatic differentiation support (e.g. for opinizing) you should use the `calculate_separables_nokw` function, since it supports AD, whereas the keyword-based version does not.
+
+#Arguments
++ `AT`:     optional type signfying the array result type. You can for example use `CuArray{Float32}` using `CUDA` to create the views on the GPU.
++ `fct`:    the function to calculate for each axis index (no need for broadcasting!) of this iterable of seperable axes. Note that the first arguments of `fct` have to be the index of this coordinate and the size of this axis. Any further `args` and `nargs` can follow. Often the second argument is not used but it still needs to be present.
++ `sz`:     the size of the result array (when appying the one-D axes)
++ `offset`: specifying the center (zero-position) of the result array in one-based coordinates. The default corresponds to the Fourier-center.
++ `scale`:  multiplies the index before passing it to `fct`
++ `args`:   further arguments which are passed over to the function `fct`.
++ `all_axes`: if provided, this memory is used instead of allocating a new one. This can be useful if you want to use the same memory for multiple calculations.
+
+#Example:
+```julia
+julia> fct = (r, sz, sigma)-> exp(-r^2/(2*sigma^2))
+julia> gauss_sep = SeparableFunctions.calculate_separables_nokw(Array{Float32}, fct, (6,5), (0.5,1.0), 1.0, 1.0)
+2-element Vector{Array{Float32}}:
+ [4.4963495f-9; 0.00014774836; … ; 0.1978987; 0.0007318024;;]
+ [0.088921614 0.48675224 … 0.726149 0.1978987]
+ julia> my_gaussian = .*(gauss_sep...) # this is how to broadcast it
+ 6×5 Matrix{Float32}:
+ 3.99823f-10  2.18861f-9   4.40732f-9   3.26502f-9   8.89822f-10
+ 1.3138f-5    7.19168f-5   0.000144823  0.000107287  2.92392f-5
+ 0.00790705   0.0432828    0.0871608    0.0645703    0.0175975
+ 0.0871608    0.477114     0.960789     0.71177      0.19398
+ 0.0175975    0.0963276    0.19398      0.143704     0.0391639
+ 6.50731f-5   0.000356206  0.000717312  0.000531398  0.000144823
+```
+"""
 function calculate_separables(::Type{AT}, fct, sz::NTuple{N, Int}, 
     args...; 
     all_axes = (similar_arr_type(AT, eltype(AT), Val(1)))(undef, sum(sz)),