comparison src/Grids/mapped_grid.jl @ 1506:535f32316637 feature/grids/curvilinear

Rename from curvilinear to mapped
author Jonatan Werpers <jonatan@werpers.com>
date Fri, 16 Feb 2024 15:28:19 +0100
parents src/Grids/curvilinear_grid.jl@63101a1cd0e6
children 69790e9d1652
comparison
equal deleted inserted replaced
1505:63101a1cd0e6 1506:535f32316637
1 struct MappedGrid{T,D, GT<:Grid{<:Any,D}, CT<:AbstractArray{T,D}, JT<:AbstractArray{<:AbstractArray{<:Any, 2}, D}} <: Grid{T,D}
2 logicalgrid::GT
3 physicalcoordinates::CT
4 jacobian::JT
5 end
6
7 jacobian(g::MappedGrid) = g.jacobian
8 logicalgrid(g::MappedGrid) = g.logicalgrid
9
10
11 # Indexing interface
12 Base.getindex(g::MappedGrid, I::Vararg{Int}) = g.physicalcoordinates[I...]
13 Base.eachindex(g::MappedGrid) = eachindex(g.logicalgrid)
14
15 Base.firstindex(g::MappedGrid, d) = firstindex(g.logicalgrid, d)
16 Base.lastindex(g::MappedGrid, d) = lastindex(g.logicalgrid, d)
17
18 # Iteration interface
19
20 Base.iterate(g::MappedGrid) = iterate(g.physicalcoordinates)
21 Base.iterate(g::MappedGrid, state) = iterate(g.physicalcoordinates, state)
22
23 Base.IteratorSize(::Type{<:MappedGrid{<:Any, D}}) where D = Base.HasShape{D}()
24 Base.length(g::MappedGrid) = length(g.logicalgrid)
25 Base.size(g::MappedGrid) = size(g.logicalgrid)
26 Base.size(g::MappedGrid, d) = size(g.logicalgrid, d)
27
28 boundary_identifiers(g::MappedGrid) = boundary_identifiers(g.logicalgrid)
29 boundary_indices(g::MappedGrid, id::TensorGridBoundary) = boundary_indices(g.logicalgrid, id)
30
31 function boundary_grid(g::MappedGrid, id::TensorGridBoundary)
32 b_indices = boundary_indices(g.logicalgrid, id)
33
34 # Calculate indices of needed jacobian components
35 D = ndims(g)
36 all_indices = SVector{D}(1:D)
37 free_variable_indices = deleteat(all_indices, grid_id(id))
38 jacobian_components = (:, free_variable_indices)
39
40 # Create grid function for boundary grid jacobian
41 boundary_jacobian = componentview((@view g.jacobian[b_indices...]) , jacobian_components...)
42 boundary_physicalcoordinates = @view g.physicalcoordinates[b_indices...]
43
44 return MappedGrid(
45 boundary_grid(g.logicalgrid, id),
46 boundary_physicalcoordinates,
47 boundary_jacobian,
48 )
49 end
50
51 function mapped_grid(x, J, size...)
52 D = length(size)
53 lg = equidistant_grid(size, ntuple(i->0., D), ntuple(i->1., D))
54 return MappedGrid(
55 lg,
56 map(x,lg),
57 map(J,lg),
58 )
59 end