diff src/Grids/grid.jl @ 1531:9da4ab4fb85e bugfix/sbp_operators/stencil_return_type

Merge default
author Jonatan Werpers <jonatan@werpers.com>
date Thu, 11 Apr 2024 22:50:42 +0200
parents d641798539c2
children 9113f437431d 611ae2308aa1
line wrap: on
line diff
--- a/src/Grids/grid.jl	Sun Nov 26 22:30:44 2023 +0100
+++ b/src/Grids/grid.jl	Thu Apr 11 22:50:42 2024 +0200
@@ -27,18 +27,60 @@
 """
     coordinate_size(g)
 
-The lenght of the coordinate vector of `Grid` `g`.
+The length of the coordinate vector of `Grid` `g`.
 """
 coordinate_size(::Type{<:Grid{T}}) where T = _ncomponents(T)
 coordinate_size(g::Grid) = coordinate_size(typeof(g)) # TBD: Name of this function?!
 
 """
-    component_type(g)
+    component_type(gf)
+
+The type of the components of the elements of `gf`. I.e if `gf` is a vector
+valued grid function, `component_view(gf)` is the element type of the vectors
+at each grid point.
+
+# Examples
+```julia-repl
+julia> component_type([[1,2], [2,3], [3,4]])
+Int64
+```
+"""
+component_type(T::Type) = eltype(eltype(T))
+component_type(t) = component_type(typeof(t))
+
+"""
+    componentview(gf, component_index...)
+
+A view of `gf` with only the components specified by `component_index...`.
 
-The type of the components of the coordinate vector of `Grid` `g`.
+# Examples
+```julia-repl
+julia> componentview([[1,2], [2,3], [3,4]],2)
+3-element ArrayComponentView{Int64, Vector{Int64}, 1, Vector{Vector{Int64}}, Tuple{Int64}}:
+ 2
+ 3
+ 4
+```
 """
-component_type(::Type{<:Grid{T}}) where T = eltype(T)
-component_type(g::Grid) = component_type(typeof(g))
+componentview(gf, component_index...) = ArrayComponentView(gf, component_index)
+
+struct ArrayComponentView{CT,T,D,AT <: AbstractArray{T,D}, IT} <: AbstractArray{CT,D}
+    v::AT
+    component_index::IT
+
+    function ArrayComponentView(v, component_index)
+        CT = typeof(first(v)[component_index...])
+        return new{CT, eltype(v), ndims(v), typeof(v), typeof(component_index)}(v,component_index)
+    end
+end
+
+Base.size(cv) = size(cv.v)
+Base.getindex(cv::ArrayComponentView, i::Int) = cv.v[i][cv.component_index...]
+Base.getindex(cv::ArrayComponentView, I::Vararg{Int}) = cv.v[I...][cv.component_index...]
+IndexStyle(::Type{<:ArrayComponentView{<:Any,<:Any,AT}}) where AT = IndexStyle(AT)
+
+# TODO: Implement `setindex!`?
+# TODO: Implement a more general ComponentView that can handle non-AbstractArrays.
 
 """
     refine(g::Grid, r)
@@ -74,6 +116,18 @@
 # TBD: Can we implement a version here that accepts multiple ids and grouped boundaries? Maybe we need multiblock stuff?
 
 """
+    boundary_indices(g::Grid, id::BoundaryIdentifier)
+
+A collection of indices corresponding to the boundary with given id. For grids
+with Cartesian indexing these collections will be tuples with elements of type
+``Union{Int,Colon}``.
+
+When implementing this method it is expected that the returned collection can
+be used to index grid functions to obtain grid functions on the boundary grid.
+"""
+function boundary_indices end
+
+"""
     eval_on(g::Grid, f)
 
 Lazy evaluation of `f` on the grid. `f` can either be on the form `f(x,y,...)`
@@ -87,6 +141,8 @@
     if hasmethod(f, (Any,))
         return LazyTensors.LazyFunctionArray((I...)->f(g[I...]), size(g))
     else
+        # TBD This branch can be removed if we accept the trade off that we define f with the syntax f((x,y)) instead if we don't want to handle the vector in the body of f. (Add an example in the docs)
+        # Also see Notes.md
         return LazyTensors.LazyFunctionArray((I...)->f(g[I...]...), size(g))
     end
 end
@@ -100,3 +156,5 @@
 
 _ncomponents(::Type{<:Number}) = 1
 _ncomponents(T::Type{<:SVector}) = length(T)
+
+