From cfb240775af173c705ddef5d7ece6dc75d44bfbc Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva Date: Sat, 1 Nov 2025 20:45:25 -0300 Subject: [PATCH 01/27] Add fast path for indices on non-invertible TransformedGrid --- src/indices.jl | 14 ++++++++++++++ test/indices.jl | 14 ++++++++++++++ 2 files changed, 28 insertions(+) diff --git a/src/indices.jl b/src/indices.jl index 722379bd7..53ba867f1 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -74,6 +74,20 @@ function indices(grid::OrthoRectilinearGrid, box::Box) LinearIndices(size(grid))[range] |> vec end +function indices(grid::TransformedGrid, geometry::Geometry) + g = grid.mesh + t = grid.transform + + # find the invertible part (e.g., Rotate, Affine) + invertible_part = only(x for x in t if isinvertible(x)) + + # un-transform the geometry using only that part + p = revert(invertible_part, geometry, nothing) + + # find indices on the non-deformed mesh + indices(g, p) +end + # ---------------- # CARTESIAN RANGE # ---------------- diff --git a/test/indices.jl b/test/indices.jl index 6d5f2ccc4..82caebfc5 100644 --- a/test/indices.jl +++ b/test/indices.jl @@ -202,4 +202,18 @@ @test linds[2, 2] == only(indices(grid, p4)) @test linds[10, 10] == only(indices(grid, p5)) @test isempty(indices(grid, p6)) + + # transformed grid + g = cartgrid(20, 20) + p_orig = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) + expected_inds = indices(g, p_orig) + linds = LinearIndices(size(g)) + @test linds[12, 6] ∈ expected_inds + t_aff = Rotate(Angle2d(T(π / 2))) + t_morph = Morphological(c -> Cartesian(c.x + c.y, c.y, c.x - c.y)) + t_seq = t_aff → t_morph + grid = TransformedGrid(g, t_seq) + poly = t_aff(p_orig) + actual_inds = indices(grid, poly) + @test actual_inds == expected_inds end From dceddeb84dc9a1c0c0057e5aadd0ae0888b09383 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:15:00 -0300 Subject: [PATCH 02/27] Update src/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- src/indices.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/indices.jl b/src/indices.jl index 53ba867f1..3d32053d2 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -75,8 +75,8 @@ function indices(grid::OrthoRectilinearGrid, box::Box) end function indices(grid::TransformedGrid, geometry::Geometry) - g = grid.mesh - t = grid.transform + g = parent(grid) + t = transform(grid) # find the invertible part (e.g., Rotate, Affine) invertible_part = only(x for x in t if isinvertible(x)) From 6f0baf069d6ad4e4ca2101b5f85f872d99a77597 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:15:26 -0300 Subject: [PATCH 03/27] Update src/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- src/indices.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/indices.jl b/src/indices.jl index 3d32053d2..f9b929fdf 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -79,7 +79,7 @@ function indices(grid::TransformedGrid, geometry::Geometry) t = transform(grid) # find the invertible part (e.g., Rotate, Affine) - invertible_part = only(x for x in t if isinvertible(x)) + r = reduce(→, reverse(filter(isrevertible, t))) # un-transform the geometry using only that part p = revert(invertible_part, geometry, nothing) From c24448cb0b9e960f7c23adbf0b6eb101ff37ced7 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:15:36 -0300 Subject: [PATCH 04/27] Update src/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- src/indices.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/indices.jl b/src/indices.jl index f9b929fdf..d299f5485 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -78,7 +78,7 @@ function indices(grid::TransformedGrid, geometry::Geometry) g = parent(grid) t = transform(grid) - # find the invertible part (e.g., Rotate, Affine) + # construct reverse transform from revertible steps r = reduce(→, reverse(filter(isrevertible, t))) # un-transform the geometry using only that part From fcb44c6b34b21695c9b4616b2810207aa4ee744c Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:16:04 -0300 Subject: [PATCH 05/27] Update src/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- src/indices.jl | 7 ++----- 1 file changed, 2 insertions(+), 5 deletions(-) diff --git a/src/indices.jl b/src/indices.jl index d299f5485..f46e962eb 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -81,11 +81,8 @@ function indices(grid::TransformedGrid, geometry::Geometry) # construct reverse transform from revertible steps r = reduce(→, reverse(filter(isrevertible, t))) - # un-transform the geometry using only that part - p = revert(invertible_part, geometry, nothing) - - # find indices on the non-deformed mesh - indices(g, p) + # find indices before transformation + indices(g, r(geometry)) end # ---------------- From 5d092014991ffecf4688d6822faaa1e9138ef578 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:30:30 -0300 Subject: [PATCH 06/27] Update src/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- src/indices.jl | 9 +++------ 1 file changed, 3 insertions(+), 6 deletions(-) diff --git a/src/indices.jl b/src/indices.jl index f46e962eb..0724a62ab 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -75,14 +75,11 @@ function indices(grid::OrthoRectilinearGrid, box::Box) end function indices(grid::TransformedGrid, geometry::Geometry) - g = parent(grid) - t = transform(grid) - # construct reverse transform from revertible steps - r = reduce(→, reverse(filter(isrevertible, t))) + revtrans = reduce(→, reverse(filter(isrevertible, transform(t)))) - # find indices before transformation - indices(g, r(geometry)) + # find indices in non-transformed space + indices(parent(grid), revtrans(geometry)) end # ---------------- From 0af0ce50443865567fa3beea49a01ea09c6ad767 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva Date: Sat, 1 Nov 2025 21:31:03 -0300 Subject: [PATCH 07/27] fix test names --- test/indices.jl | 20 ++++++++++---------- 1 file changed, 10 insertions(+), 10 deletions(-) diff --git a/test/indices.jl b/test/indices.jl index 82caebfc5..44f7da98c 100644 --- a/test/indices.jl +++ b/test/indices.jl @@ -205,15 +205,15 @@ # transformed grid g = cartgrid(20, 20) - p_orig = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) - expected_inds = indices(g, p_orig) + poly1 = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) + expectedinds = indices(g, poly1) linds = LinearIndices(size(g)) - @test linds[12, 6] ∈ expected_inds - t_aff = Rotate(Angle2d(T(π / 2))) - t_morph = Morphological(c -> Cartesian(c.x + c.y, c.y, c.x - c.y)) - t_seq = t_aff → t_morph - grid = TransformedGrid(g, t_seq) - poly = t_aff(p_orig) - actual_inds = indices(grid, poly) - @test actual_inds == expected_inds + @test linds[12, 6] ∈ expectedinds + taff = Rotate(Angle2d(T(π / 2))) + tmorph = Morphological(c -> Cartesian(c.x + c.y, c.y, c.x - c.y)) + tseq = taff → tmorph + grid = TransformedGrid(g, tseq) + poly2 = taff(poly1) + actualinds = indices(grid, poly2) + @test actualinds == expectedinds end From 67f4573e8a3f921fc8774d1a9fc26017000a9360 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:51:14 -0300 Subject: [PATCH 08/27] Update test/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- test/indices.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/indices.jl b/test/indices.jl index 44f7da98c..4e78faa9f 100644 --- a/test/indices.jl +++ b/test/indices.jl @@ -204,7 +204,7 @@ @test isempty(indices(grid, p6)) # transformed grid - g = cartgrid(20, 20) + grid = cartgrid(20, 20) poly1 = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) expectedinds = indices(g, poly1) linds = LinearIndices(size(g)) From a9325df6775aad61d65def7e7a11edf2d78105f0 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:51:21 -0300 Subject: [PATCH 09/27] Update test/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- test/indices.jl | 19 ++++++++----------- 1 file changed, 8 insertions(+), 11 deletions(-) diff --git a/test/indices.jl b/test/indices.jl index 4e78faa9f..519eb4d8e 100644 --- a/test/indices.jl +++ b/test/indices.jl @@ -205,15 +205,12 @@ # transformed grid grid = cartgrid(20, 20) - poly1 = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) - expectedinds = indices(g, poly1) - linds = LinearIndices(size(g)) - @test linds[12, 6] ∈ expectedinds - taff = Rotate(Angle2d(T(π / 2))) - tmorph = Morphological(c -> Cartesian(c.x + c.y, c.y, c.x - c.y)) - tseq = taff → tmorph - grid = TransformedGrid(g, tseq) - poly2 = taff(poly1) - actualinds = indices(grid, poly2) - @test actualinds == expectedinds + poly = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) + inds = indices(grid, poly) + linds = LinearIndices(size(grid)) + @test linds[12, 6] ∈ indices(grid, poly) + trans = Rotate(Angle2d(T(π / 2))) → Morphological(c -> Cartesian(c.x + c.y, c.y, c.x - c.y)) + tgrid = trans(grid) + tpoly = trans(poly) + @test indices(tgrid, tpoly) == indices(grid, poly) end From e0f3f4237550d8d80ede0e46bc204bbdf5680ded Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva Date: Sat, 1 Nov 2025 21:53:21 -0300 Subject: [PATCH 10/27] fix indices transform --- src/indices.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/indices.jl b/src/indices.jl index 0724a62ab..6bd058cb6 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -74,9 +74,9 @@ function indices(grid::OrthoRectilinearGrid, box::Box) LinearIndices(size(grid))[range] |> vec end -function indices(grid::TransformedGrid, geometry::Geometry) +function indice(grid::TransformedGrid, geometry::Geometry) # construct reverse transform from revertible steps - revtrans = reduce(→, reverse(filter(isrevertible, transform(t)))) + revtrans = reduce(→, reverse(filter(isrevertible, transform(grid)))) # find indices in non-transformed space indices(parent(grid), revtrans(geometry)) From 472f37d2d111811e369e0f096ef5859c0e849f66 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva Date: Sat, 1 Nov 2025 21:53:59 -0300 Subject: [PATCH 11/27] fix transform in indices --- src/indices.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/indices.jl b/src/indices.jl index 6bd058cb6..0b468211d 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -74,7 +74,7 @@ function indices(grid::OrthoRectilinearGrid, box::Box) LinearIndices(size(grid))[range] |> vec end -function indice(grid::TransformedGrid, geometry::Geometry) +function indices(grid::TransformedGrid, geometry::Geometry) # construct reverse transform from revertible steps revtrans = reduce(→, reverse(filter(isrevertible, transform(grid)))) From c7c91f967fcb868cac7b2a34b9aeecfad09490a8 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva Date: Wed, 19 Nov 2025 17:45:02 -0300 Subject: [PATCH 12/27] _reversetransform and indices --- src/Meshes.jl | 882 ++++++++++++++++++++++++------------------------ src/indices.jl | 420 ++++++++++++----------- test/indices.jl | 418 +++++++++++------------ 3 files changed, 864 insertions(+), 856 deletions(-) diff --git a/src/Meshes.jl b/src/Meshes.jl index 5ce0a66e0..d5d690e07 100644 --- a/src/Meshes.jl +++ b/src/Meshes.jl @@ -32,7 +32,7 @@ using ScopedValues: ScopedValue using Base.Cartesian: @nloops, @nref, @ntuple using Base: @propagate_inbounds -import Random +using Random: Random import Base: sort import Base: ==, ! import Base: +, -, * @@ -41,7 +41,7 @@ import Distances: evaluate import NearestNeighbors: MinkowskiMetric # Transforms API -import TransformsBase: Transform, → +import TransformsBase: Transform, →, SequentialTransform, Identity import TransformsBase: isrevertible, isinvertible import TransformsBase: apply, revert, reapply, inverse import TransformsBase: parameters, preprocess @@ -78,6 +78,9 @@ include("domains.jl") # utilities include("utils.jl") +# transforms +include("transforms.jl") + # domain indices include("indices.jl") @@ -120,8 +123,7 @@ include("discretization.jl") include("refinement.jl") include("coarsening.jl") -# transforms -include("transforms.jl") + # miscellaneous include("rand.jl") @@ -134,441 +136,441 @@ include("projecting.jl") include("viz.jl") export - # vectors - Vec, - ∠, - ⋅, - ×, - - # manifolds - 𝔼, - 🌐, - - # geometries - Geometry, - embeddim, - paramdim, - crs, - manifold, - - # primitives - Primitive, - Point, - Ray, - Line, - BezierCurve, - ParametrizedCurve, - Plane, - Box, - Ball, - Sphere, - Ellipsoid, - Disk, - Circle, - Cylinder, - CylinderSurface, - ParaboloidSurface, - Cone, - ConeSurface, - Frustum, - FrustumSurface, - Torus, - controls, - ncontrols, - degree, - Horner, - DeCasteljau, - coords, - to, - center, - radius, - radii, - plane, - bottom, - top, - axis, - base, - apex, - isright, - sides, - diagonal, - focallength, - direction, - - # polytopes - Polytope, - Chain, - Segment, - Rope, - Ring, - Polygon, - Ngon, - Triangle, - Quadrangle, - Pentagon, - Hexagon, - Heptagon, - Octagon, - Nonagon, - Decagon, - PolyArea, - Polyhedron, - Tetrahedron, - Hexahedron, - Pyramid, - Wedge, - vertex, - vertices, - nvertices, - eachvertex, - rings, - segments, - angles, - innerangles, - normal, - ≗, - - # multi-geometries - Multi, - MultiPoint, - MultiSegment, - MultiRope, - MultiRing, - MultiPolygon, - MultiPolyhedron, - - # transformed geometry - TransformedGeometry, - - # connectivities - Connectivity, - paramdim, - indices, - connect, - materialize, - pltype, - - # topologies - Topology, - vertex, - vertices, - nvertices, - element, - elements, - nelements, - facet, - facets, - nfacets, - faces, - nfaces, - GridTopology, - HalfEdgeTopology, - HalfEdge, - SimpleTopology, - elem2cart, - cart2elem, - corner2cart, - cart2corner, - elem2corner, - corner2elem, - elementtype, - facettype, - half4elem, - half4vert, - half4edge, - half4pair, - edge4pair, - connec4elem, - - # topological relations - TopologicalRelation, - Boundary, - Coboundary, - Adjacency, - - # domain traits - Domain, - SubDomain, - TransformedDomain, - embeddim, - paramdim, - crs, - manifold, - element, - nelements, - - # sets - GeometrySet, - PointSet, - - # meshes - Mesh, - SimpleMesh, - TransformedMesh, - Grid, - SubGrid, - RegularGrid, - CartesianGrid, - RectilinearGrid, - StructuredGrid, - TransformedGrid, - vertex, - vertices, - nvertices, - eachvertex, - element, - elements, - nelements, - facet, - facets, - nfacets, - topology, - topoconvert, - spacing, - offset, - - # trajectories - CylindricalTrajectory, - - # indices - indices, - - # partitions - Partition, - indices, - metadata, - - # partitioning - PartitionMethod, - IndexPredicatePartitionMethod, - PointPredicatePartitionMethod, - UniformPartition, - FractionPartition, - BlockPartition, - BisectPointPartition, - BisectFractionPartition, - BallPartition, - PlanePartition, - DirectionPartition, - IndexPredicatePartition, - PointPredicatePartition, - ProductPartition, - HierarchicalPartition, - partition, - split, - - # sorting - SortingMethod, - DirectionSort, - sortinds, - - # traversing - Path, - LinearPath, - RandomPath, - ShiftedPath, - SourcePath, - MultiGridPath, - traverse, - - # neighborhoods - Neighborhood, - MetricBall, - hasequalradii, - radii, - rotation, - metric, - radius, - - # neighborhood search - NeighborSearchMethod, - BoundedNeighborSearchMethod, - BallSearch, - KNearestSearch, - KBallSearch, - search!, - searchdists!, - search, - searchdists, - maxneighbors, - - # predicates - isparametrized, - iscurve, - issurface, - issolid, - isperiodic, - issimplex, - isclosed, - isconvex, - issimple, - hasholes, - intersects, - iscollinear, - iscoplanar, - ≺, - ≻, - ⪯, - ⪰, - - # centroids - centroid, - - # measures - measure, - area, - volume, - perimeter, - - # boundary - boundary, - - # winding number - winding, - - # sideof - sideof, - SideType, - IN, - OUT, - ON, - LEFT, - RIGHT, - - # orientation - orientation, - OrientationType, - CW, - CCW, - - # clipping - ClippingMethod, - SutherlandHodgmanClipping, - clip, - - # intersections - IntersectionType, - Crossing, - CornerCrossing, - EdgeCrossing, - Touching, - CornerTouching, - EdgeTouching, - Overlapping, - PosOverlapping, - NegOverlapping, - NotIntersecting, - Intersecting, - - # intersecting - Intersection, - intersection, - type, - - # simplification - SimplificationMethod, - SelingerSimplification, - DouglasPeuckerSimplification, - MinMaxSimplification, - simplify, - - # bounding boxes - boundingbox, - - # hulls - HullMethod, - GrahamScan, - JarvisMarch, - hull, - convexhull, - - # sampling - SamplingMethod, - DiscreteSamplingMethod, - ContinuousSamplingMethod, - UniformSampling, - WeightedSampling, - BallSampling, - BlockSampling, - RegularSampling, - HomogeneousSampling, - MinDistanceSampling, - FibonacciSampling, - sampleinds, - sample, - - # pointification - pointify, - - # tesselation - TesselationMethod, - DelaunayTesselation, - VoronoiTesselation, - tesselate, - - # discretization - DiscretizationMethod, - BoundaryTriangulationMethod, - FanTriangulation, - DehnTriangulation, - HeldTriangulation, - DelaunayTriangulation, - ManualSimplexification, - RegularDiscretization, - MaxLengthDiscretization, - discretize, - discretizewithin, - simplexify, - - # refinement - RefinementMethod, - TriRefinement, - QuadRefinement, - RegularRefinement, - CatmullClarkRefinement, - TriSubdivision, - refine, - - # coarsening - CoarseningMethod, - RegularCoarsening, - coarsen, - - # transforms - GeometricTransform, - CoordinateTransform, - Rotate, - Translate, - Scale, - Affine, - Stretch, - StdCoords, - Proj, - Morphological, - LengthUnit, - ValidCoords, - Shadow, - Slice, - Repair, - Bridge, - LambdaMuSmoothing, - LaplaceSmoothing, - TaubinSmoothing, - isaffine, - isinvertible, - inverse, - - # miscellaneous - signarea, - householderbasis, - supportfun, - laplacematrix, - measurematrix, - adjacencymatrix, - atol, - - # visualization - viz, - viz! + # vectors + Vec, + ∠, + ⋅, + ×, + + # manifolds + 𝔼, + 🌐, + + # geometries + Geometry, + embeddim, + paramdim, + crs, + manifold, + + # primitives + Primitive, + Point, + Ray, + Line, + BezierCurve, + ParametrizedCurve, + Plane, + Box, + Ball, + Sphere, + Ellipsoid, + Disk, + Circle, + Cylinder, + CylinderSurface, + ParaboloidSurface, + Cone, + ConeSurface, + Frustum, + FrustumSurface, + Torus, + controls, + ncontrols, + degree, + Horner, + DeCasteljau, + coords, + to, + center, + radius, + radii, + plane, + bottom, + top, + axis, + base, + apex, + isright, + sides, + diagonal, + focallength, + direction, + + # polytopes + Polytope, + Chain, + Segment, + Rope, + Ring, + Polygon, + Ngon, + Triangle, + Quadrangle, + Pentagon, + Hexagon, + Heptagon, + Octagon, + Nonagon, + Decagon, + PolyArea, + Polyhedron, + Tetrahedron, + Hexahedron, + Pyramid, + Wedge, + vertex, + vertices, + nvertices, + eachvertex, + rings, + segments, + angles, + innerangles, + normal, + ≗, + + # multi-geometries + Multi, + MultiPoint, + MultiSegment, + MultiRope, + MultiRing, + MultiPolygon, + MultiPolyhedron, + + # transformed geometry + TransformedGeometry, + + # connectivities + Connectivity, + paramdim, + indices, + connect, + materialize, + pltype, + + # topologies + Topology, + vertex, + vertices, + nvertices, + element, + elements, + nelements, + facet, + facets, + nfacets, + faces, + nfaces, + GridTopology, + HalfEdgeTopology, + HalfEdge, + SimpleTopology, + elem2cart, + cart2elem, + corner2cart, + cart2corner, + elem2corner, + corner2elem, + elementtype, + facettype, + half4elem, + half4vert, + half4edge, + half4pair, + edge4pair, + connec4elem, + + # topological relations + TopologicalRelation, + Boundary, + Coboundary, + Adjacency, + + # domain traits + Domain, + SubDomain, + TransformedDomain, + embeddim, + paramdim, + crs, + manifold, + element, + nelements, + + # sets + GeometrySet, + PointSet, + + # meshes + Mesh, + SimpleMesh, + TransformedMesh, + Grid, + SubGrid, + RegularGrid, + CartesianGrid, + RectilinearGrid, + StructuredGrid, + TransformedGrid, + vertex, + vertices, + nvertices, + eachvertex, + element, + elements, + nelements, + facet, + facets, + nfacets, + topology, + topoconvert, + spacing, + offset, + + # trajectories + CylindricalTrajectory, + + # indices + indices, + + # partitions + Partition, + indices, + metadata, + + # partitioning + PartitionMethod, + IndexPredicatePartitionMethod, + PointPredicatePartitionMethod, + UniformPartition, + FractionPartition, + BlockPartition, + BisectPointPartition, + BisectFractionPartition, + BallPartition, + PlanePartition, + DirectionPartition, + IndexPredicatePartition, + PointPredicatePartition, + ProductPartition, + HierarchicalPartition, + partition, + split, + + # sorting + SortingMethod, + DirectionSort, + sortinds, + + # traversing + Path, + LinearPath, + RandomPath, + ShiftedPath, + SourcePath, + MultiGridPath, + traverse, + + # neighborhoods + Neighborhood, + MetricBall, + hasequalradii, + radii, + rotation, + metric, + radius, + + # neighborhood search + NeighborSearchMethod, + BoundedNeighborSearchMethod, + BallSearch, + KNearestSearch, + KBallSearch, + search!, + searchdists!, + search, + searchdists, + maxneighbors, + + # predicates + isparametrized, + iscurve, + issurface, + issolid, + isperiodic, + issimplex, + isclosed, + isconvex, + issimple, + hasholes, + intersects, + iscollinear, + iscoplanar, + ≺, + ≻, + ⪯, + ⪰, + + # centroids + centroid, + + # measures + measure, + area, + volume, + perimeter, + + # boundary + boundary, + + # winding number + winding, + + # sideof + sideof, + SideType, + IN, + OUT, + ON, + LEFT, + RIGHT, + + # orientation + orientation, + OrientationType, + CW, + CCW, + + # clipping + ClippingMethod, + SutherlandHodgmanClipping, + clip, + + # intersections + IntersectionType, + Crossing, + CornerCrossing, + EdgeCrossing, + Touching, + CornerTouching, + EdgeTouching, + Overlapping, + PosOverlapping, + NegOverlapping, + NotIntersecting, + Intersecting, + + # intersecting + Intersection, + intersection, + type, + + # simplification + SimplificationMethod, + SelingerSimplification, + DouglasPeuckerSimplification, + MinMaxSimplification, + simplify, + + # bounding boxes + boundingbox, + + # hulls + HullMethod, + GrahamScan, + JarvisMarch, + hull, + convexhull, + + # sampling + SamplingMethod, + DiscreteSamplingMethod, + ContinuousSamplingMethod, + UniformSampling, + WeightedSampling, + BallSampling, + BlockSampling, + RegularSampling, + HomogeneousSampling, + MinDistanceSampling, + FibonacciSampling, + sampleinds, + sample, + + # pointification + pointify, + + # tesselation + TesselationMethod, + DelaunayTesselation, + VoronoiTesselation, + tesselate, + + # discretization + DiscretizationMethod, + BoundaryTriangulationMethod, + FanTriangulation, + DehnTriangulation, + HeldTriangulation, + DelaunayTriangulation, + ManualSimplexification, + RegularDiscretization, + MaxLengthDiscretization, + discretize, + discretizewithin, + simplexify, + + # refinement + RefinementMethod, + TriRefinement, + QuadRefinement, + RegularRefinement, + CatmullClarkRefinement, + TriSubdivision, + refine, + + # coarsening + CoarseningMethod, + RegularCoarsening, + coarsen, + + # transforms + GeometricTransform, + CoordinateTransform, + Rotate, + Translate, + Scale, + Affine, + Stretch, + StdCoords, + Proj, + Morphological, + LengthUnit, + ValidCoords, + Shadow, + Slice, + Repair, + Bridge, + LambdaMuSmoothing, + LaplaceSmoothing, + TaubinSmoothing, + isaffine, + isinvertible, + inverse, + + # miscellaneous + signarea, + householderbasis, + supportfun, + laplacematrix, + measurematrix, + adjacencymatrix, + atol, + + # visualization + viz, + viz! end # module diff --git a/src/indices.jl b/src/indices.jl index 0b468211d..ae5bae4bf 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -7,79 +7,95 @@ # --------------- """ - indices(domain, geometry) + indices(domain, geometry) Return the indices of the elements of the `domain` that intersect with the `geometry`. """ indices(domain::Domain, geometry::Geometry) = findall(intersects(geometry), domain) function indices(grid::OrthoRegularGrid, point::Point) - point ∉ grid && return Int[] + point ∉ grid && return Int[] - # grid properties - orig = minimum(grid) - spac = spacing(grid) - dims = size(grid) + # grid properties + orig = minimum(grid) + spac = spacing(grid) + dims = size(grid) - # integer coordinates - coords = ceil.(Int, (point - orig) ./ spac) + # integer coordinates + coords = ceil.(Int, (point - orig) ./ spac) - # fix coordinates that are on the grid border - coords = clamp.(coords, 1, dims) + # fix coordinates that are on the grid border + coords = clamp.(coords, 1, dims) - # convert to linear index - [LinearIndices(dims)[coords...]] + # convert to linear index + [LinearIndices(dims)[coords...]] end function indices(grid::OrthoRegularGrid, chain::Chain) - dims = size(grid) - mask = falses(dims) + dims = size(grid) + mask = falses(dims) - for segment in segments(chain) - p₁, p₂ = vertices(segment) - _bresenham!(mask, grid, true, p₁, p₂) - end + for segment in segments(chain) + p₁, p₂ = vertices(segment) + _bresenham!(mask, grid, true, p₁, p₂) + end - LinearIndices(dims)[mask] + LinearIndices(dims)[mask] end function indices(grid::OrthoRegularGrid, poly::Polygon) - dims = size(grid) - mask = zeros(Int, dims) - cpoly = poly ∩ boundingbox(grid) - isnothing(cpoly) && return Int[] + dims = size(grid) + mask = zeros(Int, dims) + cpoly = poly ∩ boundingbox(grid) + isnothing(cpoly) && return Int[] - for (i, triangle) in enumerate(simplexify(cpoly)) - _fill!(mask, grid, i, triangle) - end + for (i, triangle) in enumerate(simplexify(cpoly)) + _fill!(mask, grid, i, triangle) + end - LinearIndices(dims)[mask .> 0] + LinearIndices(dims)[mask .> 0] end function indices(grid::OrthoRegularGrid, box::Box) - # cartesian range - range = cartesianrange(grid, box) + # cartesian range + range = cartesianrange(grid, box) - # convert to linear indices - LinearIndices(size(grid))[range] |> vec + # convert to linear indices + LinearIndices(size(grid))[range] |> vec end indices(grid::OrthoRegularGrid, multi::Multi) = mapreduce(geom -> indices(grid, geom), vcat, parent(multi)) |> unique function indices(grid::OrthoRectilinearGrid, box::Box) - # cartesian range - range = cartesianrange(grid, box) + # cartesian range + range = cartesianrange(grid, box) - # convert to linear indices - LinearIndices(size(grid))[range] |> vec + # convert to linear indices + LinearIndices(size(grid))[range] |> vec +end + +# Fallback: inverse if possible, otherwise error +function _reversetransform(t::Transform) + isrevertible(t) ? inverse(t) : throw(ArgumentError("Irreversible transform: $(typeof(t))")) +end + +# Treat specific irreversible transforms as Identity (no-op) +_reversetransform(::Morphological) = Identity() +_reversetransform(::Repair) = Identity() +_reversetransform(::Bridge) = Identity() + +function _reversetransform(st::SequentialTransform) + # Recursively reverse steps, replacing exceptions with Identity + newsteps = map(_reversetransform, st.transforms) + reduce(→, reverse(newsteps)) end function indices(grid::TransformedGrid, geometry::Geometry) - # construct reverse transform from revertible steps - revtrans = reduce(→, reverse(filter(isrevertible, transform(grid)))) + # construct reverse transform from revertible steps + revtrans = _reversetransform(transform(grid)) - # find indices in non-transformed space - indices(parent(grid), revtrans(geometry)) + # find indices in non-transformed space + indices(parent(grid), revtrans(geometry)) end # ---------------- @@ -87,7 +103,7 @@ end # ---------------- """ - cartesianrange(grid, box) + cartesianrange(grid, box) Return the Cartesian range of the elements of the `grid` that intersect with the `box`. """ @@ -98,97 +114,97 @@ _manifoldrange(::Type{<:𝔼}, grid::Grid, box::Box) = _euclideanrange(grid, box _manifoldrange(::Type{<:🌐}, grid::Grid, box::Box) = _geodesicrange(grid, box) function _euclideanrange(grid::OrthoRegularGrid, box::Box) - # grid properties - or = minimum(grid) - sp = spacing(grid) - sz = size(grid) + # grid properties + or = minimum(grid) + sp = spacing(grid) + sz = size(grid) - # intersection of boxes - lo, up = extrema(boundingbox(grid) ∩ box) + # intersection of boxes + lo, up = extrema(boundingbox(grid) ∩ box) - # Cartesian indices of new corners - ijkₛ = max.(ceil.(Int, (lo - or) ./ sp), 1) - ijkₑ = min.(floor.(Int, (up - or) ./ sp) .+ 1, sz) + # Cartesian indices of new corners + ijkₛ = max.(ceil.(Int, (lo - or) ./ sp), 1) + ijkₑ = min.(floor.(Int, (up - or) ./ sp) .+ 1, sz) - # Cartesian range from corner to corner - CartesianIndex(Tuple(ijkₛ)):CartesianIndex(Tuple(ijkₑ)) + # Cartesian range from corner to corner + CartesianIndex(Tuple(ijkₛ)):CartesianIndex(Tuple(ijkₑ)) end function _euclideanrange(grid::OrthoRectilinearGrid, box::Box) - # grid properties - nd = paramdim(grid) + # grid properties + nd = paramdim(grid) - # intersection of boxes - lo, up = to.(extrema(boundingbox(grid) ∩ box)) + # intersection of boxes + lo, up = to.(extrema(boundingbox(grid) ∩ box)) - # integer coordinates of lower point - ijkₛ = ntuple(nd) do i - findlast(x -> x ≤ lo[i], xyz(grid)[i]) - end + # integer coordinates of lower point + ijkₛ = ntuple(nd) do i + findlast(x -> x ≤ lo[i], xyz(grid)[i]) + end - # integer coordinates of upper point - ijkₑ = ntuple(nd) do i - findfirst(x -> x ≥ up[i], xyz(grid)[i]) - end + # integer coordinates of upper point + ijkₑ = ntuple(nd) do i + findfirst(x -> x ≥ up[i], xyz(grid)[i]) + end - # integer coordinates of elements - CartesianIndex(ijkₛ):CartesianIndex(ijkₑ .- 1) + # integer coordinates of elements + CartesianIndex(ijkₛ):CartesianIndex(ijkₑ .- 1) end function _geodesicrange(grid::Grid, box::Box) - nlon, nlat = vsize(grid) - - boxmin = convert(LatLon, coords(minimum(box))) - boxmax = convert(LatLon, coords(maximum(box))) - - a = convert(LatLon, coords(vertex(grid, (1, 1)))) - b = convert(LatLon, coords(vertex(grid, (nlon, 1)))) - c = convert(LatLon, coords(vertex(grid, (1, nlat)))) - - swaplon = a.lon > b.lon - swaplat = a.lat > c.lat - - loninds = swaplon ? (nlon:-1:1) : (1:1:nlon) - latinds = swaplat ? (nlat:-1:1) : (1:1:nlat) - - gridlonₛ, gridlonₑ = swaplon ? (b.lon, a.lon) : (a.lon, b.lon) - gridlatₛ, gridlatₑ = swaplat ? (c.lat, a.lat) : (a.lat, c.lat) - - lonmin = max(boxmin.lon, gridlonₛ) - latmin = max(boxmin.lat, gridlatₛ) - lonmax = min(boxmax.lon, gridlonₑ) - latmax = min(boxmax.lat, gridlatₑ) - - iₛ = findlast(loninds) do i - p = vertex(grid, (i, 1)) - c = convert(LatLon, coords(p)) - c.lon ≤ lonmin - end - iₑ = findfirst(loninds) do i - p = vertex(grid, (i, 1)) - c = convert(LatLon, coords(p)) - c.lon ≥ lonmax - end - - jₛ = findlast(latinds) do i - p = vertex(grid, (1, i)) - c = convert(LatLon, coords(p)) - c.lat ≤ latmin - end - jₑ = findfirst(latinds) do i - p = vertex(grid, (1, i)) - c = convert(LatLon, coords(p)) - c.lat ≥ latmax - end - - if iₛ == iₑ || jₛ == jₑ - throw(ArgumentError("the passed limits are not valid for the grid")) - end - - iₛ, iₑ = swaplon ? (iₑ, iₛ) : (iₛ, iₑ) - jₛ, jₑ = swaplat ? (jₑ, jₛ) : (jₛ, jₑ) - - CartesianIndex(loninds[iₛ], latinds[jₛ]):CartesianIndex(loninds[iₑ] - 1, latinds[jₑ] - 1) + nlon, nlat = vsize(grid) + + boxmin = convert(LatLon, coords(minimum(box))) + boxmax = convert(LatLon, coords(maximum(box))) + + a = convert(LatLon, coords(vertex(grid, (1, 1)))) + b = convert(LatLon, coords(vertex(grid, (nlon, 1)))) + c = convert(LatLon, coords(vertex(grid, (1, nlat)))) + + swaplon = a.lon > b.lon + swaplat = a.lat > c.lat + + loninds = swaplon ? (nlon:-1:1) : (1:1:nlon) + latinds = swaplat ? (nlat:-1:1) : (1:1:nlat) + + gridlonₛ, gridlonₑ = swaplon ? (b.lon, a.lon) : (a.lon, b.lon) + gridlatₛ, gridlatₑ = swaplat ? (c.lat, a.lat) : (a.lat, c.lat) + + lonmin = max(boxmin.lon, gridlonₛ) + latmin = max(boxmin.lat, gridlatₛ) + lonmax = min(boxmax.lon, gridlonₑ) + latmax = min(boxmax.lat, gridlatₑ) + + iₛ = findlast(loninds) do i + p = vertex(grid, (i, 1)) + c = convert(LatLon, coords(p)) + c.lon ≤ lonmin + end + iₑ = findfirst(loninds) do i + p = vertex(grid, (i, 1)) + c = convert(LatLon, coords(p)) + c.lon ≥ lonmax + end + + jₛ = findlast(latinds) do i + p = vertex(grid, (1, i)) + c = convert(LatLon, coords(p)) + c.lat ≤ latmin + end + jₑ = findfirst(latinds) do i + p = vertex(grid, (1, i)) + c = convert(LatLon, coords(p)) + c.lat ≥ latmax + end + + if iₛ == iₑ || jₛ == jₑ + throw(ArgumentError("the passed limits are not valid for the grid")) + end + + iₛ, iₑ = swaplon ? (iₑ, iₛ) : (iₛ, iₑ) + jₛ, jₑ = swaplat ? (jₑ, jₛ) : (jₛ, jₑ) + + CartesianIndex(loninds[iₛ], latinds[jₛ]):CartesianIndex(loninds[iₑ]-1, latinds[jₑ]-1) end # ----------------- @@ -196,98 +212,98 @@ end # ----------------- function _fill!(mask, grid, val, triangle) - v = vertices(triangle) - - # fill edges of triangle - _bresenham!(mask, grid, val, v[1], v[2]) - _bresenham!(mask, grid, val, v[2], v[3]) - _bresenham!(mask, grid, val, v[3], v[1]) - - # fill interior of triangle - j₁ = findfirst(==(val), mask).I[2] - j₂ = findlast(==(val), mask).I[2] - for j in j₁:j₂ - i₁ = findfirst(==(val), @view(mask[:, j])) - i₂ = findlast(==(val), @view(mask[:, j])) - mask[i₁:i₂, j] .= val - end + v = vertices(triangle) + + # fill edges of triangle + _bresenham!(mask, grid, val, v[1], v[2]) + _bresenham!(mask, grid, val, v[2], v[3]) + _bresenham!(mask, grid, val, v[3], v[1]) + + # fill interior of triangle + j₁ = findfirst(==(val), mask).I[2] + j₂ = findlast(==(val), mask).I[2] + for j in j₁:j₂ + i₁ = findfirst(==(val), @view(mask[:, j])) + i₂ = findlast(==(val), @view(mask[:, j])) + mask[i₁:i₂, j] .= val + end end # Bresenham's line algorithm: https://en.wikipedia.org/wiki/Bresenham's_line_algorithm function _bresenham!(mask, grid, val, p₁, p₂) - o = minimum(grid) - s = spacing(grid) - - # integer coordinates - x₁, y₁ = ceil.(Int, (p₁ - o) ./ s) - x₂, y₂ = ceil.(Int, (p₂ - o) ./ s) - - # fix coordinates of points that are on the grid border - xmax, ymax = size(grid) - x₁ = clamp(x₁, 1, xmax) - y₁ = clamp(y₁, 1, ymax) - x₂ = clamp(x₂, 1, xmax) - y₂ = clamp(y₂, 1, ymax) - - if abs(y₂ - y₁) < abs(x₂ - x₁) - if x₁ > x₂ - _bresenhamlow!(mask, val, x₂, y₂, x₁, y₁) - else - _bresenhamlow!(mask, val, x₁, y₁, x₂, y₂) - end - else - if y₁ > y₂ - _bresenhamhigh!(mask, val, x₂, y₂, x₁, y₁) - else - _bresenhamhigh!(mask, val, x₁, y₁, x₂, y₂) - end - end + o = minimum(grid) + s = spacing(grid) + + # integer coordinates + x₁, y₁ = ceil.(Int, (p₁ - o) ./ s) + x₂, y₂ = ceil.(Int, (p₂ - o) ./ s) + + # fix coordinates of points that are on the grid border + xmax, ymax = size(grid) + x₁ = clamp(x₁, 1, xmax) + y₁ = clamp(y₁, 1, ymax) + x₂ = clamp(x₂, 1, xmax) + y₂ = clamp(y₂, 1, ymax) + + if abs(y₂ - y₁) < abs(x₂ - x₁) + if x₁ > x₂ + _bresenhamlow!(mask, val, x₂, y₂, x₁, y₁) + else + _bresenhamlow!(mask, val, x₁, y₁, x₂, y₂) + end + else + if y₁ > y₂ + _bresenhamhigh!(mask, val, x₂, y₂, x₁, y₁) + else + _bresenhamhigh!(mask, val, x₁, y₁, x₂, y₂) + end + end end function _bresenhamlow!(mask, val, x₁, y₁, x₂, y₂) - dx = x₂ - x₁ - dy = y₂ - y₁ - yi = 1 - if dy < 0 - yi = -1 - dy = -dy - end - - D = 2dy - dx - y = y₁ - - for x in x₁:x₂ - mask[x, y] = val - - if D > 0 - y = y + yi - D = D + 2dy - 2dx - else - D = D + 2dy - end - end + dx = x₂ - x₁ + dy = y₂ - y₁ + yi = 1 + if dy < 0 + yi = -1 + dy = -dy + end + + D = 2dy - dx + y = y₁ + + for x in x₁:x₂ + mask[x, y] = val + + if D > 0 + y = y + yi + D = D + 2dy - 2dx + else + D = D + 2dy + end + end end function _bresenhamhigh!(mask, val, x₁, y₁, x₂, y₂) - dx = x₂ - x₁ - dy = y₂ - y₁ - xi = 1 - if dx < 0 - xi = -1 - dx = -dx - end - - D = 2dx - dy - x = x₁ - - for y in y₁:y₂ - mask[x, y] = val - - if D > 0 - x = x + xi - D = D + 2dx - 2dy - else - D = D + 2dx - end - end + dx = x₂ - x₁ + dy = y₂ - y₁ + xi = 1 + if dx < 0 + xi = -1 + dx = -dx + end + + D = 2dx - dy + x = x₁ + + for y in y₁:y₂ + mask[x, y] = val + + if D > 0 + x = x + xi + D = D + 2dx - 2dy + else + D = D + 2dx + end + end end diff --git a/test/indices.jl b/test/indices.jl index 519eb4d8e..ba4e6fad0 100644 --- a/test/indices.jl +++ b/test/indices.jl @@ -1,216 +1,206 @@ @testitem "Indices" setup = [Setup] begin - g = cartgrid(10, 10) - b = Box(cart(1, 1), cart(5, 5)) - v = view(g, b) - @test v == CartesianGrid(cart(0, 0), cart(6, 6), dims=(6, 6)) - - p = PointSet(vertices(g)) - v = view(p, b) - @test centroid(v, 1) == cart(1, 1) - @test centroid(v, nelements(v)) == cart(5, 5) - - # boxes - g = cartgrid(10, 10) - p = PointSet(vertices(g)) - b = Ball(cart(0, 0), T(2)) - v = view(g, b) - @test nelements(v) == 4 - @test v[1] == g[1] - v = view(p, b) - @test nelements(v) == 6 - @test to.(v) == vector.([(0, 0), (1, 0), (2, 0), (0, 1), (1, 1), (0, 2)]) - - b1 = Box(cart(0.6, 0.7), cart(1.0, 1.0)) - b2 = Box(cart(0.6, 0.7), cart(1.2, 1.2)) - b3 = Box(cart(0.0, 0.0), cart(0.6, 0.3)) - b4 = Box(cart(-0.2, -0.2), cart(0.6, 0.3)) - b5 = Box(cart(0.25, 0.15), cart(0.95, 0.85)) - b6 = Box(cart(0.35, 0.25), cart(0.85, 0.75)) - b7 = Box(cart(0.05, 0.05), cart(0.65, 0.35)) - b8 = Box(cart(0.55, 0.65), cart(0.95, 0.95)) - x = range(zero(T), stop=one(T), length=6) - y = T[0.0, 0.1, 0.3, 0.7, 0.9, 1.0] - g = RectilinearGrid(x, y) - linds = LinearIndices(size(g)) - @test issetequal(indices(g, b1), [linds[4, 4], linds[5, 4], linds[4, 5], linds[5, 5]]) - @test issetequal(indices(g, b2), [linds[4, 4], linds[5, 4], linds[4, 5], linds[5, 5]]) - @test issetequal(indices(g, b3), [linds[1, 1], linds[2, 1], linds[3, 1], linds[1, 2], linds[2, 2], linds[3, 2]]) - @test issetequal(indices(g, b4), [linds[1, 1], linds[2, 1], linds[3, 1], linds[1, 2], linds[2, 2], linds[3, 2]]) - @test issetequal( - indices(g, b5), - [ - linds[2, 2], - linds[3, 2], - linds[4, 2], - linds[5, 2], - linds[2, 3], - linds[3, 3], - linds[4, 3], - linds[5, 3], - linds[2, 4], - linds[3, 4], - linds[4, 4], - linds[5, 4] - ] - ) - @test issetequal( - indices(g, b6), - [ - linds[2, 2], - linds[3, 2], - linds[4, 2], - linds[5, 2], - linds[2, 3], - linds[3, 3], - linds[4, 3], - linds[5, 3], - linds[2, 4], - linds[3, 4], - linds[4, 4], - linds[5, 4] - ] - ) - @test issetequal( - indices(g, b7), - [ - linds[1, 1], - linds[2, 1], - linds[3, 1], - linds[4, 1], - linds[1, 2], - linds[2, 2], - linds[3, 2], - linds[4, 2], - linds[1, 3], - linds[2, 3], - linds[3, 3], - linds[4, 3] - ] - ) - @test issetequal( - indices(g, b8), - [ - linds[3, 3], - linds[4, 3], - linds[5, 3], - linds[3, 4], - linds[4, 4], - linds[5, 4], - linds[3, 5], - linds[4, 5], - linds[5, 5] - ] - ) - - # convex polygons - tri = Triangle(cart(5, 7), cart(10, 12), cart(15, 7)) - pent = Pentagon(cart(6, 1), cart(2, 10), cart(10, 16), cart(18, 10), cart(14, 1)) - - grid = cartgrid(20, 20) - linds = LinearIndices(size(grid)) - @test linds[10, 10] ∈ indices(grid, tri) - @test linds[10, 6] ∈ indices(grid, pent) - - grid = CartesianGrid(cart(-2, -2), cart(20, 20), T.((0.5, 1.5))) - linds = LinearIndices(size(grid)) - @test linds[21, 7] ∈ indices(grid, tri) - @test linds[21, 4] ∈ indices(grid, pent) - - grid = CartesianGrid(cart(-100, -100), cart(20, 20), T.((2, 2))) - linds = LinearIndices(size(grid)) - @test linds[57, 54] ∈ indices(grid, tri) - @test linds[55, 53] ∈ indices(grid, pent) - - # non-convex polygons - poly1 = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) - poly2 = PolyArea([pointify(pent), pointify(tri)]) - - grid = cartgrid(20, 20) - linds = LinearIndices(size(grid)) - @test linds[12, 6] ∈ indices(grid, poly1) - @test linds[10, 3] ∈ indices(grid, poly2) - - grid = CartesianGrid(cart(-2, -2), cart(20, 20), T.((0.5, 1.5))) - linds = LinearIndices(size(grid)) - @test linds[22, 6] ∈ indices(grid, poly1) - @test linds[17, 4] ∈ indices(grid, poly2) - - grid = CartesianGrid(cart(-100, -100), cart(20, 20), T.((2, 2))) - linds = LinearIndices(size(grid)) - @test linds[57, 54] ∈ indices(grid, poly1) - @test linds[55, 53] ∈ indices(grid, poly2) - - # rotate - poly1 = poly1 |> Rotate(Angle2d(T(π / 2))) - poly2 = poly2 |> Rotate(Angle2d(T(π / 2))) - - grid = CartesianGrid(cart(-20, 0), cart(0, 20), T.((1, 1))) - linds = LinearIndices(size(grid)) - @test linds[12, 12] ∈ indices(grid, poly1) - @test linds[16, 11] ∈ indices(grid, poly2) - - grid = CartesianGrid(cart(-22, -2), cart(0, 20), T.((0.5, 1.5))) - linds = LinearIndices(size(grid)) - @test linds[26, 8] ∈ indices(grid, poly1) - @test linds[36, 9] ∈ indices(grid, poly2) - - grid = CartesianGrid(cart(-100, -100), cart(20, 20), T.((2, 2))) - linds = LinearIndices(size(grid)) - @test linds[46, 57] ∈ indices(grid, poly1) - @test linds[48, 55] ∈ indices(grid, poly2) - - # multi - multi = Multi([tri, pent]) - grid = cartgrid(20, 20) - linds = LinearIndices(size(grid)) - @test linds[10, 10] ∈ indices(grid, multi) - @test linds[10, 6] ∈ indices(grid, multi) - - # clipping - tri = Triangle(cart(-4, 10), cart(5, 19), cart(5, 1)) - grid = cartgrid(20, 20) - linds = LinearIndices(size(grid)) - @test linds[3, 10] ∈ indices(grid, tri) - - # out of grid - tri = Triangle(cart(-12, 8), cart(-8, 14), cart(-4, 8)) - grid = cartgrid(20, 20) - @test isempty(indices(grid, tri)) - - # chain - seg = Segment(cart(2, 12), cart(16, 18)) - rope = Rope(cart(8, 1), cart(5, 9), cart(9, 13), cart(17, 10)) - ring = Ring(cart(8, 1), cart(5, 9), cart(9, 13), cart(17, 10)) - grid = cartgrid(20, 20) - linds = LinearIndices(size(grid)) - @test linds[9, 15] ∈ indices(grid, seg) - @test linds[7, 11] ∈ indices(grid, rope) - @test linds[12, 5] ∈ indices(grid, ring) - - # points - p1 = cart(0, 0) - p2 = cart(0.5, 0.5) - p3 = cart(1, 1) - p4 = cart(2, 2) - p5 = cart(10, 10) - p6 = cart(11, 11) - grid = cartgrid(10, 10) - linds = LinearIndices(size(grid)) - @test linds[1, 1] == only(indices(grid, p1)) - @test linds[1, 1] == only(indices(grid, p2)) - @test linds[1, 1] == only(indices(grid, p3)) - @test linds[2, 2] == only(indices(grid, p4)) - @test linds[10, 10] == only(indices(grid, p5)) - @test isempty(indices(grid, p6)) - - # transformed grid - grid = cartgrid(20, 20) - poly = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) - inds = indices(grid, poly) - linds = LinearIndices(size(grid)) - @test linds[12, 6] ∈ indices(grid, poly) - trans = Rotate(Angle2d(T(π / 2))) → Morphological(c -> Cartesian(c.x + c.y, c.y, c.x - c.y)) - tgrid = trans(grid) - tpoly = trans(poly) - @test indices(tgrid, tpoly) == indices(grid, poly) + g = cartgrid(10, 10) + b = Box(cart(1, 1), cart(5, 5)) + v = view(g, b) + @test v == CartesianGrid(cart(0, 0), cart(6, 6), dims=(6, 6)) + + p = PointSet(vertices(g)) + v = view(p, b) + @test centroid(v, 1) == cart(1, 1) + @test centroid(v, nelements(v)) == cart(5, 5) + + # boxes + g = cartgrid(10, 10) + p = PointSet(vertices(g)) + b = Ball(cart(0, 0), T(2)) + v = view(g, b) + @test nelements(v) == 4 + @test v[1] == g[1] + v = view(p, b) + @test nelements(v) == 6 + @test to.(v) == vector.([(0, 0), (1, 0), (2, 0), (0, 1), (1, 1), (0, 2)]) + + b1 = Box(cart(0.6, 0.7), cart(1.0, 1.0)) + b2 = Box(cart(0.6, 0.7), cart(1.2, 1.2)) + b3 = Box(cart(0.0, 0.0), cart(0.6, 0.3)) + b4 = Box(cart(-0.2, -0.2), cart(0.6, 0.3)) + b5 = Box(cart(0.25, 0.15), cart(0.95, 0.85)) + b6 = Box(cart(0.35, 0.25), cart(0.85, 0.75)) + b7 = Box(cart(0.05, 0.05), cart(0.65, 0.35)) + b8 = Box(cart(0.55, 0.65), cart(0.95, 0.95)) + x = range(zero(T), stop=one(T), length=6) + y = T[0.0, 0.1, 0.3, 0.7, 0.9, 1.0] + g = RectilinearGrid(x, y) + linds = LinearIndices(size(g)) + @test issetequal(indices(g, b1), [linds[4, 4], linds[5, 4], linds[4, 5], linds[5, 5]]) + @test issetequal(indices(g, b2), [linds[4, 4], linds[5, 4], linds[4, 5], linds[5, 5]]) + @test issetequal(indices(g, b3), [linds[1, 1], linds[2, 1], linds[3, 1], linds[1, 2], linds[2, 2], linds[3, 2]]) + @test issetequal(indices(g, b4), [linds[1, 1], linds[2, 1], linds[3, 1], linds[1, 2], linds[2, 2], linds[3, 2]]) + @test issetequal( + indices(g, b5), + [ + linds[2, 2], + linds[3, 2], + linds[4, 2], + linds[5, 2], + linds[2, 3], + linds[3, 3], + linds[4, 3], + linds[5, 3], + linds[2, 4], + linds[3, 4], + linds[4, 4], + linds[5, 4], + ], + ) + @test issetequal( + indices(g, b6), + [ + linds[2, 2], + linds[3, 2], + linds[4, 2], + linds[5, 2], + linds[2, 3], + linds[3, 3], + linds[4, 3], + linds[5, 3], + linds[2, 4], + linds[3, 4], + linds[4, 4], + linds[5, 4], + ], + ) + @test issetequal( + indices(g, b7), + [ + linds[1, 1], + linds[2, 1], + linds[3, 1], + linds[4, 1], + linds[1, 2], + linds[2, 2], + linds[3, 2], + linds[4, 2], + linds[1, 3], + linds[2, 3], + linds[3, 3], + linds[4, 3], + ], + ) + @test issetequal( + indices(g, b8), + [ + linds[3, 3], + linds[4, 3], + linds[5, 3], + linds[3, 4], + linds[4, 4], + linds[5, 4], + linds[3, 5], + linds[4, 5], + linds[5, 5], + ], + ) + + # convex polygons + tri = Triangle(cart(5, 7), cart(10, 12), cart(15, 7)) + pent = Pentagon(cart(6, 1), cart(2, 10), cart(10, 16), cart(18, 10), cart(14, 1)) + + grid = cartgrid(20, 20) + linds = LinearIndices(size(grid)) + @test linds[10, 10] ∈ indices(grid, tri) + @test linds[10, 6] ∈ indices(grid, pent) + + grid = CartesianGrid(cart(-2, -2), cart(20, 20), T.((0.5, 1.5))) + linds = LinearIndices(size(grid)) + @test linds[21, 7] ∈ indices(grid, tri) + @test linds[21, 4] ∈ indices(grid, pent) + + grid = CartesianGrid(cart(-100, -100), cart(20, 20), T.((2, 2))) + linds = LinearIndices(size(grid)) + @test linds[57, 54] ∈ indices(grid, tri) + @test linds[55, 53] ∈ indices(grid, pent) + + # non-convex polygons + poly1 = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) + poly2 = PolyArea([pointify(pent), pointify(tri)]) + + grid = cartgrid(20, 20) + linds = LinearIndices(size(grid)) + @test linds[12, 6] ∈ indices(grid, poly1) + @test linds[10, 3] ∈ indices(grid, poly2) + + grid = CartesianGrid(cart(-2, -2), cart(20, 20), T.((0.5, 1.5))) + linds = LinearIndices(size(grid)) + @test linds[22, 6] ∈ indices(grid, poly1) + @test linds[17, 4] ∈ indices(grid, poly2) + + grid = CartesianGrid(cart(-100, -100), cart(20, 20), T.((2, 2))) + linds = LinearIndices(size(grid)) + @test linds[57, 54] ∈ indices(grid, poly1) + @test linds[55, 53] ∈ indices(grid, poly2) + + # rotate + poly1 = poly1 |> Rotate(Angle2d(T(π / 2))) + poly2 = poly2 |> Rotate(Angle2d(T(π / 2))) + + grid = CartesianGrid(cart(-20, 0), cart(0, 20), T.((1, 1))) + linds = LinearIndices(size(grid)) + @test linds[12, 12] ∈ indices(grid, poly1) + @test linds[16, 11] ∈ indices(grid, poly2) + + grid = CartesianGrid(cart(-22, -2), cart(0, 20), T.((0.5, 1.5))) + linds = LinearIndices(size(grid)) + @test linds[26, 8] ∈ indices(grid, poly1) + @test linds[36, 9] ∈ indices(grid, poly2) + + grid = CartesianGrid(cart(-100, -100), cart(20, 20), T.((2, 2))) + linds = LinearIndices(size(grid)) + @test linds[46, 57] ∈ indices(grid, poly1) + @test linds[48, 55] ∈ indices(grid, poly2) + + # multi + multi = Multi([tri, pent]) + grid = cartgrid(20, 20) + linds = LinearIndices(size(grid)) + @test linds[10, 10] ∈ indices(grid, multi) + @test linds[10, 6] ∈ indices(grid, multi) + + # clipping + tri = Triangle(cart(-4, 10), cart(5, 19), cart(5, 1)) + grid = cartgrid(20, 20) + linds = LinearIndices(size(grid)) + @test linds[3, 10] ∈ indices(grid, tri) + + # out of grid + tri = Triangle(cart(-12, 8), cart(-8, 14), cart(-4, 8)) + grid = cartgrid(20, 20) + @test isempty(indices(grid, tri)) + + # chain + seg = Segment(cart(2, 12), cart(16, 18)) + rope = Rope(cart(8, 1), cart(5, 9), cart(9, 13), cart(17, 10)) + ring = Ring(cart(8, 1), cart(5, 9), cart(9, 13), cart(17, 10)) + grid = cartgrid(20, 20) + linds = LinearIndices(size(grid)) + @test linds[9, 15] ∈ indices(grid, seg) + @test linds[7, 11] ∈ indices(grid, rope) + @test linds[12, 5] ∈ indices(grid, ring) + + # points + p1 = cart(0, 0) + p2 = cart(0.5, 0.5) + p3 = cart(1, 1) + p4 = cart(2, 2) + p5 = cart(10, 10) + p6 = cart(11, 11) + grid = cartgrid(10, 10) + linds = LinearIndices(size(grid)) + @test linds[1, 1] == only(indices(grid, p1)) + @test linds[1, 1] == only(indices(grid, p2)) + @test linds[1, 1] == only(indices(grid, p3)) + @test linds[2, 2] == only(indices(grid, p4)) + @test linds[10, 10] == only(indices(grid, p5)) + @test isempty(indices(grid, p6)) + end From 47b064bd53a36e9e4f5092e811c80c3df0767d35 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva Date: Sat, 1 Nov 2025 20:45:25 -0300 Subject: [PATCH 13/27] Add fast path for indices on non-invertible TransformedGrid --- src/indices.jl | 14 ++++++++++++++ test/indices.jl | 14 ++++++++++++++ 2 files changed, 28 insertions(+) diff --git a/src/indices.jl b/src/indices.jl index 722379bd7..53ba867f1 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -74,6 +74,20 @@ function indices(grid::OrthoRectilinearGrid, box::Box) LinearIndices(size(grid))[range] |> vec end +function indices(grid::TransformedGrid, geometry::Geometry) + g = grid.mesh + t = grid.transform + + # find the invertible part (e.g., Rotate, Affine) + invertible_part = only(x for x in t if isinvertible(x)) + + # un-transform the geometry using only that part + p = revert(invertible_part, geometry, nothing) + + # find indices on the non-deformed mesh + indices(g, p) +end + # ---------------- # CARTESIAN RANGE # ---------------- diff --git a/test/indices.jl b/test/indices.jl index 6d5f2ccc4..82caebfc5 100644 --- a/test/indices.jl +++ b/test/indices.jl @@ -202,4 +202,18 @@ @test linds[2, 2] == only(indices(grid, p4)) @test linds[10, 10] == only(indices(grid, p5)) @test isempty(indices(grid, p6)) + + # transformed grid + g = cartgrid(20, 20) + p_orig = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) + expected_inds = indices(g, p_orig) + linds = LinearIndices(size(g)) + @test linds[12, 6] ∈ expected_inds + t_aff = Rotate(Angle2d(T(π / 2))) + t_morph = Morphological(c -> Cartesian(c.x + c.y, c.y, c.x - c.y)) + t_seq = t_aff → t_morph + grid = TransformedGrid(g, t_seq) + poly = t_aff(p_orig) + actual_inds = indices(grid, poly) + @test actual_inds == expected_inds end From 882db7a0bded69f7f8a6fcb9a953b0868b5cd514 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:15:00 -0300 Subject: [PATCH 14/27] Update src/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- src/indices.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/indices.jl b/src/indices.jl index 53ba867f1..3d32053d2 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -75,8 +75,8 @@ function indices(grid::OrthoRectilinearGrid, box::Box) end function indices(grid::TransformedGrid, geometry::Geometry) - g = grid.mesh - t = grid.transform + g = parent(grid) + t = transform(grid) # find the invertible part (e.g., Rotate, Affine) invertible_part = only(x for x in t if isinvertible(x)) From 4201d1acf2d89e8c681f726a2ad51f37e799d314 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:15:26 -0300 Subject: [PATCH 15/27] Update src/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- src/indices.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/indices.jl b/src/indices.jl index 3d32053d2..f9b929fdf 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -79,7 +79,7 @@ function indices(grid::TransformedGrid, geometry::Geometry) t = transform(grid) # find the invertible part (e.g., Rotate, Affine) - invertible_part = only(x for x in t if isinvertible(x)) + r = reduce(→, reverse(filter(isrevertible, t))) # un-transform the geometry using only that part p = revert(invertible_part, geometry, nothing) From 14f0e65ac992005d53f4869e0f540bff195253f8 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:15:36 -0300 Subject: [PATCH 16/27] Update src/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- src/indices.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/indices.jl b/src/indices.jl index f9b929fdf..d299f5485 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -78,7 +78,7 @@ function indices(grid::TransformedGrid, geometry::Geometry) g = parent(grid) t = transform(grid) - # find the invertible part (e.g., Rotate, Affine) + # construct reverse transform from revertible steps r = reduce(→, reverse(filter(isrevertible, t))) # un-transform the geometry using only that part From 04b9b9d5a4ea21763ee08fb49fa9aab79d4b57ca Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:16:04 -0300 Subject: [PATCH 17/27] Update src/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- src/indices.jl | 7 ++----- 1 file changed, 2 insertions(+), 5 deletions(-) diff --git a/src/indices.jl b/src/indices.jl index d299f5485..f46e962eb 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -81,11 +81,8 @@ function indices(grid::TransformedGrid, geometry::Geometry) # construct reverse transform from revertible steps r = reduce(→, reverse(filter(isrevertible, t))) - # un-transform the geometry using only that part - p = revert(invertible_part, geometry, nothing) - - # find indices on the non-deformed mesh - indices(g, p) + # find indices before transformation + indices(g, r(geometry)) end # ---------------- From c49779f984954f12382f73d7630aab61f9a47ed3 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva Date: Sat, 1 Nov 2025 21:31:03 -0300 Subject: [PATCH 18/27] fix test names --- test/indices.jl | 20 ++++++++++---------- 1 file changed, 10 insertions(+), 10 deletions(-) diff --git a/test/indices.jl b/test/indices.jl index 82caebfc5..44f7da98c 100644 --- a/test/indices.jl +++ b/test/indices.jl @@ -205,15 +205,15 @@ # transformed grid g = cartgrid(20, 20) - p_orig = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) - expected_inds = indices(g, p_orig) + poly1 = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) + expectedinds = indices(g, poly1) linds = LinearIndices(size(g)) - @test linds[12, 6] ∈ expected_inds - t_aff = Rotate(Angle2d(T(π / 2))) - t_morph = Morphological(c -> Cartesian(c.x + c.y, c.y, c.x - c.y)) - t_seq = t_aff → t_morph - grid = TransformedGrid(g, t_seq) - poly = t_aff(p_orig) - actual_inds = indices(grid, poly) - @test actual_inds == expected_inds + @test linds[12, 6] ∈ expectedinds + taff = Rotate(Angle2d(T(π / 2))) + tmorph = Morphological(c -> Cartesian(c.x + c.y, c.y, c.x - c.y)) + tseq = taff → tmorph + grid = TransformedGrid(g, tseq) + poly2 = taff(poly1) + actualinds = indices(grid, poly2) + @test actualinds == expectedinds end From 8cb03f5948de5d91e6ab7c5defc52d86c90fcfcf Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:30:30 -0300 Subject: [PATCH 19/27] Update src/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- src/indices.jl | 9 +++------ 1 file changed, 3 insertions(+), 6 deletions(-) diff --git a/src/indices.jl b/src/indices.jl index f46e962eb..0724a62ab 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -75,14 +75,11 @@ function indices(grid::OrthoRectilinearGrid, box::Box) end function indices(grid::TransformedGrid, geometry::Geometry) - g = parent(grid) - t = transform(grid) - # construct reverse transform from revertible steps - r = reduce(→, reverse(filter(isrevertible, t))) + revtrans = reduce(→, reverse(filter(isrevertible, transform(t)))) - # find indices before transformation - indices(g, r(geometry)) + # find indices in non-transformed space + indices(parent(grid), revtrans(geometry)) end # ---------------- From 2cdb2672999ce0ff46a7fa91a6729b393f9d4398 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva Date: Sat, 1 Nov 2025 21:53:21 -0300 Subject: [PATCH 20/27] fix indices transform --- src/indices.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/indices.jl b/src/indices.jl index 0724a62ab..6bd058cb6 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -74,9 +74,9 @@ function indices(grid::OrthoRectilinearGrid, box::Box) LinearIndices(size(grid))[range] |> vec end -function indices(grid::TransformedGrid, geometry::Geometry) +function indice(grid::TransformedGrid, geometry::Geometry) # construct reverse transform from revertible steps - revtrans = reduce(→, reverse(filter(isrevertible, transform(t)))) + revtrans = reduce(→, reverse(filter(isrevertible, transform(grid)))) # find indices in non-transformed space indices(parent(grid), revtrans(geometry)) From 9894eeac7cb6046cc16e7c330ccc88f345b71f12 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva Date: Sat, 1 Nov 2025 21:53:59 -0300 Subject: [PATCH 21/27] fix transform in indices --- src/indices.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/indices.jl b/src/indices.jl index 6bd058cb6..0b468211d 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -74,7 +74,7 @@ function indices(grid::OrthoRectilinearGrid, box::Box) LinearIndices(size(grid))[range] |> vec end -function indice(grid::TransformedGrid, geometry::Geometry) +function indices(grid::TransformedGrid, geometry::Geometry) # construct reverse transform from revertible steps revtrans = reduce(→, reverse(filter(isrevertible, transform(grid)))) From e5fcb803e6a6cc284ea20e9b59d6afcf7ee32c91 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:51:14 -0300 Subject: [PATCH 22/27] Update test/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- test/indices.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/indices.jl b/test/indices.jl index 44f7da98c..4e78faa9f 100644 --- a/test/indices.jl +++ b/test/indices.jl @@ -204,7 +204,7 @@ @test isempty(indices(grid, p6)) # transformed grid - g = cartgrid(20, 20) + grid = cartgrid(20, 20) poly1 = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) expectedinds = indices(g, poly1) linds = LinearIndices(size(g)) From 52144ea6f6a6c65615115c959c17181f9088b470 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Sat, 1 Nov 2025 21:51:21 -0300 Subject: [PATCH 23/27] Update test/indices.jl MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-authored-by: Júlio Hoffimann --- test/indices.jl | 19 ++++++++----------- 1 file changed, 8 insertions(+), 11 deletions(-) diff --git a/test/indices.jl b/test/indices.jl index 4e78faa9f..519eb4d8e 100644 --- a/test/indices.jl +++ b/test/indices.jl @@ -205,15 +205,12 @@ # transformed grid grid = cartgrid(20, 20) - poly1 = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) - expectedinds = indices(g, poly1) - linds = LinearIndices(size(g)) - @test linds[12, 6] ∈ expectedinds - taff = Rotate(Angle2d(T(π / 2))) - tmorph = Morphological(c -> Cartesian(c.x + c.y, c.y, c.x - c.y)) - tseq = taff → tmorph - grid = TransformedGrid(g, tseq) - poly2 = taff(poly1) - actualinds = indices(grid, poly2) - @test actualinds == expectedinds + poly = PolyArea(cart.([(3, 3), (9, 9), (3, 15), (17, 15), (17, 3)])) + inds = indices(grid, poly) + linds = LinearIndices(size(grid)) + @test linds[12, 6] ∈ indices(grid, poly) + trans = Rotate(Angle2d(T(π / 2))) → Morphological(c -> Cartesian(c.x + c.y, c.y, c.x - c.y)) + tgrid = trans(grid) + tpoly = trans(poly) + @test indices(tgrid, tpoly) == indices(grid, poly) end From 5627ca48769a95886f473eaa7225f00ffa2e0cbd Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva Date: Wed, 19 Nov 2025 17:45:02 -0300 Subject: [PATCH 24/27] _reversetransform and indices --- src/Meshes.jl | 936 +++++++++++++++++++++++++------------------------ src/indices.jl | 420 +++++++++++----------- 2 files changed, 687 insertions(+), 669 deletions(-) diff --git a/src/Meshes.jl b/src/Meshes.jl index 0dc6a8ff8..fa9b22d27 100644 --- a/src/Meshes.jl +++ b/src/Meshes.jl @@ -32,7 +32,7 @@ using ScopedValues: ScopedValue using Base.Cartesian: @nloops, @nref, @ntuple using Base: @propagate_inbounds -import Random +using Random: Random import Base: sort import Base: ==, ! import Base: +, -, * @@ -41,7 +41,7 @@ import Distances: evaluate import NearestNeighbors: MinkowskiMetric # Transforms API -import TransformsBase: Transform, → +import TransformsBase: Transform, →, SequentialTransform, Identity import TransformsBase: isrevertible, isinvertible import TransformsBase: apply, revert, reapply, inverse import TransformsBase: parameters, preprocess @@ -78,6 +78,9 @@ include("domains.jl") # utilities include("utils.jl") +# transforms +include("transforms.jl") + # domain indices include("indices.jl") @@ -120,8 +123,7 @@ include("discretization.jl") include("refinement.jl") include("coarsening.jl") -# transforms -include("transforms.jl") + # miscellaneous include("rand.jl") @@ -134,471 +136,471 @@ include("projecting.jl") include("viz.jl") function __init__() - # register error hint for visualization functions - # since this is a recurring issue for new users - Base.Experimental.register_error_hint(MethodError) do io, exc, argtypes, kwargs - if exc.f == viz || exc.f == viz! - if isnothing(Base.get_extension(Meshes, :MeshesMakieExt)) - print( - io, - """ - - Did you import a Makie.jl backend (e.g., GLMakie.jl, CairoMakie.jl) for visualization? - - """ - ) - printstyled( - io, - """ - julia> using Meshes - julia> import GLMakie # or CairoMakie, WGLMakie, etc. - """, - color=:cyan, - bold=true - ) - end - end - end + # register error hint for visualization functions + # since this is a recurring issue for new users + Base.Experimental.register_error_hint(MethodError) do io, exc, argtypes, kwargs + if exc.f == viz || exc.f == viz! + if isnothing(Base.get_extension(Meshes, :MeshesMakieExt)) + print( + io, + """ + + Did you import a Makie.jl backend (e.g., GLMakie.jl, CairoMakie.jl) for visualization? + + """, + ) + printstyled( + io, + """ + julia> using Meshes + julia> import GLMakie # or CairoMakie, WGLMakie, etc. + """, + color=:cyan, + bold=true, + ) + end + end + end end export - # vectors - Vec, - ∠, - ⋅, - ×, - - # manifolds - 𝔼, - 🌐, - - # geometries - Geometry, - embeddim, - paramdim, - crs, - manifold, - - # primitives - Primitive, - Point, - Ray, - Line, - BezierCurve, - ParametrizedCurve, - Plane, - Box, - Ball, - Sphere, - Ellipsoid, - Disk, - Circle, - Cylinder, - CylinderSurface, - ParaboloidSurface, - Cone, - ConeSurface, - Frustum, - FrustumSurface, - Torus, - controls, - ncontrols, - degree, - Horner, - DeCasteljau, - coords, - to, - center, - radius, - radii, - plane, - bottom, - top, - axis, - base, - apex, - isright, - sides, - diagonal, - focallength, - direction, - - # polytopes - Polytope, - Chain, - Segment, - Rope, - Ring, - Polygon, - Ngon, - Triangle, - Quadrangle, - Pentagon, - Hexagon, - Heptagon, - Octagon, - Nonagon, - Decagon, - PolyArea, - Polyhedron, - Tetrahedron, - Hexahedron, - Pyramid, - Wedge, - vertex, - vertices, - nvertices, - eachvertex, - rings, - segments, - angles, - innerangles, - normal, - base, - apex, - ≗, - - # multi-geometries - Multi, - MultiPoint, - MultiSegment, - MultiRope, - MultiRing, - MultiPolygon, - MultiPolyhedron, - - # transformed geometry - TransformedGeometry, - - # connectivities - Connectivity, - paramdim, - indices, - connect, - materialize, - pltype, - - # topologies - Topology, - vertex, - vertices, - nvertices, - element, - elements, - nelements, - facet, - facets, - nfacets, - faces, - nfaces, - GridTopology, - HalfEdgeTopology, - HalfEdge, - SimpleTopology, - elem2cart, - cart2elem, - corner2cart, - cart2corner, - elem2corner, - corner2elem, - elementtype, - facettype, - half4elem, - half4vert, - half4edge, - half4pair, - edge4pair, - connec4elem, - - # topological relations - TopologicalRelation, - Boundary, - Coboundary, - Adjacency, - - # domain traits - Domain, - SubDomain, - TransformedDomain, - embeddim, - paramdim, - crs, - manifold, - element, - nelements, - - # sets - GeometrySet, - PointSet, - - # meshes - Mesh, - SimpleMesh, - TransformedMesh, - Grid, - SubGrid, - RegularGrid, - CartesianGrid, - RectilinearGrid, - StructuredGrid, - TransformedGrid, - vertex, - vertices, - nvertices, - eachvertex, - element, - elements, - nelements, - facet, - facets, - nfacets, - topology, - topoconvert, - spacing, - offset, - - # trajectories - CylindricalTrajectory, - - # indices - indices, - - # partitions - Partition, - indices, - metadata, - - # partitioning - PartitionMethod, - IndexPredicatePartitionMethod, - PointPredicatePartitionMethod, - UniformPartition, - FractionPartition, - BlockPartition, - BisectPointPartition, - BisectFractionPartition, - BallPartition, - PlanePartition, - DirectionPartition, - IndexPredicatePartition, - PointPredicatePartition, - ProductPartition, - HierarchicalPartition, - partition, - split, - - # sorting - SortingMethod, - DirectionSort, - sortinds, - - # traversing - Path, - LinearPath, - RandomPath, - ShiftedPath, - SourcePath, - MultiGridPath, - traverse, - - # neighborhoods - Neighborhood, - MetricBall, - hasequalradii, - radii, - rotation, - metric, - radius, - - # neighborhood search - NeighborSearchMethod, - BoundedNeighborSearchMethod, - BallSearch, - KNearestSearch, - KBallSearch, - search!, - searchdists!, - search, - searchdists, - maxneighbors, - - # predicates - isparametrized, - iscurve, - issurface, - issolid, - isperiodic, - issimplex, - isclosed, - isconvex, - issimple, - hasholes, - intersects, - iscollinear, - iscoplanar, - ≺, - ≻, - ⪯, - ⪰, - - # centroids - centroid, - - # measures - measure, - area, - volume, - perimeter, - - # boundary - boundary, - - # winding number - winding, - - # sideof - sideof, - SideType, - IN, - OUT, - ON, - LEFT, - RIGHT, - - # orientation - orientation, - OrientationType, - CW, - CCW, - - # clipping - ClippingMethod, - SutherlandHodgmanClipping, - clip, - - # intersections - IntersectionType, - Crossing, - CornerCrossing, - EdgeCrossing, - Touching, - CornerTouching, - EdgeTouching, - Overlapping, - PosOverlapping, - NegOverlapping, - NotIntersecting, - Intersecting, - - # intersecting - Intersection, - intersection, - type, - - # simplification - SimplificationMethod, - SelingerSimplification, - DouglasPeuckerSimplification, - MinMaxSimplification, - simplify, - - # bounding boxes - boundingbox, - - # hulls - HullMethod, - GrahamScan, - JarvisMarch, - hull, - convexhull, - - # sampling - SamplingMethod, - DiscreteSamplingMethod, - ContinuousSamplingMethod, - UniformSampling, - WeightedSampling, - BallSampling, - BlockSampling, - RegularSampling, - HomogeneousSampling, - MinDistanceSampling, - FibonacciSampling, - sampleinds, - sample, - - # pointification - pointify, - - # tesselation - TesselationMethod, - DelaunayTesselation, - VoronoiTesselation, - tesselate, - - # discretization - DiscretizationMethod, - BoundaryTriangulationMethod, - FanTriangulation, - DehnTriangulation, - HeldTriangulation, - DelaunayTriangulation, - ManualSimplexification, - RegularDiscretization, - MaxLengthDiscretization, - discretize, - discretizewithin, - simplexify, - - # refinement - RefinementMethod, - TriRefinement, - QuadRefinement, - RegularRefinement, - CatmullClarkRefinement, - TriSubdivision, - refine, - - # coarsening - CoarseningMethod, - RegularCoarsening, - coarsen, - - # transforms - GeometricTransform, - CoordinateTransform, - Rotate, - Translate, - Scale, - Affine, - Stretch, - StdCoords, - Proj, - Morphological, - LengthUnit, - ValidCoords, - Shadow, - Slice, - Repair, - Bridge, - LambdaMuSmoothing, - LaplaceSmoothing, - TaubinSmoothing, - isaffine, - isinvertible, - inverse, - - # miscellaneous - signarea, - householderbasis, - supportfun, - laplacematrix, - measurematrix, - adjacencymatrix, - atol, - - # visualization - viz, - viz! + # vectors + Vec, + ∠, + ⋅, + ×, + + # manifolds + 𝔼, + 🌐, + + # geometries + Geometry, + embeddim, + paramdim, + crs, + manifold, + + # primitives + Primitive, + Point, + Ray, + Line, + BezierCurve, + ParametrizedCurve, + Plane, + Box, + Ball, + Sphere, + Ellipsoid, + Disk, + Circle, + Cylinder, + CylinderSurface, + ParaboloidSurface, + Cone, + ConeSurface, + Frustum, + FrustumSurface, + Torus, + controls, + ncontrols, + degree, + Horner, + DeCasteljau, + coords, + to, + center, + radius, + radii, + plane, + bottom, + top, + axis, + base, + apex, + isright, + sides, + diagonal, + focallength, + direction, + + # polytopes + Polytope, + Chain, + Segment, + Rope, + Ring, + Polygon, + Ngon, + Triangle, + Quadrangle, + Pentagon, + Hexagon, + Heptagon, + Octagon, + Nonagon, + Decagon, + PolyArea, + Polyhedron, + Tetrahedron, + Hexahedron, + Pyramid, + Wedge, + vertex, + vertices, + nvertices, + eachvertex, + rings, + segments, + angles, + innerangles, + normal, + base, + apex, + ≗, + + # multi-geometries + Multi, + MultiPoint, + MultiSegment, + MultiRope, + MultiRing, + MultiPolygon, + MultiPolyhedron, + + # transformed geometry + TransformedGeometry, + + # connectivities + Connectivity, + paramdim, + indices, + connect, + materialize, + pltype, + + # topologies + Topology, + vertex, + vertices, + nvertices, + element, + elements, + nelements, + facet, + facets, + nfacets, + faces, + nfaces, + GridTopology, + HalfEdgeTopology, + HalfEdge, + SimpleTopology, + elem2cart, + cart2elem, + corner2cart, + cart2corner, + elem2corner, + corner2elem, + elementtype, + facettype, + half4elem, + half4vert, + half4edge, + half4pair, + edge4pair, + connec4elem, + + # topological relations + TopologicalRelation, + Boundary, + Coboundary, + Adjacency, + + # domain traits + Domain, + SubDomain, + TransformedDomain, + embeddim, + paramdim, + crs, + manifold, + element, + nelements, + + # sets + GeometrySet, + PointSet, + + # meshes + Mesh, + SimpleMesh, + TransformedMesh, + Grid, + SubGrid, + RegularGrid, + CartesianGrid, + RectilinearGrid, + StructuredGrid, + TransformedGrid, + vertex, + vertices, + nvertices, + eachvertex, + element, + elements, + nelements, + facet, + facets, + nfacets, + topology, + topoconvert, + spacing, + offset, + + # trajectories + CylindricalTrajectory, + + # indices + indices, + + # partitions + Partition, + indices, + metadata, + + # partitioning + PartitionMethod, + IndexPredicatePartitionMethod, + PointPredicatePartitionMethod, + UniformPartition, + FractionPartition, + BlockPartition, + BisectPointPartition, + BisectFractionPartition, + BallPartition, + PlanePartition, + DirectionPartition, + IndexPredicatePartition, + PointPredicatePartition, + ProductPartition, + HierarchicalPartition, + partition, + split, + + # sorting + SortingMethod, + DirectionSort, + sortinds, + + # traversing + Path, + LinearPath, + RandomPath, + ShiftedPath, + SourcePath, + MultiGridPath, + traverse, + + # neighborhoods + Neighborhood, + MetricBall, + hasequalradii, + radii, + rotation, + metric, + radius, + + # neighborhood search + NeighborSearchMethod, + BoundedNeighborSearchMethod, + BallSearch, + KNearestSearch, + KBallSearch, + search!, + searchdists!, + search, + searchdists, + maxneighbors, + + # predicates + isparametrized, + iscurve, + issurface, + issolid, + isperiodic, + issimplex, + isclosed, + isconvex, + issimple, + hasholes, + intersects, + iscollinear, + iscoplanar, + ≺, + ≻, + ⪯, + ⪰, + + # centroids + centroid, + + # measures + measure, + area, + volume, + perimeter, + + # boundary + boundary, + + # winding number + winding, + + # sideof + sideof, + SideType, + IN, + OUT, + ON, + LEFT, + RIGHT, + + # orientation + orientation, + OrientationType, + CW, + CCW, + + # clipping + ClippingMethod, + SutherlandHodgmanClipping, + clip, + + # intersections + IntersectionType, + Crossing, + CornerCrossing, + EdgeCrossing, + Touching, + CornerTouching, + EdgeTouching, + Overlapping, + PosOverlapping, + NegOverlapping, + NotIntersecting, + Intersecting, + + # intersecting + Intersection, + intersection, + type, + + # simplification + SimplificationMethod, + SelingerSimplification, + DouglasPeuckerSimplification, + MinMaxSimplification, + simplify, + + # bounding boxes + boundingbox, + + # hulls + HullMethod, + GrahamScan, + JarvisMarch, + hull, + convexhull, + + # sampling + SamplingMethod, + DiscreteSamplingMethod, + ContinuousSamplingMethod, + UniformSampling, + WeightedSampling, + BallSampling, + BlockSampling, + RegularSampling, + HomogeneousSampling, + MinDistanceSampling, + FibonacciSampling, + sampleinds, + sample, + + # pointification + pointify, + + # tesselation + TesselationMethod, + DelaunayTesselation, + VoronoiTesselation, + tesselate, + + # discretization + DiscretizationMethod, + BoundaryTriangulationMethod, + FanTriangulation, + DehnTriangulation, + HeldTriangulation, + DelaunayTriangulation, + ManualSimplexification, + RegularDiscretization, + MaxLengthDiscretization, + discretize, + discretizewithin, + simplexify, + + # refinement + RefinementMethod, + TriRefinement, + QuadRefinement, + RegularRefinement, + CatmullClarkRefinement, + TriSubdivision, + refine, + + # coarsening + CoarseningMethod, + RegularCoarsening, + coarsen, + + # transforms + GeometricTransform, + CoordinateTransform, + Rotate, + Translate, + Scale, + Affine, + Stretch, + StdCoords, + Proj, + Morphological, + LengthUnit, + ValidCoords, + Shadow, + Slice, + Repair, + Bridge, + LambdaMuSmoothing, + LaplaceSmoothing, + TaubinSmoothing, + isaffine, + isinvertible, + inverse, + + # miscellaneous + signarea, + householderbasis, + supportfun, + laplacematrix, + measurematrix, + adjacencymatrix, + atol, + + # visualization + viz, + viz! end # module diff --git a/src/indices.jl b/src/indices.jl index 0b468211d..ae5bae4bf 100644 --- a/src/indices.jl +++ b/src/indices.jl @@ -7,79 +7,95 @@ # --------------- """ - indices(domain, geometry) + indices(domain, geometry) Return the indices of the elements of the `domain` that intersect with the `geometry`. """ indices(domain::Domain, geometry::Geometry) = findall(intersects(geometry), domain) function indices(grid::OrthoRegularGrid, point::Point) - point ∉ grid && return Int[] + point ∉ grid && return Int[] - # grid properties - orig = minimum(grid) - spac = spacing(grid) - dims = size(grid) + # grid properties + orig = minimum(grid) + spac = spacing(grid) + dims = size(grid) - # integer coordinates - coords = ceil.(Int, (point - orig) ./ spac) + # integer coordinates + coords = ceil.(Int, (point - orig) ./ spac) - # fix coordinates that are on the grid border - coords = clamp.(coords, 1, dims) + # fix coordinates that are on the grid border + coords = clamp.(coords, 1, dims) - # convert to linear index - [LinearIndices(dims)[coords...]] + # convert to linear index + [LinearIndices(dims)[coords...]] end function indices(grid::OrthoRegularGrid, chain::Chain) - dims = size(grid) - mask = falses(dims) + dims = size(grid) + mask = falses(dims) - for segment in segments(chain) - p₁, p₂ = vertices(segment) - _bresenham!(mask, grid, true, p₁, p₂) - end + for segment in segments(chain) + p₁, p₂ = vertices(segment) + _bresenham!(mask, grid, true, p₁, p₂) + end - LinearIndices(dims)[mask] + LinearIndices(dims)[mask] end function indices(grid::OrthoRegularGrid, poly::Polygon) - dims = size(grid) - mask = zeros(Int, dims) - cpoly = poly ∩ boundingbox(grid) - isnothing(cpoly) && return Int[] + dims = size(grid) + mask = zeros(Int, dims) + cpoly = poly ∩ boundingbox(grid) + isnothing(cpoly) && return Int[] - for (i, triangle) in enumerate(simplexify(cpoly)) - _fill!(mask, grid, i, triangle) - end + for (i, triangle) in enumerate(simplexify(cpoly)) + _fill!(mask, grid, i, triangle) + end - LinearIndices(dims)[mask .> 0] + LinearIndices(dims)[mask .> 0] end function indices(grid::OrthoRegularGrid, box::Box) - # cartesian range - range = cartesianrange(grid, box) + # cartesian range + range = cartesianrange(grid, box) - # convert to linear indices - LinearIndices(size(grid))[range] |> vec + # convert to linear indices + LinearIndices(size(grid))[range] |> vec end indices(grid::OrthoRegularGrid, multi::Multi) = mapreduce(geom -> indices(grid, geom), vcat, parent(multi)) |> unique function indices(grid::OrthoRectilinearGrid, box::Box) - # cartesian range - range = cartesianrange(grid, box) + # cartesian range + range = cartesianrange(grid, box) - # convert to linear indices - LinearIndices(size(grid))[range] |> vec + # convert to linear indices + LinearIndices(size(grid))[range] |> vec +end + +# Fallback: inverse if possible, otherwise error +function _reversetransform(t::Transform) + isrevertible(t) ? inverse(t) : throw(ArgumentError("Irreversible transform: $(typeof(t))")) +end + +# Treat specific irreversible transforms as Identity (no-op) +_reversetransform(::Morphological) = Identity() +_reversetransform(::Repair) = Identity() +_reversetransform(::Bridge) = Identity() + +function _reversetransform(st::SequentialTransform) + # Recursively reverse steps, replacing exceptions with Identity + newsteps = map(_reversetransform, st.transforms) + reduce(→, reverse(newsteps)) end function indices(grid::TransformedGrid, geometry::Geometry) - # construct reverse transform from revertible steps - revtrans = reduce(→, reverse(filter(isrevertible, transform(grid)))) + # construct reverse transform from revertible steps + revtrans = _reversetransform(transform(grid)) - # find indices in non-transformed space - indices(parent(grid), revtrans(geometry)) + # find indices in non-transformed space + indices(parent(grid), revtrans(geometry)) end # ---------------- @@ -87,7 +103,7 @@ end # ---------------- """ - cartesianrange(grid, box) + cartesianrange(grid, box) Return the Cartesian range of the elements of the `grid` that intersect with the `box`. """ @@ -98,97 +114,97 @@ _manifoldrange(::Type{<:𝔼}, grid::Grid, box::Box) = _euclideanrange(grid, box _manifoldrange(::Type{<:🌐}, grid::Grid, box::Box) = _geodesicrange(grid, box) function _euclideanrange(grid::OrthoRegularGrid, box::Box) - # grid properties - or = minimum(grid) - sp = spacing(grid) - sz = size(grid) + # grid properties + or = minimum(grid) + sp = spacing(grid) + sz = size(grid) - # intersection of boxes - lo, up = extrema(boundingbox(grid) ∩ box) + # intersection of boxes + lo, up = extrema(boundingbox(grid) ∩ box) - # Cartesian indices of new corners - ijkₛ = max.(ceil.(Int, (lo - or) ./ sp), 1) - ijkₑ = min.(floor.(Int, (up - or) ./ sp) .+ 1, sz) + # Cartesian indices of new corners + ijkₛ = max.(ceil.(Int, (lo - or) ./ sp), 1) + ijkₑ = min.(floor.(Int, (up - or) ./ sp) .+ 1, sz) - # Cartesian range from corner to corner - CartesianIndex(Tuple(ijkₛ)):CartesianIndex(Tuple(ijkₑ)) + # Cartesian range from corner to corner + CartesianIndex(Tuple(ijkₛ)):CartesianIndex(Tuple(ijkₑ)) end function _euclideanrange(grid::OrthoRectilinearGrid, box::Box) - # grid properties - nd = paramdim(grid) + # grid properties + nd = paramdim(grid) - # intersection of boxes - lo, up = to.(extrema(boundingbox(grid) ∩ box)) + # intersection of boxes + lo, up = to.(extrema(boundingbox(grid) ∩ box)) - # integer coordinates of lower point - ijkₛ = ntuple(nd) do i - findlast(x -> x ≤ lo[i], xyz(grid)[i]) - end + # integer coordinates of lower point + ijkₛ = ntuple(nd) do i + findlast(x -> x ≤ lo[i], xyz(grid)[i]) + end - # integer coordinates of upper point - ijkₑ = ntuple(nd) do i - findfirst(x -> x ≥ up[i], xyz(grid)[i]) - end + # integer coordinates of upper point + ijkₑ = ntuple(nd) do i + findfirst(x -> x ≥ up[i], xyz(grid)[i]) + end - # integer coordinates of elements - CartesianIndex(ijkₛ):CartesianIndex(ijkₑ .- 1) + # integer coordinates of elements + CartesianIndex(ijkₛ):CartesianIndex(ijkₑ .- 1) end function _geodesicrange(grid::Grid, box::Box) - nlon, nlat = vsize(grid) - - boxmin = convert(LatLon, coords(minimum(box))) - boxmax = convert(LatLon, coords(maximum(box))) - - a = convert(LatLon, coords(vertex(grid, (1, 1)))) - b = convert(LatLon, coords(vertex(grid, (nlon, 1)))) - c = convert(LatLon, coords(vertex(grid, (1, nlat)))) - - swaplon = a.lon > b.lon - swaplat = a.lat > c.lat - - loninds = swaplon ? (nlon:-1:1) : (1:1:nlon) - latinds = swaplat ? (nlat:-1:1) : (1:1:nlat) - - gridlonₛ, gridlonₑ = swaplon ? (b.lon, a.lon) : (a.lon, b.lon) - gridlatₛ, gridlatₑ = swaplat ? (c.lat, a.lat) : (a.lat, c.lat) - - lonmin = max(boxmin.lon, gridlonₛ) - latmin = max(boxmin.lat, gridlatₛ) - lonmax = min(boxmax.lon, gridlonₑ) - latmax = min(boxmax.lat, gridlatₑ) - - iₛ = findlast(loninds) do i - p = vertex(grid, (i, 1)) - c = convert(LatLon, coords(p)) - c.lon ≤ lonmin - end - iₑ = findfirst(loninds) do i - p = vertex(grid, (i, 1)) - c = convert(LatLon, coords(p)) - c.lon ≥ lonmax - end - - jₛ = findlast(latinds) do i - p = vertex(grid, (1, i)) - c = convert(LatLon, coords(p)) - c.lat ≤ latmin - end - jₑ = findfirst(latinds) do i - p = vertex(grid, (1, i)) - c = convert(LatLon, coords(p)) - c.lat ≥ latmax - end - - if iₛ == iₑ || jₛ == jₑ - throw(ArgumentError("the passed limits are not valid for the grid")) - end - - iₛ, iₑ = swaplon ? (iₑ, iₛ) : (iₛ, iₑ) - jₛ, jₑ = swaplat ? (jₑ, jₛ) : (jₛ, jₑ) - - CartesianIndex(loninds[iₛ], latinds[jₛ]):CartesianIndex(loninds[iₑ] - 1, latinds[jₑ] - 1) + nlon, nlat = vsize(grid) + + boxmin = convert(LatLon, coords(minimum(box))) + boxmax = convert(LatLon, coords(maximum(box))) + + a = convert(LatLon, coords(vertex(grid, (1, 1)))) + b = convert(LatLon, coords(vertex(grid, (nlon, 1)))) + c = convert(LatLon, coords(vertex(grid, (1, nlat)))) + + swaplon = a.lon > b.lon + swaplat = a.lat > c.lat + + loninds = swaplon ? (nlon:-1:1) : (1:1:nlon) + latinds = swaplat ? (nlat:-1:1) : (1:1:nlat) + + gridlonₛ, gridlonₑ = swaplon ? (b.lon, a.lon) : (a.lon, b.lon) + gridlatₛ, gridlatₑ = swaplat ? (c.lat, a.lat) : (a.lat, c.lat) + + lonmin = max(boxmin.lon, gridlonₛ) + latmin = max(boxmin.lat, gridlatₛ) + lonmax = min(boxmax.lon, gridlonₑ) + latmax = min(boxmax.lat, gridlatₑ) + + iₛ = findlast(loninds) do i + p = vertex(grid, (i, 1)) + c = convert(LatLon, coords(p)) + c.lon ≤ lonmin + end + iₑ = findfirst(loninds) do i + p = vertex(grid, (i, 1)) + c = convert(LatLon, coords(p)) + c.lon ≥ lonmax + end + + jₛ = findlast(latinds) do i + p = vertex(grid, (1, i)) + c = convert(LatLon, coords(p)) + c.lat ≤ latmin + end + jₑ = findfirst(latinds) do i + p = vertex(grid, (1, i)) + c = convert(LatLon, coords(p)) + c.lat ≥ latmax + end + + if iₛ == iₑ || jₛ == jₑ + throw(ArgumentError("the passed limits are not valid for the grid")) + end + + iₛ, iₑ = swaplon ? (iₑ, iₛ) : (iₛ, iₑ) + jₛ, jₑ = swaplat ? (jₑ, jₛ) : (jₛ, jₑ) + + CartesianIndex(loninds[iₛ], latinds[jₛ]):CartesianIndex(loninds[iₑ]-1, latinds[jₑ]-1) end # ----------------- @@ -196,98 +212,98 @@ end # ----------------- function _fill!(mask, grid, val, triangle) - v = vertices(triangle) - - # fill edges of triangle - _bresenham!(mask, grid, val, v[1], v[2]) - _bresenham!(mask, grid, val, v[2], v[3]) - _bresenham!(mask, grid, val, v[3], v[1]) - - # fill interior of triangle - j₁ = findfirst(==(val), mask).I[2] - j₂ = findlast(==(val), mask).I[2] - for j in j₁:j₂ - i₁ = findfirst(==(val), @view(mask[:, j])) - i₂ = findlast(==(val), @view(mask[:, j])) - mask[i₁:i₂, j] .= val - end + v = vertices(triangle) + + # fill edges of triangle + _bresenham!(mask, grid, val, v[1], v[2]) + _bresenham!(mask, grid, val, v[2], v[3]) + _bresenham!(mask, grid, val, v[3], v[1]) + + # fill interior of triangle + j₁ = findfirst(==(val), mask).I[2] + j₂ = findlast(==(val), mask).I[2] + for j in j₁:j₂ + i₁ = findfirst(==(val), @view(mask[:, j])) + i₂ = findlast(==(val), @view(mask[:, j])) + mask[i₁:i₂, j] .= val + end end # Bresenham's line algorithm: https://en.wikipedia.org/wiki/Bresenham's_line_algorithm function _bresenham!(mask, grid, val, p₁, p₂) - o = minimum(grid) - s = spacing(grid) - - # integer coordinates - x₁, y₁ = ceil.(Int, (p₁ - o) ./ s) - x₂, y₂ = ceil.(Int, (p₂ - o) ./ s) - - # fix coordinates of points that are on the grid border - xmax, ymax = size(grid) - x₁ = clamp(x₁, 1, xmax) - y₁ = clamp(y₁, 1, ymax) - x₂ = clamp(x₂, 1, xmax) - y₂ = clamp(y₂, 1, ymax) - - if abs(y₂ - y₁) < abs(x₂ - x₁) - if x₁ > x₂ - _bresenhamlow!(mask, val, x₂, y₂, x₁, y₁) - else - _bresenhamlow!(mask, val, x₁, y₁, x₂, y₂) - end - else - if y₁ > y₂ - _bresenhamhigh!(mask, val, x₂, y₂, x₁, y₁) - else - _bresenhamhigh!(mask, val, x₁, y₁, x₂, y₂) - end - end + o = minimum(grid) + s = spacing(grid) + + # integer coordinates + x₁, y₁ = ceil.(Int, (p₁ - o) ./ s) + x₂, y₂ = ceil.(Int, (p₂ - o) ./ s) + + # fix coordinates of points that are on the grid border + xmax, ymax = size(grid) + x₁ = clamp(x₁, 1, xmax) + y₁ = clamp(y₁, 1, ymax) + x₂ = clamp(x₂, 1, xmax) + y₂ = clamp(y₂, 1, ymax) + + if abs(y₂ - y₁) < abs(x₂ - x₁) + if x₁ > x₂ + _bresenhamlow!(mask, val, x₂, y₂, x₁, y₁) + else + _bresenhamlow!(mask, val, x₁, y₁, x₂, y₂) + end + else + if y₁ > y₂ + _bresenhamhigh!(mask, val, x₂, y₂, x₁, y₁) + else + _bresenhamhigh!(mask, val, x₁, y₁, x₂, y₂) + end + end end function _bresenhamlow!(mask, val, x₁, y₁, x₂, y₂) - dx = x₂ - x₁ - dy = y₂ - y₁ - yi = 1 - if dy < 0 - yi = -1 - dy = -dy - end - - D = 2dy - dx - y = y₁ - - for x in x₁:x₂ - mask[x, y] = val - - if D > 0 - y = y + yi - D = D + 2dy - 2dx - else - D = D + 2dy - end - end + dx = x₂ - x₁ + dy = y₂ - y₁ + yi = 1 + if dy < 0 + yi = -1 + dy = -dy + end + + D = 2dy - dx + y = y₁ + + for x in x₁:x₂ + mask[x, y] = val + + if D > 0 + y = y + yi + D = D + 2dy - 2dx + else + D = D + 2dy + end + end end function _bresenhamhigh!(mask, val, x₁, y₁, x₂, y₂) - dx = x₂ - x₁ - dy = y₂ - y₁ - xi = 1 - if dx < 0 - xi = -1 - dx = -dx - end - - D = 2dx - dy - x = x₁ - - for y in y₁:y₂ - mask[x, y] = val - - if D > 0 - x = x + xi - D = D + 2dx - 2dy - else - D = D + 2dx - end - end + dx = x₂ - x₁ + dy = y₂ - y₁ + xi = 1 + if dx < 0 + xi = -1 + dx = -dx + end + + D = 2dx - dy + x = x₁ + + for y in y₁:y₂ + mask[x, y] = val + + if D > 0 + x = x + xi + D = D + 2dx - 2dy + else + D = D + 2dx + end + end end From 308b316b82862ad29ce6ebf919c7054b0faadada Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva Date: Fri, 21 Nov 2025 08:21:01 -0300 Subject: [PATCH 25/27] update --- src/Meshes.jl | 3 --- 1 file changed, 3 deletions(-) diff --git a/src/Meshes.jl b/src/Meshes.jl index 7802ce259..fa9b22d27 100644 --- a/src/Meshes.jl +++ b/src/Meshes.jl @@ -256,11 +256,8 @@ export angles, innerangles, normal, -<<<<<<< HEAD base, apex, -======= ->>>>>>> c7c91f967fcb868cac7b2a34b9aeecfad09490a8 ≗, # multi-geometries From 281ff849cac022d74283164b915a7a3aceab4704 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva <135757954+marcosdanieldasilva@users.noreply.github.com> Date: Fri, 21 Nov 2025 08:28:02 -0300 Subject: [PATCH 26/27] tab in Meshes.jl --- src/Meshes.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/Meshes.jl b/src/Meshes.jl index fa9b22d27..2a17af51d 100644 --- a/src/Meshes.jl +++ b/src/Meshes.jl @@ -32,7 +32,7 @@ using ScopedValues: ScopedValue using Base.Cartesian: @nloops, @nref, @ntuple using Base: @propagate_inbounds -using Random: Random +import Random import Base: sort import Base: ==, ! import Base: +, -, * From cc2d07c0b31faffbf58a4f47bbfd44c80c1601d9 Mon Sep 17 00:00:00 2001 From: Marcos Daniel da Silva Date: Fri, 21 Nov 2025 08:36:15 -0300 Subject: [PATCH 27/27] formatter --- src/Meshes.jl | 928 +++++++++++++++++++++++++------------------------- 1 file changed, 463 insertions(+), 465 deletions(-) diff --git a/src/Meshes.jl b/src/Meshes.jl index 2a17af51d..69f1efc0c 100644 --- a/src/Meshes.jl +++ b/src/Meshes.jl @@ -123,8 +123,6 @@ include("discretization.jl") include("refinement.jl") include("coarsening.jl") - - # miscellaneous include("rand.jl") include("distances.jl") @@ -136,471 +134,471 @@ include("projecting.jl") include("viz.jl") function __init__() - # register error hint for visualization functions - # since this is a recurring issue for new users - Base.Experimental.register_error_hint(MethodError) do io, exc, argtypes, kwargs - if exc.f == viz || exc.f == viz! - if isnothing(Base.get_extension(Meshes, :MeshesMakieExt)) - print( - io, - """ - - Did you import a Makie.jl backend (e.g., GLMakie.jl, CairoMakie.jl) for visualization? - - """, - ) - printstyled( - io, - """ - julia> using Meshes - julia> import GLMakie # or CairoMakie, WGLMakie, etc. - """, - color=:cyan, - bold=true, - ) - end - end - end + # register error hint for visualization functions + # since this is a recurring issue for new users + Base.Experimental.register_error_hint(MethodError) do io, exc, argtypes, kwargs + if exc.f == viz || exc.f == viz! + if isnothing(Base.get_extension(Meshes, :MeshesMakieExt)) + print( + io, + """ + + Did you import a Makie.jl backend (e.g., GLMakie.jl, CairoMakie.jl) for visualization? + + """ + ) + printstyled( + io, + """ + julia> using Meshes + julia> import GLMakie # or CairoMakie, WGLMakie, etc. + """, + color=:cyan, + bold=true + ) + end + end + end end export - # vectors - Vec, - ∠, - ⋅, - ×, - - # manifolds - 𝔼, - 🌐, - - # geometries - Geometry, - embeddim, - paramdim, - crs, - manifold, - - # primitives - Primitive, - Point, - Ray, - Line, - BezierCurve, - ParametrizedCurve, - Plane, - Box, - Ball, - Sphere, - Ellipsoid, - Disk, - Circle, - Cylinder, - CylinderSurface, - ParaboloidSurface, - Cone, - ConeSurface, - Frustum, - FrustumSurface, - Torus, - controls, - ncontrols, - degree, - Horner, - DeCasteljau, - coords, - to, - center, - radius, - radii, - plane, - bottom, - top, - axis, - base, - apex, - isright, - sides, - diagonal, - focallength, - direction, - - # polytopes - Polytope, - Chain, - Segment, - Rope, - Ring, - Polygon, - Ngon, - Triangle, - Quadrangle, - Pentagon, - Hexagon, - Heptagon, - Octagon, - Nonagon, - Decagon, - PolyArea, - Polyhedron, - Tetrahedron, - Hexahedron, - Pyramid, - Wedge, - vertex, - vertices, - nvertices, - eachvertex, - rings, - segments, - angles, - innerangles, - normal, - base, - apex, - ≗, - - # multi-geometries - Multi, - MultiPoint, - MultiSegment, - MultiRope, - MultiRing, - MultiPolygon, - MultiPolyhedron, - - # transformed geometry - TransformedGeometry, - - # connectivities - Connectivity, - paramdim, - indices, - connect, - materialize, - pltype, - - # topologies - Topology, - vertex, - vertices, - nvertices, - element, - elements, - nelements, - facet, - facets, - nfacets, - faces, - nfaces, - GridTopology, - HalfEdgeTopology, - HalfEdge, - SimpleTopology, - elem2cart, - cart2elem, - corner2cart, - cart2corner, - elem2corner, - corner2elem, - elementtype, - facettype, - half4elem, - half4vert, - half4edge, - half4pair, - edge4pair, - connec4elem, - - # topological relations - TopologicalRelation, - Boundary, - Coboundary, - Adjacency, - - # domain traits - Domain, - SubDomain, - TransformedDomain, - embeddim, - paramdim, - crs, - manifold, - element, - nelements, - - # sets - GeometrySet, - PointSet, - - # meshes - Mesh, - SimpleMesh, - TransformedMesh, - Grid, - SubGrid, - RegularGrid, - CartesianGrid, - RectilinearGrid, - StructuredGrid, - TransformedGrid, - vertex, - vertices, - nvertices, - eachvertex, - element, - elements, - nelements, - facet, - facets, - nfacets, - topology, - topoconvert, - spacing, - offset, - - # trajectories - CylindricalTrajectory, - - # indices - indices, - - # partitions - Partition, - indices, - metadata, - - # partitioning - PartitionMethod, - IndexPredicatePartitionMethod, - PointPredicatePartitionMethod, - UniformPartition, - FractionPartition, - BlockPartition, - BisectPointPartition, - BisectFractionPartition, - BallPartition, - PlanePartition, - DirectionPartition, - IndexPredicatePartition, - PointPredicatePartition, - ProductPartition, - HierarchicalPartition, - partition, - split, - - # sorting - SortingMethod, - DirectionSort, - sortinds, - - # traversing - Path, - LinearPath, - RandomPath, - ShiftedPath, - SourcePath, - MultiGridPath, - traverse, - - # neighborhoods - Neighborhood, - MetricBall, - hasequalradii, - radii, - rotation, - metric, - radius, - - # neighborhood search - NeighborSearchMethod, - BoundedNeighborSearchMethod, - BallSearch, - KNearestSearch, - KBallSearch, - search!, - searchdists!, - search, - searchdists, - maxneighbors, - - # predicates - isparametrized, - iscurve, - issurface, - issolid, - isperiodic, - issimplex, - isclosed, - isconvex, - issimple, - hasholes, - intersects, - iscollinear, - iscoplanar, - ≺, - ≻, - ⪯, - ⪰, - - # centroids - centroid, - - # measures - measure, - area, - volume, - perimeter, - - # boundary - boundary, - - # winding number - winding, - - # sideof - sideof, - SideType, - IN, - OUT, - ON, - LEFT, - RIGHT, - - # orientation - orientation, - OrientationType, - CW, - CCW, - - # clipping - ClippingMethod, - SutherlandHodgmanClipping, - clip, - - # intersections - IntersectionType, - Crossing, - CornerCrossing, - EdgeCrossing, - Touching, - CornerTouching, - EdgeTouching, - Overlapping, - PosOverlapping, - NegOverlapping, - NotIntersecting, - Intersecting, - - # intersecting - Intersection, - intersection, - type, - - # simplification - SimplificationMethod, - SelingerSimplification, - DouglasPeuckerSimplification, - MinMaxSimplification, - simplify, - - # bounding boxes - boundingbox, - - # hulls - HullMethod, - GrahamScan, - JarvisMarch, - hull, - convexhull, - - # sampling - SamplingMethod, - DiscreteSamplingMethod, - ContinuousSamplingMethod, - UniformSampling, - WeightedSampling, - BallSampling, - BlockSampling, - RegularSampling, - HomogeneousSampling, - MinDistanceSampling, - FibonacciSampling, - sampleinds, - sample, - - # pointification - pointify, - - # tesselation - TesselationMethod, - DelaunayTesselation, - VoronoiTesselation, - tesselate, - - # discretization - DiscretizationMethod, - BoundaryTriangulationMethod, - FanTriangulation, - DehnTriangulation, - HeldTriangulation, - DelaunayTriangulation, - ManualSimplexification, - RegularDiscretization, - MaxLengthDiscretization, - discretize, - discretizewithin, - simplexify, - - # refinement - RefinementMethod, - TriRefinement, - QuadRefinement, - RegularRefinement, - CatmullClarkRefinement, - TriSubdivision, - refine, - - # coarsening - CoarseningMethod, - RegularCoarsening, - coarsen, - - # transforms - GeometricTransform, - CoordinateTransform, - Rotate, - Translate, - Scale, - Affine, - Stretch, - StdCoords, - Proj, - Morphological, - LengthUnit, - ValidCoords, - Shadow, - Slice, - Repair, - Bridge, - LambdaMuSmoothing, - LaplaceSmoothing, - TaubinSmoothing, - isaffine, - isinvertible, - inverse, - - # miscellaneous - signarea, - householderbasis, - supportfun, - laplacematrix, - measurematrix, - adjacencymatrix, - atol, - - # visualization - viz, - viz! + # vectors + Vec, + ∠, + ⋅, + ×, + + # manifolds + 𝔼, + 🌐, + + # geometries + Geometry, + embeddim, + paramdim, + crs, + manifold, + + # primitives + Primitive, + Point, + Ray, + Line, + BezierCurve, + ParametrizedCurve, + Plane, + Box, + Ball, + Sphere, + Ellipsoid, + Disk, + Circle, + Cylinder, + CylinderSurface, + ParaboloidSurface, + Cone, + ConeSurface, + Frustum, + FrustumSurface, + Torus, + controls, + ncontrols, + degree, + Horner, + DeCasteljau, + coords, + to, + center, + radius, + radii, + plane, + bottom, + top, + axis, + base, + apex, + isright, + sides, + diagonal, + focallength, + direction, + + # polytopes + Polytope, + Chain, + Segment, + Rope, + Ring, + Polygon, + Ngon, + Triangle, + Quadrangle, + Pentagon, + Hexagon, + Heptagon, + Octagon, + Nonagon, + Decagon, + PolyArea, + Polyhedron, + Tetrahedron, + Hexahedron, + Pyramid, + Wedge, + vertex, + vertices, + nvertices, + eachvertex, + rings, + segments, + angles, + innerangles, + normal, + base, + apex, + ≗, + + # multi-geometries + Multi, + MultiPoint, + MultiSegment, + MultiRope, + MultiRing, + MultiPolygon, + MultiPolyhedron, + + # transformed geometry + TransformedGeometry, + + # connectivities + Connectivity, + paramdim, + indices, + connect, + materialize, + pltype, + + # topologies + Topology, + vertex, + vertices, + nvertices, + element, + elements, + nelements, + facet, + facets, + nfacets, + faces, + nfaces, + GridTopology, + HalfEdgeTopology, + HalfEdge, + SimpleTopology, + elem2cart, + cart2elem, + corner2cart, + cart2corner, + elem2corner, + corner2elem, + elementtype, + facettype, + half4elem, + half4vert, + half4edge, + half4pair, + edge4pair, + connec4elem, + + # topological relations + TopologicalRelation, + Boundary, + Coboundary, + Adjacency, + + # domain traits + Domain, + SubDomain, + TransformedDomain, + embeddim, + paramdim, + crs, + manifold, + element, + nelements, + + # sets + GeometrySet, + PointSet, + + # meshes + Mesh, + SimpleMesh, + TransformedMesh, + Grid, + SubGrid, + RegularGrid, + CartesianGrid, + RectilinearGrid, + StructuredGrid, + TransformedGrid, + vertex, + vertices, + nvertices, + eachvertex, + element, + elements, + nelements, + facet, + facets, + nfacets, + topology, + topoconvert, + spacing, + offset, + + # trajectories + CylindricalTrajectory, + + # indices + indices, + + # partitions + Partition, + indices, + metadata, + + # partitioning + PartitionMethod, + IndexPredicatePartitionMethod, + PointPredicatePartitionMethod, + UniformPartition, + FractionPartition, + BlockPartition, + BisectPointPartition, + BisectFractionPartition, + BallPartition, + PlanePartition, + DirectionPartition, + IndexPredicatePartition, + PointPredicatePartition, + ProductPartition, + HierarchicalPartition, + partition, + split, + + # sorting + SortingMethod, + DirectionSort, + sortinds, + + # traversing + Path, + LinearPath, + RandomPath, + ShiftedPath, + SourcePath, + MultiGridPath, + traverse, + + # neighborhoods + Neighborhood, + MetricBall, + hasequalradii, + radii, + rotation, + metric, + radius, + + # neighborhood search + NeighborSearchMethod, + BoundedNeighborSearchMethod, + BallSearch, + KNearestSearch, + KBallSearch, + search!, + searchdists!, + search, + searchdists, + maxneighbors, + + # predicates + isparametrized, + iscurve, + issurface, + issolid, + isperiodic, + issimplex, + isclosed, + isconvex, + issimple, + hasholes, + intersects, + iscollinear, + iscoplanar, + ≺, + ≻, + ⪯, + ⪰, + + # centroids + centroid, + + # measures + measure, + area, + volume, + perimeter, + + # boundary + boundary, + + # winding number + winding, + + # sideof + sideof, + SideType, + IN, + OUT, + ON, + LEFT, + RIGHT, + + # orientation + orientation, + OrientationType, + CW, + CCW, + + # clipping + ClippingMethod, + SutherlandHodgmanClipping, + clip, + + # intersections + IntersectionType, + Crossing, + CornerCrossing, + EdgeCrossing, + Touching, + CornerTouching, + EdgeTouching, + Overlapping, + PosOverlapping, + NegOverlapping, + NotIntersecting, + Intersecting, + + # intersecting + Intersection, + intersection, + type, + + # simplification + SimplificationMethod, + SelingerSimplification, + DouglasPeuckerSimplification, + MinMaxSimplification, + simplify, + + # bounding boxes + boundingbox, + + # hulls + HullMethod, + GrahamScan, + JarvisMarch, + hull, + convexhull, + + # sampling + SamplingMethod, + DiscreteSamplingMethod, + ContinuousSamplingMethod, + UniformSampling, + WeightedSampling, + BallSampling, + BlockSampling, + RegularSampling, + HomogeneousSampling, + MinDistanceSampling, + FibonacciSampling, + sampleinds, + sample, + + # pointification + pointify, + + # tesselation + TesselationMethod, + DelaunayTesselation, + VoronoiTesselation, + tesselate, + + # discretization + DiscretizationMethod, + BoundaryTriangulationMethod, + FanTriangulation, + DehnTriangulation, + HeldTriangulation, + DelaunayTriangulation, + ManualSimplexification, + RegularDiscretization, + MaxLengthDiscretization, + discretize, + discretizewithin, + simplexify, + + # refinement + RefinementMethod, + TriRefinement, + QuadRefinement, + RegularRefinement, + CatmullClarkRefinement, + TriSubdivision, + refine, + + # coarsening + CoarseningMethod, + RegularCoarsening, + coarsen, + + # transforms + GeometricTransform, + CoordinateTransform, + Rotate, + Translate, + Scale, + Affine, + Stretch, + StdCoords, + Proj, + Morphological, + LengthUnit, + ValidCoords, + Shadow, + Slice, + Repair, + Bridge, + LambdaMuSmoothing, + LaplaceSmoothing, + TaubinSmoothing, + isaffine, + isinvertible, + inverse, + + # miscellaneous + signarea, + householderbasis, + supportfun, + laplacematrix, + measurematrix, + adjacencymatrix, + atol, + + # visualization + viz, + viz! end # module