Mercurial > repos > public > sbplib_julia

comparison src/Grids/grid.jl @ 2057:8a2a0d678d6f feature/lazy_tensors/pretty_printing

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Merge default
author Jonatan Werpers <jonatan@werpers.com>
date Tue, 10 Feb 2026 22:41:19 +0100
parents 491887181b8c
children
comparison
equal deleted inserted replaced
1110:c0bff9f6e0fb 2057:8a2a0d678d6f
1 """
2 Grid{T,D}
3
4 A grid with coordinates of type `T`, e.g. `SVector{3,Float64}`, and dimension
5 `D`. The grid can be embedded in a higher dimension in which case the number
6 of indices and the number of components of the coordinate vectors will be
7 different.
8
9 All grids are expected to behave as a grid function for the coordinates.
10
11 `Grids` is top level abstract type for grids. A grid should implement Julia's interfaces for
12 indexing and iteration.
13
14 ## Note
15
16 Importantly a grid does not have to be an `AbstractArray`. The reason is to
17 allow flexible handling of special types of grids like multi-block grids, or
18 grids with special indexing.
19 """
20 abstract type Grid{T,D} end
21
22 Base.ndims(::Grid{T,D}) where {T,D} = D
23 Base.eltype(::Type{<:Grid{T}}) where T = T
24
25 Base.getindex(g::Grid, I::CartesianIndex) = g[Tuple(I)...]
26
27 """
28 coordinate_size(g)
29
30 The length of the coordinate vector of `Grid` `g`.
31 """
32 coordinate_size(::Type{<:Grid{T}}) where T = _ncomponents(T)
33 coordinate_size(g::Grid) = coordinate_size(typeof(g)) # TBD: Name of this function?!
34
35 """
36 component_type(gf)
37
38 The type of the components of the elements of `gf`. I.e if `gf` is a vector
39 valued grid function, `component_view(gf)` is the element type of the vectors
40 at each grid point.
41
42 # Examples
43 ```julia-repl
44 julia> component_type([[1,2], [2,3], [3,4]])
45 Int64
46 ```
47 """
48 component_type(T::Type) = eltype(eltype(T))
49 component_type(t) = component_type(typeof(t))
50
51 """
52 componentview(gf, component_index...)
53
54 A view of `gf` with only the components specified by `component_index...`.
55
56 # Examples
57 ```julia-repl
58 julia> componentview([[1,2], [2,3], [3,4]],2)
59 3-element ArrayComponentView{Int64, Vector{Int64}, 1, Vector{Vector{Int64}}, Tuple{Int64}}:
60 2
61 3
62 4
63 ```
64 """
65 componentview(gf, component_index...) = ArrayComponentView(gf, component_index)
66
67 struct ArrayComponentView{CT,T,D,AT <: AbstractArray{T,D}, IT} <: AbstractArray{CT,D}
68 v::AT
69 component_index::IT
70
71 function ArrayComponentView(v, component_index)
72 CT = typeof(first(v)[component_index...])
73 return new{CT, eltype(v), ndims(v), typeof(v), typeof(component_index)}(v,component_index)
74 end
75 end
76
77 _array_type(v::ArrayComponentView) = _array_type(typeof(v))
78 _array_type(::Type{ArrayComponentView{CT,T,D,AT,IT}}) where {CT,T,D,AT,IT} = AT
79
80 Base.size(cv::ArrayComponentView) = size(cv.v)
81 Base.getindex(cv::ArrayComponentView, i::Int) = cv.v[i][cv.component_index...]
82 Base.getindex(cv::ArrayComponentView, I::Vararg{Int}) = cv.v[I...][cv.component_index...]
83 Base.IndexStyle(ACT::Type{<:ArrayComponentView}) = IndexStyle(_array_type(ACT))
84
85 # TODO: Implement `setindex!`?
86 # TODO: Implement a more general ComponentView that can handle non-AbstractArrays.
87
88
89 """
90 min_spacing(g::Grid)
91
92 The smallest distance between any pair of grid points in `g`.
93 """
94 function min_spacing end
95
96 """
97 refine(g::Grid, r)
98
99 The grid where `g` is refined by the factor `r`.
100
101 See also: [`coarsen`](@ref).
102 """
103 function refine end
104
105 """
106 coarsen(g::Grid, r)
107
108 The grid where `g` is coarsened by the factor `r`.
109
110 See also: [`refine`](@ref).
111 """
112 function coarsen end
113
114 """
115 BoundaryIdentifier
116
117 An identifier for a boundary of a grid.
118 """
119 abstract type BoundaryIdentifier end
120
121 """
122 boundary_identifiers(g::Grid)
123
124 Identifiers for all the boundaries of `g`.
125 """
126 function boundary_identifiers end
127
128 """
129 boundary_grid(g::Grid, id::BoundaryIdentifier)
130
131 The grid for the boundary specified by `id`.
132 """
133 function boundary_grid end
134 # TBD: Can we implement a version here that accepts multiple ids and grouped boundaries? Maybe we need multiblock stuff?
135
136 """
137 boundary_indices(g::Grid, id::BoundaryIdentifier)
138
139 A collection of indices corresponding to the boundary with given id. For grids
140 with Cartesian indexing these collections will be tuples with elements of type
141 ``Union{Int,Colon}``.
142
143 When implementing this method it is expected that the returned collection can
144 be used to index grid functions to obtain grid functions on the boundary grid.
145 """
146 function boundary_indices end
147
148 """
149 eval_on(g::Grid, f)
150
151 Lazy evaluation of `f` on the grid. `f` can either be on the form `f(x,y,...)`
152 with each coordinate as an argument, or on the form `f(x̄)` taking a
153 coordinate vector.
154
155 For concrete array grid functions `map(f,g)` can be used instead.
156 """
157 eval_on(g::Grid, f) = eval_on(g, f, Base.IteratorSize(g))
158 function eval_on(g::Grid, f, ::Base.HasShape)
159 if hasmethod(f, (Any,))
160 return LazyTensors.LazyFunctionArray((I...)->f(g[I...]), size(g))
161 else
162 # 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)
163 # Also see Notes.md
164 return LazyTensors.LazyFunctionArray((I...)->f(g[I...]...), size(g))
165 end
166 end
167
168 """
169 eval_on(g::Grid, f::Number)
170
171 Lazy evaluation of a scalar `f` on the grid.
172 """
173 eval_on(g::Grid, f::Number) = return LazyTensors.LazyConstantArray(f, size(g))
174
175 _ncomponents(::Type{<:Number}) = 1
176 _ncomponents(T::Type{<:SVector}) = length(T)
177
178