From f56acab9d4c6c214263f14acdbc8565a2d052686 Mon Sep 17 00:00:00 2001 From: t-bltg Date: Mon, 26 Jan 2026 12:36:22 +0100 Subject: [PATCH 1/3] update actions --- .JuliaFormatter | 14 -------- .github/dependabot.yml | 15 ++++++++ .github/workflows/format-check.yml | 36 -------------------- .github/workflows/format.yml | 21 ++++++++++++ .github/workflows/{TagBot.yml => tagbot.yml} | 9 +++-- .pre-commit-config.yaml | 15 ++++++++ 6 files changed, 58 insertions(+), 52 deletions(-) delete mode 100644 .JuliaFormatter create mode 100644 .github/dependabot.yml delete mode 100644 .github/workflows/format-check.yml create mode 100644 .github/workflows/format.yml rename .github/workflows/{TagBot.yml => tagbot.yml} (74%) create mode 100644 .pre-commit-config.yaml diff --git a/.JuliaFormatter b/.JuliaFormatter deleted file mode 100644 index 3992a04..0000000 --- a/.JuliaFormatter +++ /dev/null @@ -1,14 +0,0 @@ -annotate_untyped_fields_with_any = false -short_to_long_function_def = false -long_to_short_function_def = true -whitespace_ops_in_indices = true -remove_extra_newlines = true -whitespace_in_kwargs = true -always_use_return = false -conditional_to_if = false -align_struct_field = true -align_conditional = true -align_assignment = true -align_pair_arrow = true -import_to_using = false -always_for_in = true diff --git a/.github/dependabot.yml b/.github/dependabot.yml new file mode 100644 index 0000000..14eb10f --- /dev/null +++ b/.github/dependabot.yml @@ -0,0 +1,15 @@ +# To get started with Dependabot version updates, you'll need to specify which +# package ecosystems to update and where the package manifests are located. +# Please see the documentation for all configuration options: +# https://docs.github.com/github/administering-a-repository/configuration-options-for-dependency-updates + +version: 2 +updates: + - package-ecosystem: "github-actions" # See documentation for possible values + directory: "/" # Location of package manifests + schedule: + interval: "weekly" + - package-ecosystem: "julia" + directory: "/" + schedule: + interval: "weekly" diff --git a/.github/workflows/format-check.yml b/.github/workflows/format-check.yml deleted file mode 100644 index 41fbef0..0000000 --- a/.github/workflows/format-check.yml +++ /dev/null @@ -1,36 +0,0 @@ -name: format - -on: - pull_request: - push: - branches: [main] - -concurrency: - group: ${{ github.workflow }}-${{ github.head_ref || github.run_id }} - cancel-in-progress: true - -jobs: - check: - runs-on: ubuntu-latest - steps: - - uses: julia-actions/setup-julia@latest - with: - version: '1' - - uses: actions/checkout@v5 - - name: Install JuliaFormatter and format - run: | - julia -e 'using Pkg; Pkg.add(PackageSpec(name="JuliaFormatter"))' - julia -e 'using JuliaFormatter; format(["src", "test"], verbose=true)' - - name: Format check - run: | - julia -e ' - out = Cmd(`git diff --name-only`) |> read |> String - if out == "" - exit(0) - else - @error "Some files have not been formatted !!!" - write(stdout, out) - run(`git diff`) - exit(1) - end - ' diff --git a/.github/workflows/format.yml b/.github/workflows/format.yml new file mode 100644 index 0000000..058754a --- /dev/null +++ b/.github/workflows/format.yml @@ -0,0 +1,21 @@ +name: format + +on: + push: + branches: ['release-', 'main'] + tags: + - '*' + pull_request: + +concurrency: + group: ${{ github.workflow }}-${{ github.head_ref || github.run_id }} + cancel-in-progress: true + +jobs: + runic: + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@v6 + - uses: julia-actions/setup-julia@latest + - uses: julia-actions/cache@v2 + - uses: fredrikekre/runic-action@v1 diff --git a/.github/workflows/TagBot.yml b/.github/workflows/tagbot.yml similarity index 74% rename from .github/workflows/TagBot.yml rename to .github/workflows/tagbot.yml index 623860f..4a5ae6c 100644 --- a/.github/workflows/TagBot.yml +++ b/.github/workflows/tagbot.yml @@ -1,15 +1,20 @@ name: TagBot + on: issue_comment: types: - created workflow_dispatch: + inputs: + lookback: + default: 10 + jobs: - TagBot: + tagbot: if: github.event_name == 'workflow_dispatch' || github.actor == 'JuliaTagBot' runs-on: ubuntu-latest steps: - uses: JuliaRegistries/TagBot@v1 with: token: ${{ secrets.GITHUB_TOKEN }} - ssh: ${{ secrets.DOCUMENTER_KEY }} \ No newline at end of file + ssh: ${{ secrets.DOCUMENTER_KEY }} diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml new file mode 100644 index 0000000..429943e --- /dev/null +++ b/.pre-commit-config.yaml @@ -0,0 +1,15 @@ +# See https://pre-commit.com for more information +# See https://pre-commit.com/hooks.html for more hooks +repos: +- repo: https://github.com/pre-commit/pre-commit-hooks + rev: v5.0.0 + hooks: + - id: trailing-whitespace + - id: end-of-file-fixer + - id: check-toml + - id: check-yaml + - id: check-added-large-files +- repo: https://github.com/fredrikekre/runic-pre-commit + rev: v2.0.1 + hooks: + - id: runic From 46c2d4bff4c1761bcbc440e704323a219cbeb300 Mon Sep 17 00:00:00 2001 From: t-bltg Date: Mon, 26 Jan 2026 12:37:40 +0100 Subject: [PATCH 2/3] format --- ref/check.jl | 184 ++++--- src/MarchingCubes.jl | 1246 +++++++++++++++++++++--------------------- src/example.jl | 48 +- test/runtests.jl | 22 +- 4 files changed, 753 insertions(+), 747 deletions(-) diff --git a/ref/check.jl b/ref/check.jl index b5de72b..42fe191 100644 --- a/ref/check.jl +++ b/ref/check.jl @@ -1,103 +1,109 @@ using MarchingCubes main() = begin - run(`gcc main.c`) + run(`gcc main.c`) - println("== 2d ==") - for (l, sym) ∈ enumerate(( - :cases, - :tiling1, - :tiling2, - :tiling3_1, - :tiling3_2, - :tiling4_1, - :tiling4_2, - :tiling5, - :test6, - :tiling6_1_1, - :tiling6_1_2, - :tiling6_2, - :test7, - :tiling7_1, - :tiling7_4_1, - :tiling7_4_2, - :tiling8, - :tiling9, - :test10, - :tiling10_1_1, - :tiling10_1_1_, - :tiling10_1_2, - :tiling10_2, - :tiling10_2_, - :tiling11, - :test12, - :tiling12_1_1, - :tiling12_1_1_, - :tiling12_1_2, - :tiling12_2, - :tiling12_2_, - :test13, - :tiling13_1, - :tiling13_1_, - :tiling14, - :casesClassic, - )) - @show sym, l - lut = getfield(MarchingCubes, sym) - o = sym ∈ (:test6, :test7, :test10, :test12, :test13) ? 0 : 1 - for i ∈ 1:length(lut) - for j ∈ 1:length(first(lut)) - io = PipeBuffer() - pipeline(`./a.out $l $(i - 1) $(j - 1)`; stdout=io, stderr=devnull) |> run - val = parse(Int, read(io, String) |> rstrip) + o - @assert val == lut[i][j] ((i, j), val, lut[i][j]) - end + println("== 2d ==") + for (l, sym) in enumerate( + ( + :cases, + :tiling1, + :tiling2, + :tiling3_1, + :tiling3_2, + :tiling4_1, + :tiling4_2, + :tiling5, + :test6, + :tiling6_1_1, + :tiling6_1_2, + :tiling6_2, + :test7, + :tiling7_1, + :tiling7_4_1, + :tiling7_4_2, + :tiling8, + :tiling9, + :test10, + :tiling10_1_1, + :tiling10_1_1_, + :tiling10_1_2, + :tiling10_2, + :tiling10_2_, + :tiling11, + :test12, + :tiling12_1_1, + :tiling12_1_1_, + :tiling12_1_2, + :tiling12_2, + :tiling12_2_, + :test13, + :tiling13_1, + :tiling13_1_, + :tiling14, + :casesClassic, + ) + ) + @show sym, l + lut = getfield(MarchingCubes, sym) + o = sym ∈ (:test6, :test7, :test10, :test12, :test13) ? 0 : 1 + for i in 1:length(lut) + for j in 1:length(first(lut)) + io = PipeBuffer() + pipeline(`./a.out $l $(i - 1) $(j - 1)`; stdout = io, stderr = devnull) |> run + val = parse(Int, read(io, String) |> rstrip) + o + @assert val == lut[i][j] ((i, j), val, lut[i][j]) + end + end end - end - println("== 1d ==") - for (l, sym) ∈ enumerate(( - :test3, - :test4, - :subcfg13 - )) - @show sym, l - lut = getfield(MarchingCubes, sym) - for i ∈ 1:length(lut) - io = PipeBuffer() - pipeline(`./a.out $(100 + l) $(i - 1)`; stdout=io, stderr=devnull) |> run - val = parse(Int, read(io, String) |> rstrip) - @assert val == lut[i] (i, val, lut[i]) + println("== 1d ==") + for (l, sym) in enumerate( + ( + :test3, + :test4, + :subcfg13, + ) + ) + @show sym, l + lut = getfield(MarchingCubes, sym) + for i in 1:length(lut) + io = PipeBuffer() + pipeline(`./a.out $(100 + l) $(i - 1)`; stdout = io, stderr = devnull) |> run + val = parse(Int, read(io, String) |> rstrip) + @assert val == lut[i] (i, val, lut[i]) + end end - end - println("== 3d ==") - for (l, sym) ∈ enumerate(( - :tiling7_2, - :tiling7_3, - :tiling13_2, - :tiling13_2_, - :tiling13_3, - :tiling13_3_, - :tiling13_4, - :tiling13_5_1, - :tiling13_5_2, - )) - @show sym, l - lut = getfield(MarchingCubes, sym) - for i ∈ 1:length(lut) - for j ∈ 1:length(first(lut)) - for k ∈ 1:length(first(first(lut))) - io = PipeBuffer() - pipeline(`./a.out $(200 + l) $(i - 1) $(j - 1) $(k - 1)`; stdout=io, stderr=devnull) |> run - val = parse(Int, read(io, String) |> rstrip) + 1 - @assert val == lut[i][j][k] ((i, j, k), val, lut[i][j][k]) + println("== 3d ==") + for (l, sym) in enumerate( + ( + :tiling7_2, + :tiling7_3, + :tiling13_2, + :tiling13_2_, + :tiling13_3, + :tiling13_3_, + :tiling13_4, + :tiling13_5_1, + :tiling13_5_2, + ) + ) + @show sym, l + lut = getfield(MarchingCubes, sym) + for i in 1:length(lut) + for j in 1:length(first(lut)) + for k in 1:length(first(first(lut))) + io = PipeBuffer() + pipeline(`./a.out $(200 + l) $(i - 1) $(j - 1) $(k - 1)`; stdout = io, stderr = devnull) |> run + val = parse(Int, read(io, String) |> rstrip) + 1 + @assert val == lut[i][j][k] ((i, j, k), val, lut[i][j][k]) + end + end end - end end - end - return + return end main() diff --git a/src/MarchingCubes.jl b/src/MarchingCubes.jl index e68a14c..347f047 100644 --- a/src/MarchingCubes.jl +++ b/src/MarchingCubes.jl @@ -19,688 +19,688 @@ Adapted to `Julia` by T Bltg (github.com/t-bltg). module MarchingCubes -using PrecompileTools -using StaticArrays - -import Base: RefValue - -export MC, march, march_legacy - -include("lut.jl") - -const Vertex{F} = SVector{3,F} -const Normal{F} = SVector{3,F} -const Triangle = SVector{3,Int} - -""" - MC(vol::Array{F,3}, I::Type = Int; x = [], y = [], z = []) - -# Description -Structure to hold temporaries and output of the marching cubes algorithm. -Vertices, normals and triangles are accessible using `m.vertices`, `m.normals` and `m.triangles`. -1-based indexing is assumed. -Optionally, pass `x`, `y`, `z` to normalize the vertices coordinates instead of `0:nx-1`, `0:ny-1` and `0:nz-1`. - -# Arguments - - `vol`: rank 3 array of floats (volume) on which the Marching Cubes algorithm is applied. - - `I`: Int32, Int64, ...: vertices / normals / triangles index `Integer` type (defaults to `Int`). - - `x` normalize vertex.x to the coordinates of `x`. - - `y` normalize vertex.y to the coordinates of `y`. - - `z` normalize vertex.z to the coordinates of `z`. -""" -struct MC{F,I} - nx::Int # height of the grid - ny::Int # width of the grid - nz::Int # depth of the grid - vol::RefValue{Array{F,3}} # implicit function values sampled on the grid - cube::MVector{8,F} # values of the implicit function on the active cube - tv::MVector{3,I} # triangle buffer - vtx::MVector{3,F} # vertex buffer - nrm::MVector{3,F} # normal buffer - vert_indices::Array{I,4} # pre-computed vertex indices - triangles::Vector{Triangle} # output triangles - vertices::Vector{Vertex{F}} # output vertex positions - normals::Vector{Normal{F}} # output vertex normals - normal_sign::Int # direction of normal vectors (+1 for outward / -1 for inward) - x::RefValue{Vector{F}} - y::RefValue{Vector{F}} - z::RefValue{Vector{F}} - MC( - vol::Array{F,3}, - I::Type{G} = Int; - normal_sign::Integer = 1, - x::AbstractVector{F} = F[], - y::AbstractVector{F} = F[], - z::AbstractVector{F} = F[], - ) where {F<:AbstractFloat,G<:Integer} = begin - isa(x, AbstractRange) && (x = collect(x)) - isa(y, AbstractRange) && (y = collect(y)) - isa(z, AbstractRange) && (z = collect(z)) - abs(normal_sign) == 1 || - throw(ArgumentError("`normal_sign` should be either -1 or +1")) - m = new{F,I}( - size(vol)..., - Ref(vol), - zeros(F, 8), - zeros(Int, 3), - zeros(F, 3), - zeros(F, 3), - zeros(I, (3, size(vol)...)), - Triangle[], - Vertex[], - Normal[], - normal_sign, - Ref(x), - Ref(y), - Ref(z), - ) - sz = size(vol) |> prod - sizehint!(m.triangles, nextpow(2, sz ÷ 6)) - sizehint!(m.vertices, nextpow(2, sz ÷ 2)) - sizehint!(m.normals, nextpow(2, sz ÷ 2)) - m + using PrecompileTools + using StaticArrays + + import Base: RefValue + + export MC, march, march_legacy + + include("lut.jl") + + const Vertex{F} = SVector{3, F} + const Normal{F} = SVector{3, F} + const Triangle = SVector{3, Int} + + """ + MC(vol::Array{F,3}, I::Type = Int; x = [], y = [], z = []) + + # Description + Structure to hold temporaries and output of the marching cubes algorithm. + Vertices, normals and triangles are accessible using `m.vertices`, `m.normals` and `m.triangles`. + 1-based indexing is assumed. + Optionally, pass `x`, `y`, `z` to normalize the vertices coordinates instead of `0:nx-1`, `0:ny-1` and `0:nz-1`. + + # Arguments + - `vol`: rank 3 array of floats (volume) on which the Marching Cubes algorithm is applied. + - `I`: Int32, Int64, ...: vertices / normals / triangles index `Integer` type (defaults to `Int`). + - `x` normalize vertex.x to the coordinates of `x`. + - `y` normalize vertex.y to the coordinates of `y`. + - `z` normalize vertex.z to the coordinates of `z`. + """ + struct MC{F, I} + nx::Int # height of the grid + ny::Int # width of the grid + nz::Int # depth of the grid + vol::RefValue{Array{F, 3}} # implicit function values sampled on the grid + cube::MVector{8, F} # values of the implicit function on the active cube + tv::MVector{3, I} # triangle buffer + vtx::MVector{3, F} # vertex buffer + nrm::MVector{3, F} # normal buffer + vert_indices::Array{I, 4} # pre-computed vertex indices + triangles::Vector{Triangle} # output triangles + vertices::Vector{Vertex{F}} # output vertex positions + normals::Vector{Normal{F}} # output vertex normals + normal_sign::Int # direction of normal vectors (+1 for outward / -1 for inward) + x::RefValue{Vector{F}} + y::RefValue{Vector{F}} + z::RefValue{Vector{F}} + MC( + vol::Array{F, 3}, + I::Type{G} = Int; + normal_sign::Integer = 1, + x::AbstractVector{F} = F[], + y::AbstractVector{F} = F[], + z::AbstractVector{F} = F[], + ) where {F <: AbstractFloat, G <: Integer} = begin + isa(x, AbstractRange) && (x = collect(x)) + isa(y, AbstractRange) && (y = collect(y)) + isa(z, AbstractRange) && (z = collect(z)) + abs(normal_sign) == 1 || + throw(ArgumentError("`normal_sign` should be either -1 or +1")) + m = new{F, I}( + size(vol)..., + Ref(vol), + zeros(F, 8), + zeros(Int, 3), + zeros(F, 3), + zeros(F, 3), + zeros(I, (3, size(vol)...)), + Triangle[], + Vertex[], + Normal[], + normal_sign, + Ref(x), + Ref(y), + Ref(z), + ) + sz = size(vol) |> prod + sizehint!(m.triangles, nextpow(2, sz ÷ 6)) + sizehint!(m.vertices, nextpow(2, sz ÷ 2)) + sizehint!(m.normals, nextpow(2, sz ÷ 2)) + m + end end -end - -""" - lut_entry(vol, cb, i, j, k, iso) -# Description -Cube sign representation. -""" -lut_entry(vol::Array{F,3}, cb, i, j, k, iso) where {F} = begin - lut = UInt8(0) - @inbounds for p ∈ 0:7 - c = vol[i+((p⊻(p>>1))&1), j+((p>>1)&1), k+((p>>2)&1)] - iso - abs(c) < eps(F) && (c = eps(F)) - c > 0 && (lut += (UInt8(1) << p)) - cb[p+1] = c + """ + lut_entry(vol, cb, i, j, k, iso) + + # Description + Cube sign representation. + """ + lut_entry(vol::Array{F, 3}, cb, i, j, k, iso) where {F} = begin + lut = UInt8(0) + @inbounds for p in 0:7 + c = vol[i + ((p ⊻ (p >> 1)) & 1), j + ((p >> 1) & 1), k + ((p >> 2) & 1)] - iso + abs(c) < eps(F) && (c = eps(F)) + c > 0 && (lut += (UInt8(1) << p)) + cb[p + 1] = c + end + Int(lut) + 1 end - Int(lut) + 1 -end -""" - denormalize(m::MC) - -# Description -Optional denormalization step of vertices to user coordinates. -""" -denormalize(m::MC) = begin - if length(m.x[]) > 0 && length(m.y[]) > 0 && length(m.z[]) > 0 - mx, Mx = extrema(m.x[]) - my, My = extrema(m.y[]) - mz, Mz = extrema(m.z[]) - scl = - @SVector([Mx - mx, My - my, Mz - mz]) ./ - @SVector([m.nx - 1, m.ny - 1, m.nz - 1]) - off = @SVector([mx, my, mz]) - @inbounds for (n, v) in enumerate(m.vertices) - m.vertices[n] = Vertex(off .+ v .* scl) + """ + denormalize(m::MC) + + # Description + Optional denormalization step of vertices to user coordinates. + """ + denormalize(m::MC) = begin + if length(m.x[]) > 0 && length(m.y[]) > 0 && length(m.z[]) > 0 + mx, Mx = extrema(m.x[]) + my, My = extrema(m.y[]) + mz, Mz = extrema(m.z[]) + scl = + @SVector([Mx - mx, My - my, Mz - mz]) ./ + @SVector([m.nx - 1, m.ny - 1, m.nz - 1]) + off = @SVector([mx, my, mz]) + @inbounds for (n, v) in enumerate(m.vertices) + m.vertices[n] = Vertex(off .+ v .* scl) + end end + nothing end - nothing -end -""" - march_legacy(m::MC, isovalue::Number = 0) + """ + march_legacy(m::MC, isovalue::Number = 0) -# Description -Original Marching Cubes algorithm. + # Description + Original Marching Cubes algorithm. -# Arguments - - `m`: Marching Cubes data structure. - - `isovalue`: isosurface value. -""" -march_legacy(m::MC{F}, isovalue::Number = F(0)) where {F} = begin - empty!(m.triangles) - empty!(m.vertices) - empty!(m.normals) + # Arguments + - `m`: Marching Cubes data structure. + - `isovalue`: isosurface value. + """ + march_legacy(m::MC{F}, isovalue::Number = F(0)) where {F} = begin + empty!(m.triangles) + empty!(m.vertices) + empty!(m.normals) - vol = m.vol[] - iso = F(isovalue) + vol = m.vol[] + iso = F(isovalue) - compute_intersection_points(m, vol, m.cube, iso) + compute_intersection_points(m, vol, m.cube, iso) - @inbounds for k ∈ 1:m.nz-1, j ∈ 1:m.ny-1, i ∈ 1:m.nx-1 - lut = lut_entry(vol, m.cube, i, j, k, iso) - nt = 0 - while casesClassic[lut][3nt+1] > 0 - nt += 1 + @inbounds for k in 1:(m.nz - 1), j in 1:(m.ny - 1), i in 1:(m.nx - 1) + lut = lut_entry(vol, m.cube, i, j, k, iso) + nt = 0 + while casesClassic[lut][3nt + 1] > 0 + nt += 1 + end + add_triangle(m, i, j, k, casesClassic[lut], nt) end - add_triangle(m, i, j, k, casesClassic[lut], nt) - end - - denormalize(m) - nothing -end -""" - march(m::MC, isovalue::Number = 0) - -# Description -Topologically controlled Marching Cubes algorithm. + denormalize(m) + nothing + end -Arguments - - `m`: Marching Cubes data structure. - - `isovalue`: isosurface value. -""" -march(m::MC{F}, isovalue::Number = F(0)) where {F} = begin - empty!(m.triangles) - empty!(m.vertices) - empty!(m.normals) - - vol = m.vol[] - iso = F(isovalue) - - compute_intersection_points(m, vol, m.cube, iso) - - @inbounds for k ∈ 1:m.nz-1, j ∈ 1:m.ny-1, i ∈ 1:m.nx-1 - lut = lut_entry(vol, m.cube, i, j, k, iso) - - case = cases[lut][1] # case of the active cube in [1; 16] - cfg = cases[lut][2] # configuration of the active cube - subcfg = 1 # subconfiguration of the active cube - - if case == 1 - nothing - elseif case == 2 - add_triangle(m, i, j, k, tiling1[cfg], 1) - elseif case == 3 - add_triangle(m, i, j, k, tiling2[cfg], 2) - elseif case == 4 - if test_face(m.cube, test3[cfg]) - add_triangle(m, i, j, k, tiling3_2[cfg], 4) # 3.2 - else - add_triangle(m, i, j, k, tiling3_1[cfg], 2) # 3.1 - end - elseif case == 5 - if test_interior(case, m.cube, cfg, subcfg, test4[cfg]) - add_triangle(m, i, j, k, tiling4_1[cfg], 2) # 4.1.1 - else - add_triangle(m, i, j, k, tiling4_2[cfg], 6) # 4.1.2 - end - elseif case == 6 - add_triangle(m, i, j, k, tiling5[cfg], 3) - elseif case == 7 - if test_face(m.cube, test6[cfg][1]) - add_triangle(m, i, j, k, tiling6_2[cfg], 5) # 6.2 - else - if test_interior(case, m.cube, cfg, subcfg, test6[cfg][2]) - add_triangle(m, i, j, k, tiling6_1_1[cfg], 3) # 6.1.1 + """ + march(m::MC, isovalue::Number = 0) + + # Description + Topologically controlled Marching Cubes algorithm. + + Arguments + - `m`: Marching Cubes data structure. + - `isovalue`: isosurface value. + """ + march(m::MC{F}, isovalue::Number = F(0)) where {F} = begin + empty!(m.triangles) + empty!(m.vertices) + empty!(m.normals) + + vol = m.vol[] + iso = F(isovalue) + + compute_intersection_points(m, vol, m.cube, iso) + + @inbounds for k in 1:(m.nz - 1), j in 1:(m.ny - 1), i in 1:(m.nx - 1) + lut = lut_entry(vol, m.cube, i, j, k, iso) + + case = cases[lut][1] # case of the active cube in [1; 16] + cfg = cases[lut][2] # configuration of the active cube + subcfg = 1 # subconfiguration of the active cube + + if case == 1 + nothing + elseif case == 2 + add_triangle(m, i, j, k, tiling1[cfg], 1) + elseif case == 3 + add_triangle(m, i, j, k, tiling2[cfg], 2) + elseif case == 4 + if test_face(m.cube, test3[cfg]) + add_triangle(m, i, j, k, tiling3_2[cfg], 4) # 3.2 else - add_triangle(m, i, j, k, tiling6_1_2[cfg], 9, add_c_vertex(m, i, j, k)) # 6.1.2 + add_triangle(m, i, j, k, tiling3_1[cfg], 2) # 3.1 end - end - elseif case == 8 - test_face(m.cube, test7[cfg][1]) && (subcfg += 1) - test_face(m.cube, test7[cfg][2]) && (subcfg += 2) - test_face(m.cube, test7[cfg][3]) && (subcfg += 4) - if subcfg == 1 - add_triangle(m, i, j, k, tiling7_1[cfg], 3) - elseif subcfg == 2 || subcfg == 3 - add_triangle(m, i, j, k, tiling7_2[cfg][subcfg-1], 5) - elseif subcfg == 4 - add_triangle(m, i, j, k, tiling7_3[cfg][1], 9, add_c_vertex(m, i, j, k)) - elseif subcfg == 5 - add_triangle(m, i, j, k, tiling7_2[cfg][3], 5) - elseif subcfg == 6 || subcfg == 7 - add_triangle(m, i, j, k, tiling7_3[cfg][subcfg-4], 9, add_c_vertex(m, i, j, k)) - elseif subcfg == 8 - if test_interior(case, m.cube, cfg, subcfg, test7[cfg][4]) - add_triangle(m, i, j, k, tiling7_4_2[cfg], 9) + elseif case == 5 + if test_interior(case, m.cube, cfg, subcfg, test4[cfg]) + add_triangle(m, i, j, k, tiling4_1[cfg], 2) # 4.1.1 else - add_triangle(m, i, j, k, tiling7_4_1[cfg], 5) + add_triangle(m, i, j, k, tiling4_2[cfg], 6) # 4.1.2 end - end - elseif case == 9 - add_triangle(m, i, j, k, tiling8[cfg], 2) - elseif case == 10 - add_triangle(m, i, j, k, tiling9[cfg], 4) - elseif case == 11 - if test_face(m.cube, test10[cfg][1]) - if test_face(m.cube, test10[cfg][2]) - add_triangle(m, i, j, k, tiling10_1_1_[cfg], 4) # 10.1.1 + elseif case == 6 + add_triangle(m, i, j, k, tiling5[cfg], 3) + elseif case == 7 + if test_face(m.cube, test6[cfg][1]) + add_triangle(m, i, j, k, tiling6_2[cfg], 5) # 6.2 else - add_triangle(m, i, j, k, tiling10_2[cfg], 8, add_c_vertex(m, i, j, k)) # 10.2 + if test_interior(case, m.cube, cfg, subcfg, test6[cfg][2]) + add_triangle(m, i, j, k, tiling6_1_1[cfg], 3) # 6.1.1 + else + add_triangle(m, i, j, k, tiling6_1_2[cfg], 9, add_c_vertex(m, i, j, k)) # 6.1.2 + end end - else - if test_face(m.cube, test10[cfg][2]) - add_triangle(m, i, j, k, tiling10_2_[cfg], 8, add_c_vertex(m, i, j, k)) # 10.2 - else - if test_interior(case, m.cube, cfg, subcfg, test10[cfg][3]) - add_triangle(m, i, j, k, tiling10_1_1[cfg], 4) # 10.1.1 + elseif case == 8 + test_face(m.cube, test7[cfg][1]) && (subcfg += 1) + test_face(m.cube, test7[cfg][2]) && (subcfg += 2) + test_face(m.cube, test7[cfg][3]) && (subcfg += 4) + if subcfg == 1 + add_triangle(m, i, j, k, tiling7_1[cfg], 3) + elseif subcfg == 2 || subcfg == 3 + add_triangle(m, i, j, k, tiling7_2[cfg][subcfg - 1], 5) + elseif subcfg == 4 + add_triangle(m, i, j, k, tiling7_3[cfg][1], 9, add_c_vertex(m, i, j, k)) + elseif subcfg == 5 + add_triangle(m, i, j, k, tiling7_2[cfg][3], 5) + elseif subcfg == 6 || subcfg == 7 + add_triangle(m, i, j, k, tiling7_3[cfg][subcfg - 4], 9, add_c_vertex(m, i, j, k)) + elseif subcfg == 8 + if test_interior(case, m.cube, cfg, subcfg, test7[cfg][4]) + add_triangle(m, i, j, k, tiling7_4_2[cfg], 9) else - add_triangle(m, i, j, k, tiling10_1_2[cfg], 8) # 10.1.2 + add_triangle(m, i, j, k, tiling7_4_1[cfg], 5) end end - end - elseif case == 12 - add_triangle(m, i, j, k, tiling11[cfg], 4) - elseif case == 13 - if test_face(m.cube, test12[cfg][1]) - if test_face(m.cube, test12[cfg][2]) - add_triangle(m, i, j, k, tiling12_1_1_[cfg], 4) # 12.1.1 + elseif case == 9 + add_triangle(m, i, j, k, tiling8[cfg], 2) + elseif case == 10 + add_triangle(m, i, j, k, tiling9[cfg], 4) + elseif case == 11 + if test_face(m.cube, test10[cfg][1]) + if test_face(m.cube, test10[cfg][2]) + add_triangle(m, i, j, k, tiling10_1_1_[cfg], 4) # 10.1.1 + else + add_triangle(m, i, j, k, tiling10_2[cfg], 8, add_c_vertex(m, i, j, k)) # 10.2 + end else - add_triangle(m, i, j, k, tiling12_2[cfg], 8, add_c_vertex(m, i, j, k)) # 12.2 + if test_face(m.cube, test10[cfg][2]) + add_triangle(m, i, j, k, tiling10_2_[cfg], 8, add_c_vertex(m, i, j, k)) # 10.2 + else + if test_interior(case, m.cube, cfg, subcfg, test10[cfg][3]) + add_triangle(m, i, j, k, tiling10_1_1[cfg], 4) # 10.1.1 + else + add_triangle(m, i, j, k, tiling10_1_2[cfg], 8) # 10.1.2 + end + end end - else - if test_face(m.cube, test12[cfg][2]) - add_triangle(m, i, j, k, tiling12_2_[cfg], 8, add_c_vertex(m, i, j, k)) # 12.2 + elseif case == 12 + add_triangle(m, i, j, k, tiling11[cfg], 4) + elseif case == 13 + if test_face(m.cube, test12[cfg][1]) + if test_face(m.cube, test12[cfg][2]) + add_triangle(m, i, j, k, tiling12_1_1_[cfg], 4) # 12.1.1 + else + add_triangle(m, i, j, k, tiling12_2[cfg], 8, add_c_vertex(m, i, j, k)) # 12.2 + end else - if test_interior(case, m.cube, cfg, subcfg, test12[cfg][3]) - add_triangle(m, i, j, k, tiling12_1_1[cfg], 4) # 12.1.1 + if test_face(m.cube, test12[cfg][2]) + add_triangle(m, i, j, k, tiling12_2_[cfg], 8, add_c_vertex(m, i, j, k)) # 12.2 else - add_triangle(m, i, j, k, tiling12_1_2[cfg], 8) # 12.1.2 + if test_interior(case, m.cube, cfg, subcfg, test12[cfg][3]) + add_triangle(m, i, j, k, tiling12_1_1[cfg], 4) # 12.1.1 + else + add_triangle(m, i, j, k, tiling12_1_2[cfg], 8) # 12.1.2 + end end end - end - elseif case == 14 - test_face(m.cube, test13[cfg][1]) && (subcfg += 1) - test_face(m.cube, test13[cfg][2]) && (subcfg += 2) - test_face(m.cube, test13[cfg][3]) && (subcfg += 4) - test_face(m.cube, test13[cfg][4]) && (subcfg += 8) - test_face(m.cube, test13[cfg][5]) && (subcfg += 16) - test_face(m.cube, test13[cfg][6]) && (subcfg += 32) - if (sc = subcfg13[subcfg]) == 0 # 13.1 - add_triangle(m, i, j, k, tiling13_1[cfg], 4) - elseif sc ≥ 1 && sc ≤ 6 # 13.2 - add_triangle(m, i, j, k, tiling13_2[cfg][sc], 6) - elseif sc ≥ 7 && sc ≤ 18 # 13.3 - add_triangle(m, i, j, k, tiling13_3[cfg][sc-6], 10, add_c_vertex(m, i, j, k)) - elseif sc ≥ 19 && sc ≤ 22 # 13.4 - add_triangle(m, i, j, k, tiling13_4[cfg][sc-18], 12, add_c_vertex(m, i, j, k)) - elseif sc ≥ 23 && sc ≤ 26 # 13.5 - subcfg = sc - 22 - if test_interior(case, m.cube, cfg, subcfg, test13[cfg][6]) - add_triangle(m, i, j, k, tiling13_5_1[cfg][subcfg], 6) - else - add_triangle(m, i, j, k, tiling13_5_2[cfg][subcfg], 10) + elseif case == 14 + test_face(m.cube, test13[cfg][1]) && (subcfg += 1) + test_face(m.cube, test13[cfg][2]) && (subcfg += 2) + test_face(m.cube, test13[cfg][3]) && (subcfg += 4) + test_face(m.cube, test13[cfg][4]) && (subcfg += 8) + test_face(m.cube, test13[cfg][5]) && (subcfg += 16) + test_face(m.cube, test13[cfg][6]) && (subcfg += 32) + if (sc = subcfg13[subcfg]) == 0 # 13.1 + add_triangle(m, i, j, k, tiling13_1[cfg], 4) + elseif sc ≥ 1 && sc ≤ 6 # 13.2 + add_triangle(m, i, j, k, tiling13_2[cfg][sc], 6) + elseif sc ≥ 7 && sc ≤ 18 # 13.3 + add_triangle(m, i, j, k, tiling13_3[cfg][sc - 6], 10, add_c_vertex(m, i, j, k)) + elseif sc ≥ 19 && sc ≤ 22 # 13.4 + add_triangle(m, i, j, k, tiling13_4[cfg][sc - 18], 12, add_c_vertex(m, i, j, k)) + elseif sc ≥ 23 && sc ≤ 26 # 13.5 + subcfg = sc - 22 + if test_interior(case, m.cube, cfg, subcfg, test13[cfg][6]) + add_triangle(m, i, j, k, tiling13_5_1[cfg][subcfg], 6) + else + add_triangle(m, i, j, k, tiling13_5_2[cfg][subcfg], 10) + end + elseif sc ≥ 27 && sc ≤ 38 # 13.3 + add_triangle(m, i, j, k, tiling13_3_[cfg][sc - 26], 10, add_c_vertex(m, i, j, k)) + elseif sc ≥ 39 && sc ≤ 44 # 13.2 + add_triangle(m, i, j, k, tiling13_2_[cfg][sc - 38], 6) + elseif sc == 45 # 13.1 + add_triangle(m, i, j, k, tiling13_1_[cfg], 4) + elseif sc > 45 + @error "Impossible case 13 ?" end - elseif sc ≥ 27 && sc ≤ 38 # 13.3 - add_triangle(m, i, j, k, tiling13_3_[cfg][sc-26], 10, add_c_vertex(m, i, j, k)) - elseif sc ≥ 39 && sc ≤ 44 # 13.2 - add_triangle(m, i, j, k, tiling13_2_[cfg][sc-38], 6) - elseif sc == 45 # 13.1 - add_triangle(m, i, j, k, tiling13_1_[cfg], 4) - elseif sc > 45 - @error "Impossible case 13 ?" + elseif case == 15 + add_triangle(m, i, j, k, tiling14[cfg], 4) end - elseif case == 15 - add_triangle(m, i, j, k, tiling14[cfg], 4) end - end - denormalize(m) - nothing -end - -""" - compute_intersection_points(m::MC, vol, cb, iso) - -# Description -Computes almost all the vertices of the mesh by interpolation along the cubes edges. -""" -compute_intersection_points(m::MC{F}, vol, cb, iso) where {F} = begin - @inbounds for k ∈ axes(vol, 3), j ∈ axes(vol, 2), i ∈ axes(vol, 1) - c0::F = vol[i, j, k] - iso - c1::F = i < m.nx ? vol[i+1, j, k] - iso : c0 - c2::F = j < m.ny ? vol[i, j+1, k] - iso : c0 - c3::F = k < m.nz ? vol[i, j, k+1] - iso : c0 - m.cube[1] = abs(c0) < eps(F) ? eps(F) : c0 - m.cube[2] = abs(c1) < eps(F) ? eps(F) : c1 - m.cube[4] = abs(c2) < eps(F) ? eps(F) : c2 - m.cube[5] = abs(c3) < eps(F) ? eps(F) : c3 - if m.cube[1] < 0 - m.cube[2] > 0 && - (m.vert_indices[1, i, j, k] = add_x_vertex(m, vol, cb, i, j, k)) - m.cube[4] > 0 && - (m.vert_indices[2, i, j, k] = add_y_vertex(m, vol, cb, i, j, k)) - m.cube[5] > 0 && - (m.vert_indices[3, i, j, k] = add_z_vertex(m, vol, cb, i, j, k)) - else - m.cube[2] < 0 && - (m.vert_indices[1, i, j, k] = add_x_vertex(m, vol, cb, i, j, k)) - m.cube[4] < 0 && - (m.vert_indices[2, i, j, k] = add_y_vertex(m, vol, cb, i, j, k)) - m.cube[5] < 0 && - (m.vert_indices[3, i, j, k] = add_z_vertex(m, vol, cb, i, j, k)) - end + denormalize(m) + nothing end - nothing -end -""" - test_interior(case, cb, cfg, subcfg, s) - -# Description -Tests if the components of the tessellation of the cube should be connected through the interior of the cube. -""" -test_interior(case, cb::MVector{N,T}, cfg, subcfg, s) where {N,T} = begin - test = 0 - At = Bt = Ct = Dt = T(0) - @inbounds begin - if case == 5 || case == 11 - a = (cb[5] - cb[1]) * (cb[7] - cb[3]) - (cb[8] - cb[4]) * (cb[6] - cb[2]) - b = ( - cb[3] * (cb[5] - cb[1]) + cb[1] * (cb[7] - cb[3]) - - cb[2] * (cb[8] - cb[4]) - cb[4] * (cb[6] - cb[2]) - ) - t = -b / 2a - (t < 0 || t > 1) && return s > 0 - - At = cb[1] + (cb[5] - cb[1]) * t - Bt = cb[4] + (cb[8] - cb[4]) * t - Ct = cb[3] + (cb[7] - cb[3]) * t - Dt = cb[2] + (cb[6] - cb[2]) * t - elseif case == 7 || case == 8 || case == 13 || case == 14 - # reference edge of the triangulation - edge = if case == 7 - test6[cfg][3] - elseif case == 8 - test7[cfg][5] - elseif case == 13 - test12[cfg][4] - elseif case == 14 - tiling13_5_1[cfg][subcfg][1] - 1 + """ + compute_intersection_points(m::MC, vol, cb, iso) + + # Description + Computes almost all the vertices of the mesh by interpolation along the cubes edges. + """ + compute_intersection_points(m::MC{F}, vol, cb, iso) where {F} = begin + @inbounds for k in axes(vol, 3), j in axes(vol, 2), i in axes(vol, 1) + c0::F = vol[i, j, k] - iso + c1::F = i < m.nx ? vol[i + 1, j, k] - iso : c0 + c2::F = j < m.ny ? vol[i, j + 1, k] - iso : c0 + c3::F = k < m.nz ? vol[i, j, k + 1] - iso : c0 + m.cube[1] = abs(c0) < eps(F) ? eps(F) : c0 + m.cube[2] = abs(c1) < eps(F) ? eps(F) : c1 + m.cube[4] = abs(c2) < eps(F) ? eps(F) : c2 + m.cube[5] = abs(c3) < eps(F) ? eps(F) : c3 + if m.cube[1] < 0 + m.cube[2] > 0 && + (m.vert_indices[1, i, j, k] = add_x_vertex(m, vol, cb, i, j, k)) + m.cube[4] > 0 && + (m.vert_indices[2, i, j, k] = add_y_vertex(m, vol, cb, i, j, k)) + m.cube[5] > 0 && + (m.vert_indices[3, i, j, k] = add_z_vertex(m, vol, cb, i, j, k)) else - -1 + m.cube[2] < 0 && + (m.vert_indices[1, i, j, k] = add_x_vertex(m, vol, cb, i, j, k)) + m.cube[4] < 0 && + (m.vert_indices[2, i, j, k] = add_y_vertex(m, vol, cb, i, j, k)) + m.cube[5] < 0 && + (m.vert_indices[3, i, j, k] = add_z_vertex(m, vol, cb, i, j, k)) end - if edge == 0 - t = cb[1] / (cb[1] - cb[2]) - Bt = cb[4] + (cb[3] - cb[4]) * t - Ct = cb[8] + (cb[7] - cb[8]) * t - Dt = cb[5] + (cb[6] - cb[5]) * t - elseif edge == 1 - t = cb[2] / (cb[2] - cb[3]) - Bt = cb[1] + (cb[4] - cb[1]) * t - Ct = cb[5] + (cb[8] - cb[5]) * t - Dt = cb[6] + (cb[7] - cb[6]) * t - elseif edge == 2 - t = cb[3] / (cb[3] - cb[4]) - Bt = cb[2] + (cb[1] - cb[2]) * t - Ct = cb[6] + (cb[5] - cb[6]) * t - Dt = cb[7] + (cb[8] - cb[7]) * t - elseif edge == 3 - t = cb[4] / (cb[4] - cb[1]) - Bt = cb[3] + (cb[2] - cb[3]) * t - Ct = cb[7] + (cb[6] - cb[7]) * t - Dt = cb[8] + (cb[5] - cb[8]) * t - elseif edge == 4 - t = cb[5] / (cb[5] - cb[6]) - Bt = cb[8] + (cb[7] - cb[8]) * t - Ct = cb[4] + (cb[3] - cb[4]) * t - Dt = cb[1] + (cb[2] - cb[1]) * t - elseif edge == 5 - t = cb[6] / (cb[6] - cb[7]) - Bt = cb[5] + (cb[8] - cb[5]) * t - Ct = cb[1] + (cb[4] - cb[1]) * t - Dt = cb[2] + (cb[3] - cb[2]) * t - elseif edge == 6 - t = cb[7] / (cb[7] - cb[8]) - Bt = cb[6] + (cb[5] - cb[6]) * t - Ct = cb[2] + (cb[1] - cb[2]) * t - Dt = cb[3] + (cb[4] - cb[3]) * t - elseif edge == 7 - t = cb[8] / (cb[8] - cb[5]) - Bt = cb[7] + (cb[6] - cb[7]) * t - Ct = cb[3] + (cb[2] - cb[3]) * t - Dt = cb[4] + (cb[1] - cb[4]) * t - elseif edge == 8 - t = cb[1] / (cb[1] - cb[5]) + end + nothing + end + + """ + test_interior(case, cb, cfg, subcfg, s) + + # Description + Tests if the components of the tessellation of the cube should be connected through the interior of the cube. + """ + test_interior(case, cb::MVector{N, T}, cfg, subcfg, s) where {N, T} = begin + test = 0 + At = Bt = Ct = Dt = T(0) + @inbounds begin + if case == 5 || case == 11 + a = (cb[5] - cb[1]) * (cb[7] - cb[3]) - (cb[8] - cb[4]) * (cb[6] - cb[2]) + b = ( + cb[3] * (cb[5] - cb[1]) + cb[1] * (cb[7] - cb[3]) - + cb[2] * (cb[8] - cb[4]) - cb[4] * (cb[6] - cb[2]) + ) + t = -b / 2a + (t < 0 || t > 1) && return s > 0 + + At = cb[1] + (cb[5] - cb[1]) * t Bt = cb[4] + (cb[8] - cb[4]) * t Ct = cb[3] + (cb[7] - cb[3]) * t Dt = cb[2] + (cb[6] - cb[2]) * t - elseif edge == 9 - t = cb[2] / (cb[2] - cb[6]) - Bt = cb[1] + (cb[5] - cb[1]) * t - Ct = cb[4] + (cb[8] - cb[4]) * t - Dt = cb[3] + (cb[7] - cb[3]) * t - elseif edge == 10 - t = cb[3] / (cb[3] - cb[7]) - Bt = cb[2] + (cb[6] - cb[2]) * t - Ct = cb[1] + (cb[5] - cb[1]) * t - Dt = cb[4] + (cb[8] - cb[4]) * t - elseif edge == 11 - t = cb[4] / (cb[4] - cb[8]) - Bt = cb[3] + (cb[7] - cb[3]) * t - Ct = cb[2] + (cb[6] - cb[2]) * t - Dt = cb[1] + (cb[5] - cb[1]) * t + elseif case == 7 || case == 8 || case == 13 || case == 14 + # reference edge of the triangulation + edge = if case == 7 + test6[cfg][3] + elseif case == 8 + test7[cfg][5] + elseif case == 13 + test12[cfg][4] + elseif case == 14 + tiling13_5_1[cfg][subcfg][1] - 1 + else + -1 + end + if edge == 0 + t = cb[1] / (cb[1] - cb[2]) + Bt = cb[4] + (cb[3] - cb[4]) * t + Ct = cb[8] + (cb[7] - cb[8]) * t + Dt = cb[5] + (cb[6] - cb[5]) * t + elseif edge == 1 + t = cb[2] / (cb[2] - cb[3]) + Bt = cb[1] + (cb[4] - cb[1]) * t + Ct = cb[5] + (cb[8] - cb[5]) * t + Dt = cb[6] + (cb[7] - cb[6]) * t + elseif edge == 2 + t = cb[3] / (cb[3] - cb[4]) + Bt = cb[2] + (cb[1] - cb[2]) * t + Ct = cb[6] + (cb[5] - cb[6]) * t + Dt = cb[7] + (cb[8] - cb[7]) * t + elseif edge == 3 + t = cb[4] / (cb[4] - cb[1]) + Bt = cb[3] + (cb[2] - cb[3]) * t + Ct = cb[7] + (cb[6] - cb[7]) * t + Dt = cb[8] + (cb[5] - cb[8]) * t + elseif edge == 4 + t = cb[5] / (cb[5] - cb[6]) + Bt = cb[8] + (cb[7] - cb[8]) * t + Ct = cb[4] + (cb[3] - cb[4]) * t + Dt = cb[1] + (cb[2] - cb[1]) * t + elseif edge == 5 + t = cb[6] / (cb[6] - cb[7]) + Bt = cb[5] + (cb[8] - cb[5]) * t + Ct = cb[1] + (cb[4] - cb[1]) * t + Dt = cb[2] + (cb[3] - cb[2]) * t + elseif edge == 6 + t = cb[7] / (cb[7] - cb[8]) + Bt = cb[6] + (cb[5] - cb[6]) * t + Ct = cb[2] + (cb[1] - cb[2]) * t + Dt = cb[3] + (cb[4] - cb[3]) * t + elseif edge == 7 + t = cb[8] / (cb[8] - cb[5]) + Bt = cb[7] + (cb[6] - cb[7]) * t + Ct = cb[3] + (cb[2] - cb[3]) * t + Dt = cb[4] + (cb[1] - cb[4]) * t + elseif edge == 8 + t = cb[1] / (cb[1] - cb[5]) + Bt = cb[4] + (cb[8] - cb[4]) * t + Ct = cb[3] + (cb[7] - cb[3]) * t + Dt = cb[2] + (cb[6] - cb[2]) * t + elseif edge == 9 + t = cb[2] / (cb[2] - cb[6]) + Bt = cb[1] + (cb[5] - cb[1]) * t + Ct = cb[4] + (cb[8] - cb[4]) * t + Dt = cb[3] + (cb[7] - cb[3]) * t + elseif edge == 10 + t = cb[3] / (cb[3] - cb[7]) + Bt = cb[2] + (cb[6] - cb[2]) * t + Ct = cb[1] + (cb[5] - cb[1]) * t + Dt = cb[4] + (cb[8] - cb[4]) * t + elseif edge == 11 + t = cb[4] / (cb[4] - cb[8]) + Bt = cb[3] + (cb[7] - cb[3]) * t + Ct = cb[2] + (cb[6] - cb[2]) * t + Dt = cb[1] + (cb[5] - cb[1]) * t + else + @error "Invalid edge $edge" + end else - @error "Invalid edge $edge" + @error "Invalid ambiguous case $case" end - else - @error "Invalid ambiguous case $case" end - end - At ≥ 0 && (test += 1) - Bt ≥ 0 && (test += 2) - Ct ≥ 0 && (test += 4) - Dt ≥ 0 && (test += 8) - - if test == 6 || test == 8 || test == 9 || test == 12 || (test ≥ 0 && test ≤ 4) - s > 0 - elseif test == 5 - (At * Ct - Bt * Dt < eps(T)) ? s > 0 : s < 0 - elseif test == 10 - (At * Ct - Bt * Dt ≥ eps(T)) ? s > 0 : s < 0 - else # test == 7 || test == 11 || test ≥ 13 - s < 0 - end::Bool -end - -""" - test_face(cb, face) + At ≥ 0 && (test += 1) + Bt ≥ 0 && (test += 2) + Ct ≥ 0 && (test += 4) + Dt ≥ 0 && (test += 8) + + if test == 6 || test == 8 || test == 9 || test == 12 || (test ≥ 0 && test ≤ 4) + s > 0 + elseif test == 5 + (At * Ct - Bt * Dt < eps(T)) ? s > 0 : s < 0 + elseif test == 10 + (At * Ct - Bt * Dt ≥ eps(T)) ? s > 0 : s < 0 + else # test == 7 || test == 11 || test ≥ 13 + s < 0 + end::Bool + end -# Description -Tests if the components of the tessellation of the cube should be connected by the interior of an ambiguous face. -""" -test_face(cb::MVector{N,T}, face) where {N,T} = begin - @inbounds begin - if (af = abs(face)) == 1 - A = cb[1] - B = cb[5] - C = cb[6] - D = cb[2] - elseif af == 2 - A = cb[2] - B = cb[6] - C = cb[7] - D = cb[3] - elseif af == 3 - A = cb[3] - B = cb[7] - C = cb[8] - D = cb[4] - elseif af == 4 - A = cb[4] - B = cb[8] - C = cb[5] - D = cb[1] - elseif af == 5 - A = cb[1] - B = cb[4] - C = cb[3] - D = cb[2] - elseif af == 6 - A = cb[5] - B = cb[8] - C = cb[7] - D = cb[6] - else - @error "Invalid face code $face" + """ + test_face(cb, face) + + # Description + Tests if the components of the tessellation of the cube should be connected by the interior of an ambiguous face. + """ + test_face(cb::MVector{N, T}, face) where {N, T} = begin + @inbounds begin + if (af = abs(face)) == 1 + A = cb[1] + B = cb[5] + C = cb[6] + D = cb[2] + elseif af == 2 + A = cb[2] + B = cb[6] + C = cb[7] + D = cb[3] + elseif af == 3 + A = cb[3] + B = cb[7] + C = cb[8] + D = cb[4] + elseif af == 4 + A = cb[4] + B = cb[8] + C = cb[5] + D = cb[1] + elseif af == 5 + A = cb[1] + B = cb[4] + C = cb[3] + D = cb[2] + elseif af == 6 + A = cb[5] + B = cb[8] + C = cb[7] + D = cb[6] + else + @error "Invalid face code $face" + end end + abs(A * C - B * D) < eps(T) ? face ≥ 0 : face * A * (A * C - B * D) ≥ 0 # face and A invert signs end - abs(A * C - B * D) < eps(T) ? face ≥ 0 : face * A * (A * C - B * D) ≥ 0 # face and A invert signs -end - -""" - add_triangle(m::MC, i, j, k, tri, n, v12 = 0) - -# Description -Routine to add a triangle to the mesh. -# Arguments - - `tri` the code for the triangle as a sequence of edges index. - - `n` the number of triangles to produce. - - `v12` the index of the interior vertex to use, if necessary. -""" -add_triangle(m::MC{F,I}, i, j, k, tri, n, v12 = I(0)) where {F,I} = begin - @inbounds for t ∈ 1:3n - tr = tri[t] - id = (t - 1) % 3 + 1 - m.tv[id] = tv = if tr == 1 - m.vert_indices[1, i, j, k] - elseif tr == 2 - m.vert_indices[2, i+1, j, k] - elseif tr == 3 - m.vert_indices[1, i, j+1, k] - elseif tr == 4 - m.vert_indices[2, i, j, k] - elseif tr == 5 - m.vert_indices[1, i, j, k+1] - elseif tr == 6 - m.vert_indices[2, i+1, j, k+1] - elseif tr == 7 - m.vert_indices[1, i, j+1, k+1] - elseif tr == 8 - m.vert_indices[2, i, j, k+1] - elseif tr == 9 - m.vert_indices[3, i, j, k] - elseif tr == 10 - m.vert_indices[3, i+1, j, k] - elseif tr == 11 - m.vert_indices[3, i+1, j+1, k] - elseif tr == 12 - m.vert_indices[3, i, j+1, k] - else - I(v12) - end::I - tv == 0 && @warn "Invalid triangle $(length(m.triangles) + 1)" - id == 3 && push!(m.triangles, Triangle(m.tv)) + """ + add_triangle(m::MC, i, j, k, tri, n, v12 = 0) + + # Description + Routine to add a triangle to the mesh. + + # Arguments + - `tri` the code for the triangle as a sequence of edges index. + - `n` the number of triangles to produce. + - `v12` the index of the interior vertex to use, if necessary. + """ + add_triangle(m::MC{F, I}, i, j, k, tri, n, v12 = I(0)) where {F, I} = begin + @inbounds for t in 1:3n + tr = tri[t] + id = (t - 1) % 3 + 1 + m.tv[id] = tv = if tr == 1 + m.vert_indices[1, i, j, k] + elseif tr == 2 + m.vert_indices[2, i + 1, j, k] + elseif tr == 3 + m.vert_indices[1, i, j + 1, k] + elseif tr == 4 + m.vert_indices[2, i, j, k] + elseif tr == 5 + m.vert_indices[1, i, j, k + 1] + elseif tr == 6 + m.vert_indices[2, i + 1, j, k + 1] + elseif tr == 7 + m.vert_indices[1, i, j + 1, k + 1] + elseif tr == 8 + m.vert_indices[2, i, j, k + 1] + elseif tr == 9 + m.vert_indices[3, i, j, k] + elseif tr == 10 + m.vert_indices[3, i + 1, j, k] + elseif tr == 11 + m.vert_indices[3, i + 1, j + 1, k] + elseif tr == 12 + m.vert_indices[3, i, j + 1, k] + else + I(v12) + end::I + tv == 0 && @warn "Invalid triangle $(length(m.triangles) + 1)" + id == 3 && push!(m.triangles, Triangle(m.tv)) + end + nothing end - nothing -end - -""" - ∇x(args...; kwargs...) - -# Description -Interpolates the horizontal gradient of the implicit function at the lower vertex of the specified cube. -""" -@inline ∇x(vol::AbstractArray{F}, i, j, k, nx) where {F} = @inbounds( - i > 1 ? - (i < nx ? (vol[i+1, j, k] - vol[i-1, j, k]) / 2 : (vol[i, j, k] - vol[i-1, j, k])) : - (vol[i+1, j, k] - vol[i, j, k]) -)::F - -@inline ∇y(vol::AbstractArray{F}, i, j, k, ny) where {F} = @inbounds( - j > 1 ? - (j < ny ? (vol[i, j+1, k] - vol[i, j-1, k]) / 2 : (vol[i, j, k] - vol[i, j-1, k])) : - (vol[i, j+1, k] - vol[i, j, k]) -)::F - -@inline ∇z(vol::AbstractArray{F}, i, j, k, nz) where {F} = @inbounds( - k > 1 ? - (k < nz ? (vol[i, j, k+1] - vol[i, j, k-1]) / 2 : (vol[i, j, k] - vol[i, j, k-1])) : - (vol[i, j, k+1] - vol[i, j, k]) -)::F - -@inline norm(x) = @inbounds √(x[1]^2 + x[2]^2 + x[3]^2) - -""" - add_x_vertex(m::MC, vol, cb, i, j, k) -# Description -Adds a vertex on the current horizontal edge. -""" -add_x_vertex(m::MC{F}, vol, cb, i, j, k) where {F} = begin - @inbounds begin - u = cb[1] / (cb[1] - cb[2]) - m.nrm[1] = (1 - u) * ∇x(vol, i, j, k, m.nx) + u * ∇x(vol, i + 1, j, k, m.nx) - m.nrm[2] = (1 - u) * ∇y(vol, i, j, k, m.ny) + u * ∇y(vol, i + 1, j, k, m.ny) - m.nrm[3] = (1 - u) * ∇z(vol, i, j, k, m.nz) + u * ∇z(vol, i + 1, j, k, m.nz) - (mag = norm(m.nrm)) > eps(F) && (m.nrm ./= mag) + """ + ∇x(args...; kwargs...) + + # Description + Interpolates the horizontal gradient of the implicit function at the lower vertex of the specified cube. + """ + @inline ∇x(vol::AbstractArray{F}, i, j, k, nx) where {F} = @inbounds( + i > 1 ? + (i < nx ? (vol[i + 1, j, k] - vol[i - 1, j, k]) / 2 : (vol[i, j, k] - vol[i - 1, j, k])) : + (vol[i + 1, j, k] - vol[i, j, k]) + )::F + + @inline ∇y(vol::AbstractArray{F}, i, j, k, ny) where {F} = @inbounds( + j > 1 ? + (j < ny ? (vol[i, j + 1, k] - vol[i, j - 1, k]) / 2 : (vol[i, j, k] - vol[i, j - 1, k])) : + (vol[i, j + 1, k] - vol[i, j, k]) + )::F + + @inline ∇z(vol::AbstractArray{F}, i, j, k, nz) where {F} = @inbounds( + k > 1 ? + (k < nz ? (vol[i, j, k + 1] - vol[i, j, k - 1]) / 2 : (vol[i, j, k] - vol[i, j, k - 1])) : + (vol[i, j, k + 1] - vol[i, j, k]) + )::F + + @inline norm(x) = @inbounds √(x[1]^2 + x[2]^2 + x[3]^2) + + """ + add_x_vertex(m::MC, vol, cb, i, j, k) + + # Description + Adds a vertex on the current horizontal edge. + """ + add_x_vertex(m::MC{F}, vol, cb, i, j, k) where {F} = begin + @inbounds begin + u = cb[1] / (cb[1] - cb[2]) + m.nrm[1] = (1 - u) * ∇x(vol, i, j, k, m.nx) + u * ∇x(vol, i + 1, j, k, m.nx) + m.nrm[2] = (1 - u) * ∇y(vol, i, j, k, m.ny) + u * ∇y(vol, i + 1, j, k, m.ny) + m.nrm[3] = (1 - u) * ∇z(vol, i, j, k, m.nz) + u * ∇z(vol, i + 1, j, k, m.nz) + (mag = norm(m.nrm)) > eps(F) && (m.nrm ./= mag) + end + push!(m.vertices, Vertex(i - 1 + u, j - 1, k - 1)) + push!(m.normals, Normal(m.normal_sign * m.nrm)) + length(m.vertices) end - push!(m.vertices, Vertex(i - 1 + u, j - 1, k - 1)) - push!(m.normals, Normal(m.normal_sign * m.nrm)) - length(m.vertices) -end -add_y_vertex(m::MC{F}, vol, cb, i, j, k) where {F} = begin - @inbounds begin - u = cb[1] / (cb[1] - cb[4]) - m.nrm[1] = (1 - u) * ∇x(vol, i, j, k, m.nx) + u * ∇x(vol, i, j + 1, k, m.nx) - m.nrm[2] = (1 - u) * ∇y(vol, i, j, k, m.ny) + u * ∇y(vol, i, j + 1, k, m.ny) - m.nrm[3] = (1 - u) * ∇z(vol, i, j, k, m.nz) + u * ∇z(vol, i, j + 1, k, m.nz) - (mag = norm(m.nrm)) > eps(F) && (m.nrm ./= mag) + add_y_vertex(m::MC{F}, vol, cb, i, j, k) where {F} = begin + @inbounds begin + u = cb[1] / (cb[1] - cb[4]) + m.nrm[1] = (1 - u) * ∇x(vol, i, j, k, m.nx) + u * ∇x(vol, i, j + 1, k, m.nx) + m.nrm[2] = (1 - u) * ∇y(vol, i, j, k, m.ny) + u * ∇y(vol, i, j + 1, k, m.ny) + m.nrm[3] = (1 - u) * ∇z(vol, i, j, k, m.nz) + u * ∇z(vol, i, j + 1, k, m.nz) + (mag = norm(m.nrm)) > eps(F) && (m.nrm ./= mag) + end + push!(m.vertices, Vertex(i - 1, j - 1 + u, k - 1)) + push!(m.normals, Normal(m.normal_sign * m.nrm)) + length(m.vertices) end - push!(m.vertices, Vertex(i - 1, j - 1 + u, k - 1)) - push!(m.normals, Normal(m.normal_sign * m.nrm)) - length(m.vertices) -end -add_z_vertex(m::MC{F}, vol, cb, i, j, k) where {F} = begin - @inbounds begin - u = cb[1] / (cb[1] - cb[5]) - m.nrm[1] = (1 - u) * ∇x(vol, i, j, k, m.nx) + u * ∇x(vol, i, j, k + 1, m.nx) - m.nrm[2] = (1 - u) * ∇y(vol, i, j, k, m.ny) + u * ∇y(vol, i, j, k + 1, m.ny) - m.nrm[3] = (1 - u) * ∇z(vol, i, j, k, m.nz) + u * ∇z(vol, i, j, k + 1, m.nz) - (mag = norm(m.nrm)) > eps(F) && (m.nrm ./= mag) + add_z_vertex(m::MC{F}, vol, cb, i, j, k) where {F} = begin + @inbounds begin + u = cb[1] / (cb[1] - cb[5]) + m.nrm[1] = (1 - u) * ∇x(vol, i, j, k, m.nx) + u * ∇x(vol, i, j, k + 1, m.nx) + m.nrm[2] = (1 - u) * ∇y(vol, i, j, k, m.ny) + u * ∇y(vol, i, j, k + 1, m.ny) + m.nrm[3] = (1 - u) * ∇z(vol, i, j, k, m.nz) + u * ∇z(vol, i, j, k + 1, m.nz) + (mag = norm(m.nrm)) > eps(F) && (m.nrm ./= mag) + end + push!(m.vertices, Vertex(i - 1, j - 1, k - 1 + u)) + push!(m.normals, Normal(m.normal_sign * m.nrm)) + length(m.vertices) end - push!(m.vertices, Vertex(i - 1, j - 1, k - 1 + u)) - push!(m.normals, Normal(m.normal_sign * m.nrm)) - length(m.vertices) -end -""" - add_c_vertex(m::MC, i, j, k) - -# Description -Add a vertex inside the current cube. -""" -add_c_vertex(m::MC{F}, i, j, k) where {F} = begin - u = 0 - @inbounds begin - m.vtx .= 0 - m.nrm .= 0 - # compute the average of the intersection points of the cube - (n = m.vert_indices[1, i, j, k]) > 0 && - (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) - (n = m.vert_indices[2, i+1, j, k]) > 0 && - (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) - (n = m.vert_indices[1, i, j+1, k]) > 0 && - (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) - (n = m.vert_indices[2, i, j, k]) > 0 && - (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) - (n = m.vert_indices[1, i, j, k+1]) > 0 && - (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) - (n = m.vert_indices[2, i+1, j, k+1]) > 0 && - (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) - (n = m.vert_indices[1, i, j+1, k+1]) > 0 && - (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) - (n = m.vert_indices[2, i, j, k+1]) > 0 && - (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) - (n = m.vert_indices[3, i, j, k]) > 0 && - (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) - (n = m.vert_indices[3, i+1, j, k]) > 0 && - (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) - (n = m.vert_indices[3, i+1, j+1, k]) > 0 && - (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) - (n = m.vert_indices[3, i, j+1, k]) > 0 && - (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) - m.vtx ./= u - (mag = norm(m.nrm)) > eps(F) && (m.nrm ./= mag) + """ + add_c_vertex(m::MC, i, j, k) + + # Description + Add a vertex inside the current cube. + """ + add_c_vertex(m::MC{F}, i, j, k) where {F} = begin + u = 0 + @inbounds begin + m.vtx .= 0 + m.nrm .= 0 + # compute the average of the intersection points of the cube + (n = m.vert_indices[1, i, j, k]) > 0 && + (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) + (n = m.vert_indices[2, i + 1, j, k]) > 0 && + (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) + (n = m.vert_indices[1, i, j + 1, k]) > 0 && + (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) + (n = m.vert_indices[2, i, j, k]) > 0 && + (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) + (n = m.vert_indices[1, i, j, k + 1]) > 0 && + (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) + (n = m.vert_indices[2, i + 1, j, k + 1]) > 0 && + (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) + (n = m.vert_indices[1, i, j + 1, k + 1]) > 0 && + (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) + (n = m.vert_indices[2, i, j, k + 1]) > 0 && + (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) + (n = m.vert_indices[3, i, j, k]) > 0 && + (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) + (n = m.vert_indices[3, i + 1, j, k]) > 0 && + (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) + (n = m.vert_indices[3, i + 1, j + 1, k]) > 0 && + (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) + (n = m.vert_indices[3, i, j + 1, k]) > 0 && + (u += 1; m.vtx .+= m.vertices[n]; m.nrm .+= m.normals[n]) + m.vtx ./= u + (mag = norm(m.nrm)) > eps(F) && (m.nrm ./= mag) + end + push!(m.vertices, Vertex(m.vtx)) + push!(m.normals, Normal(m.normal_sign * m.nrm)) + length(m.vertices) end - push!(m.vertices, Vertex(m.vtx)) - push!(m.normals, Normal(m.normal_sign * m.nrm)) - length(m.vertices) -end -include("example.jl") + include("example.jl") -@compile_workload begin - march(scenario(4, 4, 4; F = Float64, I = Int)) -end + @compile_workload begin + march(scenario(4, 4, 4; F = Float64, I = Int)) + end end diff --git a/src/example.jl b/src/example.jl index dbfd1a8..2206659 100644 --- a/src/example.jl +++ b/src/example.jl @@ -1,13 +1,13 @@ @inline cb_cushin(x, y, z) = ( z^2 * x^2 - z^4 - 2z * x^2 + 2z^3 + x^2 - z^2 - (x^2 - z)^2 - y^4 - 2x^2 * y^2 - - y^2 * z^2 + - 2y^2 * z + - y^2 + y^2 * z^2 + + 2y^2 * z + + y^2 ) @inline cb_sphere(x, y, z) = ( ((x - 2)^2 + (y - 2)^2 + (z - 2)^2 - 1) * - ((x + 2)^2 + (y - 2)^2 + (z - 2)^2 - 1) * - ((x - 2)^2 + (y + 2)^2 + (z - 2)^2 - 1) + ((x + 2)^2 + (y - 2)^2 + (z - 2)^2 - 1) * + ((x - 2)^2 + (y + 2)^2 + (z - 2)^2 - 1) ) @inline cb_plane(x, y, z) = x + y + z - 3 @inline cb_cassini(x::F, y::F, z::F) where {F} = @@ -20,17 +20,17 @@ (x^2 + y^2 + z^2 + b^2 - d^2)^2 - 4((a * x - c * d)^2 + b^2 * y^2) @inline cb_torus2(x::F, y::F, z::F, r = F(1.85), R = F(4)) where {F} = ( ((x^2 + y^2 + z^2 + R^2 - r^2)^2 - 4R^2 * (x^2 + y^2)) * - ((x^2 + (y + R)^2 + z^2 + R^2 - r^2)^2 - 4R^2 * ((y + R)^2 + z^2)) + ((x^2 + (y + R)^2 + z^2 + R^2 - r^2)^2 - 4R^2 * ((y + R)^2 + z^2)) ) @inline cb_mc_case(x::F, y::F, z::F) where {F} = ( F(-26.5298) * (1 - x) * (1 - y) * (1 - z) + - F(+81.9199) * x * (1 - y) * (1 - z) + - F(-100.680) * x * y * (1 - z) + - F(+3.54980) * (1 - x) * y * (1 - z) + - F(+24.1201) * (1 - x) * (1 - y) * z + - F(-74.4702) * x * (1 - y) * z + - F(+91.5298) * x * y * z + - F(-3.22998) * (1 - x) * y * z + F(+81.9199) * x * (1 - y) * (1 - z) + + F(-100.68) * x * y * (1 - z) + + F(+3.5498) * (1 - x) * y * (1 - z) + + F(+24.1201) * (1 - x) * (1 - y) * z + + F(-74.4702) * x * (1 - y) * z + + F(+91.5298) * x * y * z + + F(-3.22998) * (1 - x) * y * z ) @inline cb_drip(x::F, y::F, z::F, a = 2, b = 2, c = 3, d = 6) where {F} = x^2 + y^2 - F(0.5) * (F(0.995) * z^2 + F(0.005) - z^3) + F(0.0025) @@ -44,7 +44,7 @@ fill_volume!(vol::AbstractArray{F}, cb::Function) where {F} = begin ty = F(ny / 2sy + 1.5) tz = F(nz / 2sz) - @inbounds for k ∈ 1:nz, j ∈ 1:ny, i ∈ 1:nx + @inbounds for k in 1:nz, j in 1:ny, i in 1:nx vol[i, j, k] = cb(F((i - 1) / sx - tx), F((j - 1) / sy - ty), F((k - 1) / sz - tz)) end @@ -85,10 +85,10 @@ end output( PlyIO::Module, - m::MC{F,I}, + m::MC{F, I}, fn::AbstractString = "test.ply"; verbose = true, -) where {F,I} = begin +) where {F, I} = begin nv, nt = length(m.vertices), length(m.triangles) verbose && println("Writing $nv vertices and $nt triangles using `PlyIO`.") @@ -107,7 +107,7 @@ output( ) vertex_indices = PlyIO.ListProperty("vertex_indices", I, eltype(Triangle)) - for i ∈ eachindex(m.triangles) + for i in eachindex(m.triangles) push!(vertex_indices, m.triangles[i] .- 1) # 1-based indexing -> 0-based indexing end push!(ply, PlyIO.PlyElement("face", vertex_indices)) @@ -130,10 +130,10 @@ makemesh_GeometryBasics(GeometryBasics::Module, m::MC) = begin end makemesh(mod::Module, m::MC) = - if (mod_str = string(mod)) == "Meshes" - makemesh_Meshes(mod, m) - elseif mod_str == "GeometryBasics" - makemesh_GeometryBasics(mod, m) - else - throw(ArgumentError("un-supported module `$mod`")) - end +if (mod_str = string(mod)) == "Meshes" + makemesh_Meshes(mod, m) +elseif mod_str == "GeometryBasics" + makemesh_GeometryBasics(mod, m) +else + throw(ArgumentError("un-supported module `$mod`")) +end diff --git a/test/runtests.jl b/test/runtests.jl index 4c81f33..94d4ed3 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -129,7 +129,7 @@ end y = collect(Float64, range(start_y, stop_y, length = ny)), z = collect(Float64, range(start_z, stop_z, length = nz)), ) - for callable ∈ (march, march_legacy) + for callable in (march, march_legacy) callable(mc) @test all(map(v -> start_x ≤ v[1] ≤ stop_x, mc.vertices)) @test all(map(v -> start_y ≤ v[2] ≤ stop_y, mc.vertices)) @@ -138,7 +138,7 @@ end end @testset "isovalue" begin - dat = Float32[(x - 3)^2 + (y - 3)^2 + (z - 3)^2 for x ∈ 1:5, y ∈ 1:5, z ∈ 1:5] + dat = Float32[(x - 3)^2 + (y - 3)^2 + (z - 3)^2 for x in 1:5, y in 1:5, z in 1:5] m1 = MC(dat) march(m1, 5.0) @@ -152,7 +152,7 @@ end end @testset "invert normals" begin - dat = Float32[(x - 3)^2 + (y - 3)^2 + (z - 3)^2 for x ∈ 1:5, y ∈ 1:5, z ∈ 1:5] + dat = Float32[(x - 3)^2 + (y - 3)^2 + (z - 3)^2 for x in 1:5, y in 1:5, z in 1:5] m1 = MC(dat) march(m1) @@ -184,7 +184,7 @@ end end @testset "coordinate input variations" begin - atol = 1e-3 # precision level + atol = 1.0e-3 # precision level # define coordinate ranges (also creating 3 different lengths) nx, ny, nz = 55, 46, 67 # these should be high enough to reach precision level @@ -196,7 +196,7 @@ end z = range(start_z, stop_z; length = nz) # create image (simple coordinate norm leading to spherical isosurface) - A = [√(xi^2 + yi^2 + zi^2) for xi ∈ x, yi ∈ y, zi ∈ z] + A = [√(xi^2 + yi^2 + zi^2) for xi in x, yi in y, zi in z] level = 0.5 # isolevel should produce sphere with this radius @@ -214,17 +214,17 @@ end @test mc_ranged.vertices == mc_vector.vertices @test mc_ranged.triangles == mc_vector.triangles - # test if coordinate input was used appropriately geometrically as expected + # test if coordinate input was used appropriately geometrically as expected n = length(mc_ranged.vertices) - c = sum(mc_ranged.vertices) / n # mean coordinate i.e. center - r = sum(v -> √(sum(abs2, v)), mc_ranged.vertices) / n # mean radius - @test isapprox(c, zeros(3); atol) # approximately zero mean for sphere + c = sum(mc_ranged.vertices) / n # mean coordinate i.e. center + r = sum(v -> √(sum(abs2, v)), mc_ranged.vertices) / n # mean radius + @test isapprox(c, zeros(3); atol) # approximately zero mean for sphere @test isapprox(r, level; atol) # approximately radius matching level end @testset "types" begin - for F ∈ (Float16, Float32, Float64), - I ∈ (Int16, Int32, Int64, Int128, UInt16, UInt32, UInt64, UInt128) + for F in (Float16, Float32, Float64), + I in (Int16, Int32, Int64, Int128, UInt16, UInt32, UInt64, UInt128) @test march(MarchingCubes.scenario(4, 4, 4; F, I)) isa Nothing end From 7a5494b20aeffff45948f7267a343d9500da5d5a Mon Sep 17 00:00:00 2001 From: t-bltg Date: Mon, 26 Jan 2026 12:46:23 +0100 Subject: [PATCH 3/3] up --- .github/workflows/ci.yml | 18 +++++------------- 1 file changed, 5 insertions(+), 13 deletions(-) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index b91803e..6f6c0b1 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -5,7 +5,7 @@ on: push: branches: [main] -concurrency: +concurrency: group: ${{ github.workflow }}-${{ github.head_ref || github.run_id }} cancel-in-progress: true @@ -22,24 +22,16 @@ jobs: - '1' experimental: - false - os: [ubuntu-latest] + os: [ubuntu-latest, windows-latest, macos-latest] arch: [x64] - include: # spare windows/macos CI credits - - os: windows-latest - experimental: false - version: '1' - arch: x64 - - os: macos-latest - experimental: false - version: '1' - arch: x64 + include: - os: ubuntu-latest experimental: true version: 'pre' arch: x64 steps: - - uses: actions/checkout@v5 + - uses: actions/checkout@v6 - uses: julia-actions/setup-julia@latest with: version: ${{ matrix.version }} @@ -48,7 +40,7 @@ jobs: - uses: julia-actions/julia-buildpkg@latest - uses: julia-actions/julia-runtest@latest - uses: julia-actions/julia-processcoverage@latest - - uses: codecov/codecov-action@v4 + - uses: codecov/codecov-action@v5 with: token: ${{ secrets.CODECOV_TOKEN }} fail_ci_if_error: false