diff src/Grids/mapped_grid.jl @ 1954:b0915f43b122 feature/sbp_operators/laplace_curvilinear

Merge feature/grids/geometry_functions
author Jonatan Werpers <jonatan@werpers.com>
date Sat, 08 Feb 2025 09:38:58 +0100
parents e97f4352b8d0
children
line wrap: on
line diff
--- a/src/Grids/mapped_grid.jl	Sat Feb 08 09:35:13 2025 +0100
+++ b/src/Grids/mapped_grid.jl	Sat Feb 08 09:38:58 2025 +0100
@@ -82,8 +82,8 @@
     jacobian_components = (:, free_variable_indices)
 
     # Create grid function for boundary grid jacobian
-    boundary_jacobian = componentview((@view g.jacobian[b_indices...])  , jacobian_components...)
-    boundary_physicalcoordinates = @view g.physicalcoordinates[b_indices...]
+    boundary_jacobian = componentview((@view g.jacobian[b_indices])  , jacobian_components...)
+    boundary_physicalcoordinates = @view g.physicalcoordinates[b_indices]
 
     return MappedGrid(
         boundary_grid(g.logical_grid, id),
@@ -92,27 +92,25 @@
     )
 end
 
-
 """
-    mapped_grid(x, J, size::Vararg{Int})
+    mapped_grid(x, J, size...)
 
 A `MappedGrid` with a default logical grid on the D-dimensional unit hyper 
-box [0,1]ᴰ. `x` and `J`are functions to be evaluated on the logical grid
+box [0,1]ᴰ. `x` and `J` are functions to be evaluated on the logical grid
 and `size` determines the size of the logical grid.
 """
 function mapped_grid(x, J, size::Vararg{Int})
     D = length(size)
-    lg = equidistant_grid(ntuple(i->0., D), ntuple(i->1., D), size...)
-    return mapped_grid(lg, x, J)
+    return mapped_grid(x, J, unithyperbox(D), size...)
 end
 
 """
-    mapped_grid(lg::Grid, x, J)
+    mapped_grid(x, J, lg::Grid)
 
 A `MappedGrid` with logical grid `lg`. Physical coordinates and Jacobian are
 determined by the functions `x` and `J`.
 """
-function mapped_grid(lg::Grid, x, J)
+function mapped_grid(x, J, lg::Grid)
     return MappedGrid(
         lg,
         map(x,lg),
@@ -121,6 +119,17 @@
 end
 
 """
+    mapped_grid(x, J, ps::ParameterSpace, size...)
+
+A `MappedGrid` with logical grid `lg`. Physical coordinates and Jacobian are
+determined by the functions `x` and `J`.
+"""
+function mapped_grid(x, J, ps::ParameterSpace, size::Vararg{Int})
+    lg = equidistant_grid(ps, size...)
+    return mapped_grid(x, J, lg)
+end
+
+"""
     metric_tensor(g::MappedGrid)
 
 The metric tensor of `g` as a grid function.
@@ -171,15 +180,16 @@
 The outward pointing normal as a grid function on the corresponding boundary grid.
 """
 function normal(g::MappedGrid, boundary)
-    b_indices = boundary_indices(g, boundary)
-    σ = _boundary_sign(component_type(g), boundary)
-    return map(jacobian(g)[b_indices...]) do ∂x∂ξ
-        ∂ξ∂x = inv(∂x∂ξ)
-        k = grid_id(boundary)
-        σ*∂ξ∂x[k,:]/norm(∂ξ∂x[k,:])
+    return map(boundary_indices(g, boundary)) do I
+        normal(g, boundary, Tuple(I)...)
     end
 end
 
+"""
+    normal(g::MappedGrid, boundary, i...)
+
+The outward pointing normal to the specified boundary in grid point `i`.
+"""
 function normal(g::MappedGrid, boundary, i...)
     σ = _boundary_sign(component_type(g), boundary)
     ∂ξ∂x = inv(jacobian(g)[i...])