diff --git a/julia/MOLE.jl/docs/src/index.md b/julia/MOLE.jl/docs/src/index.md index 9b8a2be7..cd46bb2f 100644 --- a/julia/MOLE.jl/docs/src/index.md +++ b/julia/MOLE.jl/docs/src/index.md @@ -86,6 +86,31 @@ Once you have built the documentation (either from the REPL or the command line) ## Functions +### Grids + +The `Grids` module provides the grid data structures. A grid stores dimensionality, topology, spacing, coordinates, and boundary metadata in one object. + +Example of usage: + +```julia +using MOLE.Grids + +grid = makeGrid(m=64, n=64, dx=1/64, dy=1/64) + +grid.dim +grid.topology +grid.nodes +grid.faces +grid.centers +``` + +```@docs +MOLE.Grids.makeGrid +MOLE.Grids.validateGrid +MOLE.Grids.Grid +MOLE.Grids.BoundaryMetadata +``` + ### Operators ```@docs diff --git a/julia/MOLE.jl/src/Grids/Grids.jl b/julia/MOLE.jl/src/Grids/Grids.jl new file mode 100644 index 00000000..735c112d --- /dev/null +++ b/julia/MOLE.jl/src/Grids/Grids.jl @@ -0,0 +1,11 @@ +module Grids + +export AbstractGrid, Grid, BoundaryMetadata +export makeGrid, validateGrid + +include("topologies.jl") +include("coordinates.jl") +include("validate_grid.jl") +include("make_grid.jl") + +end # module Grids diff --git a/julia/MOLE.jl/src/Grids/coordinates.jl b/julia/MOLE.jl/src/Grids/coordinates.jl new file mode 100644 index 00000000..545aad95 --- /dev/null +++ b/julia/MOLE.jl/src/Grids/coordinates.jl @@ -0,0 +1,120 @@ +#= + SPDX-License-Identifier: GPL-3.0-or-later + © 2008-2024 San Diego State University Research Foundation (SDSURF). + See LICENSE file or https://www.gnu.org/licenses/gpl-3.0.html for details. +=# + +function _coordinates_1d(m::Int, dx) + x_node = collect(0:m) .* dx + x_center = vcat(0.0, collect(0.5:1:(m - 0.5)) .* dx, m * dx) + x_face = x_node + + nodes = (; X = x_node) + faces = (; X = x_face) + centers = (; X = x_center) + + return nodes, faces, centers +end + +function _coordinates_2d(m::Int, n::Int, dx, dy) + xn = collect(0:m) .* dx + yn = collect(0:n) .* dy + + xc = vcat(0.0, collect(0.5:1:(m - 0.5)) .* dx, m * dx) + yc = vcat(0.0, collect(0.5:1:(n - 0.5)) .* dy, n * dy) + + xu = xn + yu = collect(0.5:1:(n - 0.5)) .* dy + + xv = collect(0.5:1:(m - 0.5)) .* dx + yv = yn + + nodes = _ndgrid2_named(xn, yn) + centers = _ndgrid2_named(xc, yc) + + faces = ( + u = _ndgrid2_named(xu, yu), + v = _ndgrid2_named(xv, yv), + ) + + return nodes, faces, centers +end + +function _coordinates_3d(m::Int, n::Int, o::Int, dx, dy, dz) + xn = collect(0:m) .* dx + yn = collect(0:n) .* dy + zn = collect(0:o) .* dz + + xc = vcat(0.0, collect(0.5:1:(m - 0.5)) .* dx, m * dx) + yc = vcat(0.0, collect(0.5:1:(n - 0.5)) .* dy, n * dy) + zc = vcat(0.0, collect(0.5:1:(o - 0.5)) .* dz, o * dz) + + nodes = _ndgrid3_named(xn, yn, zn) + centers = _ndgrid3_named(xc, yc, zc) + + faces = ( + u = _ndgrid3_named( + xn, + collect(0.5:1:(n - 0.5)) .* dy, + collect(0.5:1:(o - 0.5)) .* dz, + ), + v = _ndgrid3_named( + collect(0.5:1:(m - 0.5)) .* dx, + yn, + collect(0.5:1:(o - 0.5)) .* dz, + ), + w = _ndgrid3_named( + collect(0.5:1:(m - 0.5)) .* dx, + collect(0.5:1:(n - 0.5)) .* dy, + zn, + ), + ) + + return nodes, faces, centers +end + +function _coordinates_curvilinear_2d(nodes) + NX = nodes.X + NY = nodes.Y + + faces = ( + u = ( + X = 0.5 .* (NX[:, 1:(end - 1)] .+ NX[:, 2:end]), + Y = 0.5 .* (NY[:, 1:(end - 1)] .+ NY[:, 2:end]), + ), + v = ( + X = 0.5 .* (NX[1:(end - 1), :] .+ NX[2:end, :]), + Y = 0.5 .* (NY[1:(end - 1), :] .+ NY[2:end, :]), + ), + ) + + centers = ( + X = 0.25 .* ( + NX[1:(end - 1), 1:(end - 1)] .+ + NX[2:end, 1:(end - 1)] .+ + NX[1:(end - 1), 2:end] .+ + NX[2:end, 2:end] + ), + Y = 0.25 .* ( + NY[1:(end - 1), 1:(end - 1)] .+ + NY[2:end, 1:(end - 1)] .+ + NY[1:(end - 1), 2:end] .+ + NY[2:end, 2:end] + ), + ) + + return faces, centers +end + +function _ndgrid2_named(x, y) + X = repeat(reshape(x, :, 1), 1, length(y)) + Y = repeat(reshape(y, 1, :), length(x), 1) + return (; X, Y) +end + +function _ndgrid3_named(x, y, z) + X = repeat(reshape(x, :, 1, 1), 1, length(y), length(z)) + Y = repeat(reshape(y, 1, :, 1), length(x), 1, length(z)) + Z = repeat(reshape(z, 1, 1, :), length(x), length(y), 1) + return (; X, Y, Z) +end diff --git a/julia/MOLE.jl/src/Grids/make_grid.jl b/julia/MOLE.jl/src/Grids/make_grid.jl new file mode 100644 index 00000000..eaa19e24 --- /dev/null +++ b/julia/MOLE.jl/src/Grids/make_grid.jl @@ -0,0 +1,83 @@ +#= + SPDX-License-Identifier: GPL-3.0-or-later + © 2008-2024 San Diego State University Research Foundation (SDSURF). + See LICENSE file or https://www.gnu.org/licenses/gpl-3.0.html for details. +=# + +""" + makeGrid(grid::Grid; allowPartial::Bool = true) + makeGrid(; allowPartial::Bool = true, kwargs...) + +Construct and validate a computational grid. + +`makeGrid` is the primary constructor for `Grid` objects. It accepts either an +existing `Grid` instance or keyword arguments describing the computational +domain. The resulting grid is validated, normalized, and enriched with derived +metadata such as the spatial dimension, topology, coordinate arrays, and +boundary-condition information. + +# Arguments + +- `grid::Grid`: An existing grid to validate and normalize. +- `allowPartial::Bool=true`: If `true`, allows partially specified grids to be + created without requiring all fields. This is primarily intended for + incremental grid construction or internal use. + +When called with keyword arguments, commonly used parameters include: + +- `m`, `n`, `o`: Number of cells in each spatial direction. +- `dx`, `dy`, `dz`: Grid spacing in each spatial direction. +- `topology`: Grid topology (e.g. `:uniform`, `:periodic`, `:curvilinear`). +- `nodes`: User-supplied nodal coordinates for curvilinear grids. +- `bc`: Boundary-condition metadata. + +# Returns + +A validated [`Grid`](@ref) object. + +# Examples + +Construct a one-dimensional uniform grid: + +```julia +grid = makeGrid( + m = 100, + dx = 0.01, +) +``` + +Construct a two-dimensional uniform grid: + +```julia +grid = makeGrid( + m = 64, + n = 64, + dx = 1/64, + dy = 1/64, +) +``` + +Construct a curvilinear grid from user-defined nodal coordinates: + +```julia +grid = makeGrid( + m = m, + n = n, + topology = :curvilinear, + nodes = (X = X, Y = Y), +) +``` + +# Notes + +`makeGrid` delegates the validation and normalization of all inputs to +[`validateGrid`](@ref). +""" +function makeGrid(grid::Grid; allowPartial::Bool = true) + return validateGrid(grid; allowPartial = allowPartial) +end + +function makeGrid(; allowPartial::Bool = true, kwargs...) + raw = Dict{Symbol, Any}(kwargs) + return validateGrid(raw; allowPartial = allowPartial) +end diff --git a/julia/MOLE.jl/src/Grids/topologies.jl b/julia/MOLE.jl/src/Grids/topologies.jl new file mode 100644 index 00000000..a58bf30e --- /dev/null +++ b/julia/MOLE.jl/src/Grids/topologies.jl @@ -0,0 +1,114 @@ +#= + SPDX-License-Identifier: GPL-3.0-or-later + © 2008-2024 San Diego State University Research Foundation (SDSURF). + See LICENSE file or https://www.gnu.org/licenses/gpl-3.0.html for details. +=# + +abstract type AbstractGrid end + +""" + BoundaryMetadata{T} + +Stores the boundary-condition metadata associated with a computational grid. + +`BoundaryMetadata` is used by [`Grid`](@ref) and contains the coefficients +required to define scalar boundary conditions, together with flags indicating +whether each coordinate direction is periodic. + +# Fields + +- `dc::Vector{T}`: Dirichlet coefficients for each boundary. +- `nc::Vector{T}`: Neumann coefficients for each boundary. +- `isPeriodic::Vector{Bool}`: Periodicity flags for each coordinate direction. +- `hasData::Bool`: Indicates whether boundary-condition information has been + explicitly provided. + +The expected lengths of `dc` and `nc` depend on the grid dimension: + +| Dimension | Length | +|:----------|:------:| +| 1D | 2 | +| 2D | 4 | +| 3D | 6 | + +The length of `isPeriodic` is equal to the number of spatial dimensions. + +# Examples + +```julia +bc = BoundaryMetadata( + dc = [1.0, 1.0], + nc = [0.0, 0.0], + isPeriodic = [false], + hasData = true, +) +``` + +A periodic one-dimensional grid can be represented by zero Dirichlet and +Neumann coefficients: + +```julia +bc = BoundaryMetadata( + dc = [0.0, 0.0], + nc = [0.0, 0.0], + isPeriodic = [true], + hasData = true, +) +``` +""" +Base.@kwdef struct BoundaryMetadata{T} + dc::Vector{T} = T[] + nc::Vector{T} = T[] + isPeriodic::Vector{Bool} = Bool[] + hasData::Bool = false +end + +""" + Grid{T} <: AbstractGrid + +Stores the geometric and boundary-condition metadata for a computational grid. + +A `Grid` is the central data object used by the MOLE 2.0 grid API. It records +the grid dimension, topology, cell counts, grid spacings, coordinate arrays, +and boundary-condition metadata. Operators should treat `Grid` objects as +read-only inputs. + +# Fields + +- `dim::Int`: Spatial dimension of the grid. +- `topology::Symbol`: Grid topology, such as `:uniform`, `:periodic`, or `:curvilinear`. +- `m`, `n`, `o`: Number of cells in the x-, y-, and z-directions. +- `dx`, `dy`, `dz`: Grid spacing in the x-, y-, and z-directions. +- `nodes::NamedTuple`: Nodal coordinate arrays. +- `faces::NamedTuple`: Face coordinate arrays. +- `centers::NamedTuple`: Cell-center coordinate arrays. +- `bc::BoundaryMetadata{T}`: Boundary-condition metadata. + +# Examples + +```julia +grid = makeGrid(m = 64, dx = 1/64) + +grid.dim +grid.topology +grid.nodes +``` +""" +Base.@kwdef struct Grid{T} <: AbstractGrid + dim::Int + topology::Symbol + + m::Union{Nothing, Int} = nothing + n::Union{Nothing, Int} = nothing + o::Union{Nothing, Int} = nothing + + dx::Union{Nothing, T} = nothing + dy::Union{Nothing, T} = nothing + dz::Union{Nothing, T} = nothing + + nodes::NamedTuple = NamedTuple() + faces::NamedTuple = NamedTuple() + centers::NamedTuple = NamedTuple() + + bc::BoundaryMetadata{T} = BoundaryMetadata{T}() +end diff --git a/julia/MOLE.jl/src/Grids/validate_grid.jl b/julia/MOLE.jl/src/Grids/validate_grid.jl new file mode 100644 index 00000000..25a7c3c3 --- /dev/null +++ b/julia/MOLE.jl/src/Grids/validate_grid.jl @@ -0,0 +1,450 @@ +#= + SPDX-License-Identifier: GPL-3.0-or-later + © 2008-2024 San Diego State University Research Foundation (SDSURF). + See LICENSE file or https://www.gnu.org/licenses/gpl-3.0.html for details. +=# + +""" + validateGrid(grid::Grid; allowPartial=false) + validateGrid(raw::AbstractDict{Symbol,<:Any}; allowPartial=false) + +Validate and normalize a computational grid. + +`validateGrid` verifies that a grid definition is internally consistent, +infers missing metadata when possible, and returns a fully normalized +[`Grid`](@ref) object. Depending on the supplied information, it determines the +spatial dimension and grid topology, validates boundary-condition data, and +generates the coordinate arrays associated with nodes, faces, and cell centers. + +This function is called internally by [`makeGrid`](@ref), but may also be +used directly to validate existing grid objects or dictionaries containing +grid definitions. + +# Arguments + +- `grid::Grid`: A grid object to validate and normalize. +- `raw::AbstractDict{Symbol,<:Any}`: A dictionary containing grid parameters. +- `allowPartial::Bool=false`: If `true`, incomplete grid definitions are + accepted without requiring all fields needed for coordinate generation. + +# Validation + +Depending on the inferred spatial dimension, `validateGrid` checks the +consistency of: + +- grid dimensions (`m`, `n`, `o`) +- grid spacing (`dx`, `dy`, `dz`) +- grid topology (`:uniform`, `:periodic`, `:curvilinear`) +- boundary-condition metadata +- user-supplied nodal coordinates for curvilinear grids + +When sufficient information is available, coordinate arrays for nodes, faces, +and cell centers are generated automatically. + +# Returns + +A validated and normalized [`Grid`](@ref) object. + +# Throws + +- `ArgumentError` if required fields are missing. +- `ArgumentError` if grid dimensions or spacings are invalid. +- `ArgumentError` if boundary-condition metadata is inconsistent. +- `ArgumentError` if supplied curvilinear coordinate arrays have incompatible + dimensions. + +# Examples + +Validate a one-dimensional grid definition: + +```julia +grid = validateGrid(Dict( + :m => 100, + :dx => 0.01, +)) +``` + +Validate an existing grid object: + +```julia +grid = makeGrid(m=64, dx=1/64) + +grid = validateGrid(grid) +``` + +Validate a partially specified grid: + +```julia +grid = validateGrid( + Dict(:m => 100); + allowPartial = true, +) +``` + +# Notes + +Most users should construct grids using [`makeGrid`](@ref), which calls +`validateGrid` automatically. Direct use of `validateGrid` is primarily intended +for internal routines, advanced workflows, or validation of externally +constructed grid definitions. +""" +function validateGrid(grid::Grid; allowPartial::Bool = false) + raw = Dict{Symbol, Any}() + + raw[:dim] = grid.dim + raw[:topology] = grid.topology + + for name in (:m, :n, :o, :dx, :dy, :dz) + value = getfield(grid, name) + value !== nothing && (raw[name] = value) + end + + !isempty(grid.nodes) && (raw[:nodes] = grid.nodes) + !isempty(grid.faces) && (raw[:faces] = grid.faces) + !isempty(grid.centers) && (raw[:centers] = grid.centers) + raw[:bc] = grid.bc + + return validateGrid(raw; allowPartial = allowPartial) +end + +function validateGrid(raw::AbstractDict{Symbol, <:Any}; allowPartial::Bool = false) + bc_raw = get(raw, :bc, Dict{Symbol, Any}()) + + if haskey(raw, :dc) && !_has_bc_field(bc_raw, :dc) + bc_raw = _set_bc_field(bc_raw, :dc, raw[:dc]) + end + + if haskey(raw, :nc) && !_has_bc_field(bc_raw, :nc) + bc_raw = _set_bc_field(bc_raw, :nc, raw[:nc]) + end + + dim = get(raw, :dim, _infer_dim(raw)) + topology = Symbol(get(raw, :topology, _infer_topology(raw))) + + dim == 1 && return _normalize_1d(raw, bc_raw, topology; allowPartial = allowPartial) + dim == 2 && return _normalize_2d(raw, bc_raw, topology; allowPartial = allowPartial) + dim == 3 && return _normalize_3d(raw, bc_raw, topology; allowPartial = allowPartial) + + throw(ArgumentError("grid.dim must be 1, 2, or 3")) +end + +function _infer_dim(raw) + if any(haskey(raw, k) for k in (:o, :dz, :Z)) + return 3 + elseif any(haskey(raw, k) for k in (:n, :dy, :Y)) + return 2 + else + return 1 + end +end + +function _infer_topology(raw) + if any(haskey(raw, k) for k in (:X, :Y, :Z)) || + (haskey(raw, :nodes) && all(k -> haskey(raw[:nodes], k), (:X, :Y))) + return :curvilinear + elseif any(haskey(raw, k) for k in (:x, :y, :z)) + return :nonuniform + else + return :uniform + end +end + +function _require(raw, names, allowPartial, msg) + for name in names + if !haskey(raw, name) && !allowPartial + throw(ArgumentError("$msg: missing grid.$name")) + end + end +end + +function _normalize_1d(raw, bc_raw, topology; allowPartial) + _require(raw, (:m,), allowPartial, "validateGrid:MissingField1D") + + if haskey(raw, :m) && haskey(raw, :dx) + m = _validate_positive_int(raw[:m], "grid.m") + dx = _validate_positive_spacing(raw[:dx], "grid.dx") + + bc = _normalize_bc(bc_raw, 2) + isperiodic = bc.hasData ? [all(bc.dc .^ 2 .+ bc.nc .^ 2 .== 0)] : [false] + bc = BoundaryMetadata( + dc = bc.dc, + nc = bc.nc, + isPeriodic = isperiodic, + hasData = bc.hasData, + ) + + nodes, faces, centers = _coordinates_1d(m, dx) + grid_topology = only(bc.isPeriodic) ? :periodic : :uniform + + return Grid( + dim = 1, + topology = grid_topology, + m = m, + dx = dx, + nodes = nodes, + faces = faces, + centers = centers, + bc = bc, + ) + end + + if !allowPartial && topology == :uniform + throw( + ArgumentError( + "validateGrid:MissingUniform1D: uniform 1-D grid requires m and dx", + ), + ) + end + + return Grid{Float64}(dim = 1, topology = topology) +end + +function _normalize_2d(raw, bc_raw, topology; allowPartial) + _require(raw, (:m, :n), allowPartial, "validateGrid:MissingField2D") + + if topology == :curvilinear + if haskey(raw, :m) && haskey(raw, :n) + return _normalize_curvilinear_2d(raw, bc_raw) + end + return Grid(dim = 2, topology = :curvilinear) + end + + if all(haskey(raw, k) for k in (:m, :n, :dx, :dy)) + m = _validate_positive_int(raw[:m], "grid.m") + n = _validate_positive_int(raw[:n], "grid.n") + dx = _validate_positive_spacing(raw[:dx], "grid.dx") + dy = _validate_positive_spacing(raw[:dy], "grid.dy") + + bc = _normalize_bc(bc_raw, 4) + isperiodic = + bc.hasData ? + [ + all(bc.dc[1:2] .^ 2 .+ bc.nc[1:2] .^ 2 .== 0), + all(bc.dc[3:4] .^ 2 .+ bc.nc[3:4] .^ 2 .== 0), + ] : [false, false] + + bc = BoundaryMetadata( + dc = bc.dc, + nc = bc.nc, + isPeriodic = isperiodic, + hasData = bc.hasData, + ) + nodes, faces, centers = _coordinates_2d(m, n, dx, dy) + grid_topology = any(bc.isPeriodic) ? :periodic : :uniform + + return Grid( + dim = 2, + topology = grid_topology, + m = m, + n = n, + dx = dx, + dy = dy, + nodes = nodes, + faces = faces, + centers = centers, + bc = bc, + ) + end + + if !allowPartial && topology == :uniform + throw( + ArgumentError( + "validateGrid:MissingUniform2D: uniform 2-D grid requires m, n, dx, and dy", + ), + ) + end + + return Grid{Float64}(dim = 2, topology = topology) +end + +function _normalize_3d(raw, bc_raw, topology; allowPartial) + _require(raw, (:m, :n, :o), allowPartial, "validateGrid:MissingField3D") + + if all(haskey(raw, k) for k in (:m, :n, :o, :dx, :dy, :dz)) + m = _validate_positive_int(raw[:m], "grid.m") + n = _validate_positive_int(raw[:n], "grid.n") + o = _validate_positive_int(raw[:o], "grid.o") + + dx = _validate_positive_spacing(raw[:dx], "grid.dx") + dy = _validate_positive_spacing(raw[:dy], "grid.dy") + dz = _validate_positive_spacing(raw[:dz], "grid.dz") + + bc = _normalize_bc(bc_raw, 6) + isperiodic = + bc.hasData ? + [ + all(bc.dc[1:2] .^ 2 .+ bc.nc[1:2] .^ 2 .== 0), + all(bc.dc[3:4] .^ 2 .+ bc.nc[3:4] .^ 2 .== 0), + all(bc.dc[5:6] .^ 2 .+ bc.nc[5:6] .^ 2 .== 0), + ] : [false, false, false] + + bc = BoundaryMetadata( + dc = bc.dc, + nc = bc.nc, + isPeriodic = isperiodic, + hasData = bc.hasData, + ) + nodes, faces, centers = _coordinates_3d(m, n, o, dx, dy, dz) + grid_topology = any(bc.isPeriodic) ? :periodic : :uniform + + return Grid( + dim = 3, + topology = grid_topology, + m = m, + n = n, + o = o, + dx = dx, + dy = dy, + dz = dz, + nodes = nodes, + faces = faces, + centers = centers, + bc = bc, + ) + end + + if !allowPartial && topology == :uniform + throw( + ArgumentError( + "validateGrid:MissingUniform3D: uniform 3-D grid requires m, n, o, dx, dy, and dz", + ), + ) + end + + return Grid{Float64}(dim = 3, topology = topology) +end + +function _normalize_curvilinear_2d(raw, bc_raw) + m = _validate_positive_int(raw[:m], "grid.m") + n = _validate_positive_int(raw[:n], "grid.n") + + if !haskey(raw, :nodes) + throw( + ArgumentError( + "validateGrid:CurvilinearMissingNodes: curvilinear grid requires nodes.X and nodes.Y", + ), + ) + end + + nodes = raw[:nodes] + + if !(haskey(nodes, :X) && haskey(nodes, :Y)) + throw( + ArgumentError( + "validateGrid:CurvilinearMissingNodes: curvilinear grid requires nodes.X and nodes.Y", + ), + ) + end + + expected = (m + 1, n + 1) + + if size(nodes.X) != expected || size(nodes.Y) != expected + throw( + ArgumentError( + "validateGrid:SizeMismatch: curvilinear nodes.X/Y must have size $expected", + ), + ) + end + + faces, centers = _coordinates_curvilinear_2d(nodes) + bc = _normalize_bc(bc_raw, 4) + + return Grid( + dim = 2, + topology = :curvilinear, + m = m, + n = n, + nodes = nodes, + faces = faces, + centers = centers, + bc = bc, + ) +end + +function _normalize_bc(bc_raw, expected) + has_dc = _has_bc_field(bc_raw, :dc) + has_nc = _has_bc_field(bc_raw, :nc) + + if !has_dc && !has_nc + return BoundaryMetadata{Float64}() + end + + if has_dc != has_nc + throw(ArgumentError("grid.bc.dc and grid.bc.nc must be provided together")) + end + + dc = _bc_field(bc_raw, :dc) + nc = _bc_field(bc_raw, :nc) + + if isempty(dc) && isempty(nc) + return BoundaryMetadata{Float64}() + end + + dc = _normalize_bc_vector(dc, expected, "grid.bc.dc") + nc = _normalize_bc_vector(nc, expected, "grid.bc.nc") + + return BoundaryMetadata(dc = dc, nc = nc, isPeriodic = Bool[], hasData = true) +end + +function _normalize_bc_vector(values, expected, name) + values isa Number && return fill(float(values), expected) + + if !(values isa AbstractVector) + throw(ArgumentError("$name must be a scalar or vector")) + end + + vals = float.(collect(values)) + + length(vals) == 1 && return fill(only(vals), expected) + + if length(vals) != expected + throw(ArgumentError("$name must be a scalar or a vector of length $expected")) + end + + return vals +end + +_has_bc_field(bc::BoundaryMetadata, name::Symbol) = !isempty(getfield(bc, name)) +_has_bc_field(bc::AbstractDict, name::Symbol) = haskey(bc, name) +_has_bc_field(bc::NamedTuple, name::Symbol) = haskey(bc, name) +_has_bc_field(_, _) = false + +_bc_field(bc::BoundaryMetadata, name::Symbol) = getfield(bc, name) +_bc_field(bc::AbstractDict, name::Symbol) = bc[name] +_bc_field(bc::NamedTuple, name::Symbol) = getfield(bc, name) + +function _set_bc_field(bc::AbstractDict, name::Symbol, value) + out = Dict{Symbol, Any}(bc) + out[name] = value + return out +end + +function _set_bc_field(bc::NamedTuple, name::Symbol, value) + return merge(bc, (; name => value)) +end + +function _set_bc_field(bc::BoundaryMetadata, name::Symbol, value) + values = name == :dc ? value : bc.dc + normals = name == :nc ? value : bc.nc + return BoundaryMetadata(dc = collect(float.(values)), nc = collect(float.(normals))) +end + +function _validate_positive_int(value, name) + value = Int(value) + + if value <= 0 + throw(ArgumentError("$name must be a positive integer")) + end + + return value +end + +function _validate_positive_spacing(value, name) + value = float(value) + + if value <= 0 + throw(ArgumentError("$name must be positive")) + end + + return value +end diff --git a/julia/MOLE.jl/src/MOLE.jl b/julia/MOLE.jl/src/MOLE.jl index 77ad9d07..340ff0dd 100644 --- a/julia/MOLE.jl/src/MOLE.jl +++ b/julia/MOLE.jl/src/MOLE.jl @@ -8,5 +8,6 @@ module MOLE include("Operators/Operators.jl") include("BCs/BCs.jl") +include("Grids/Grids.jl") end # module MOLE diff --git a/julia/MOLE.jl/src/Operators/interpol.jl b/julia/MOLE.jl/src/Operators/interpol.jl index 1b96a5c9..a2864724 100644 --- a/julia/MOLE.jl/src/Operators/interpol.jl +++ b/julia/MOLE.jl/src/Operators/interpol.jl @@ -1,10 +1,10 @@ -using SparseArrays #= SPDX-License-Identifier: GPL-3.0-or-later © 2008-2024 San Diego State University Research Foundation (SDSURF). See LICENSE file or https://www.gnu.org/licenses/gpl-3.0.html for details. =# +using SparseArrays """ interpol(m, c) diff --git a/julia/MOLE.jl/test/Grids/make_grid.jl b/julia/MOLE.jl/test/Grids/make_grid.jl new file mode 100644 index 00000000..954a9fa3 --- /dev/null +++ b/julia/MOLE.jl/test/Grids/make_grid.jl @@ -0,0 +1,22 @@ +using MOLE.Grids + +@testset "makeGrid" begin + grid = makeGrid(m = 4, dx = 0.25) + + @test grid isa Grid + @test grid.dim == 1 + @test grid.topology == :uniform + @test grid.m == 4 + @test grid.dx == 0.25 + @test grid.nodes.X == [0.0, 0.25, 0.5, 0.75, 1.0] + + grid2 = makeGrid(m = 2, n = 3, dx = 0.5, dy = 0.25) + + @test grid2.dim == 2 + @test grid2.topology == :uniform + @test size(grid2.nodes.X) == (3, 4) + @test size(grid2.nodes.Y) == (3, 4) + @test size(grid2.centers.X) == (4, 5) + @test size(grid2.faces.u.X) == (3, 3) + @test size(grid2.faces.v.X) == (2, 4) +end diff --git a/julia/MOLE.jl/test/Grids/validate_grid.jl b/julia/MOLE.jl/test/Grids/validate_grid.jl new file mode 100644 index 00000000..26d3af10 --- /dev/null +++ b/julia/MOLE.jl/test/Grids/validate_grid.jl @@ -0,0 +1,247 @@ +using MOLE.Grids + +@testset "validateGrid uniform 1D" begin + grid = validateGrid(Dict(:m => 3, :dx => 0.1)) + + @test grid.dim == 1 + @test grid.topology == :uniform + @test grid.bc.hasData == false + @test grid.bc.isPeriodic == [false] + @test length(grid.nodes.X) == 4 + @test length(grid.centers.X) == 5 + @test length(grid.faces.X) == 4 +end + +@testset "validateGrid uniform 2D" begin + grid = validateGrid(Dict(:m => 2, :n => 4, :dx => 0.5, :dy => 0.25)) + + @test grid.dim == 2 + @test grid.topology == :uniform + + @test size(grid.nodes.X) == (3, 5) + @test size(grid.nodes.Y) == (3, 5) + + @test size(grid.centers.X) == (4, 6) + @test size(grid.centers.Y) == (4, 6) + + @test size(grid.faces.u.X) == (3, 4) + @test size(grid.faces.u.Y) == (3, 4) + + @test size(grid.faces.v.X) == (2, 5) + @test size(grid.faces.v.Y) == (2, 5) +end + +@testset "validateGrid uniform 3D" begin + grid = + validateGrid(Dict(:m => 2, :n => 3, :o => 4, :dx => 0.5, :dy => 0.25, :dz => 0.1)) + + @test grid.dim == 3 + @test grid.topology == :uniform + + @test size(grid.nodes.X) == (3, 4, 5) + @test size(grid.centers.X) == (4, 5, 6) + + @test size(grid.faces.u.X) == (3, 3, 4) + @test size(grid.faces.v.X) == (2, 4, 4) + @test size(grid.faces.w.X) == (2, 3, 5) +end + +@testset "boundary normalization" begin + grid = makeGrid(m = 4, dx = 0.25, dc = 0.0, nc = 0.0) + + @test grid.topology == :periodic + @test grid.bc.hasData + @test grid.bc.dc == [0.0, 0.0] + @test grid.bc.nc == [0.0, 0.0] + @test grid.bc.isPeriodic == [true] + + grid2 = makeGrid( + m = 4, + n = 5, + dx = 0.25, + dy = 0.2, + bc = (dc = [0.0, 0.0, 1.0, 1.0], nc = [0.0, 0.0, 1.0, 1.0]), + ) + + @test grid2.topology == :periodic + @test grid2.bc.isPeriodic == [true, false] +end + +@testset "curvilinear 2D" begin + m = 2 + n = 3 + + x = collect(0:m) + y = collect(0:n) + + X = repeat(reshape(x, :, 1), 1, n + 1) + Y = repeat(reshape(y, 1, :), m + 1, 1) + + grid = makeGrid(m = m, n = n, topology = :curvilinear, nodes = (; X, Y)) + + @test grid.dim == 2 + @test grid.topology == :curvilinear + + @test size(grid.nodes.X) == (3, 4) + @test size(grid.faces.u.X) == (3, 3) + @test size(grid.faces.v.X) == (2, 4) + @test size(grid.centers.X) == (2, 3) + + @test grid.centers.X[1, 1] == 0.5 + @test grid.centers.Y[1, 1] == 0.5 +end + +@testset "validation errors" begin + @test_throws ArgumentError validateGrid(Dict(:dx => 0.1); allowPartial = false) + @test_throws ArgumentError validateGrid( + Dict(:m => 4, :n => 4, :topology => :curvilinear); + allowPartial = false, + ) + + bad_nodes = (; X = zeros(2, 2), Y = zeros(2, 2)) + @test_throws ArgumentError makeGrid( + m = 4, + n = 4, + topology = :curvilinear, + nodes = bad_nodes, + ) +end + +@testset "validation tests" begin + @testset "makeGrid from existing Grid" begin + grid = makeGrid(m = 4, dx = 0.25) + grid2 = makeGrid(grid) + + @test grid2.dim == 1 + @test grid2.topology == :uniform + @test grid2.m == 4 + @test grid2.dx == 0.25 + end + + @testset "partial grids" begin + grid1 = makeGrid(m = 10; allowPartial = true) + @test grid1.dim == 1 + @test grid1.topology == :uniform + + grid2 = makeGrid(m = 10, n = 20; allowPartial = true) + @test grid2.dim == 2 + @test grid2.topology == :uniform + + grid3 = makeGrid(m = 10, n = 20, o = 30; allowPartial = true) + @test grid3.dim == 3 + @test grid3.topology == :uniform + end + + @testset "missing uniform spacing errors" begin + @test_throws ArgumentError makeGrid(m = 4; allowPartial = false) + @test_throws ArgumentError makeGrid(m = 4, n = 4, dx = 0.1; allowPartial = false) + @test_throws ArgumentError makeGrid( + m = 4, + n = 4, + o = 4, + dx = 0.1, + dy = 0.1; + allowPartial = false, + ) + end + + @testset "invalid dimensions and spacings" begin + @test_throws ArgumentError makeGrid(m = 0, dx = 0.1) + @test_throws ArgumentError makeGrid(m = -4, dx = 0.1) + @test_throws ArgumentError makeGrid(m = 4, dx = 0.0) + @test_throws ArgumentError makeGrid(m = 4, dx = -0.1) + + @test_throws ArgumentError makeGrid(m = 4, n = 0, dx = 0.1, dy = 0.1) + @test_throws ArgumentError makeGrid(m = 4, n = 4, dx = 0.1, dy = 0.0) + + @test_throws ArgumentError makeGrid( + m = 4, + n = 4, + o = 0, + dx = 0.1, + dy = 0.1, + dz = 0.1, + ) + @test_throws ArgumentError makeGrid( + m = 4, + n = 4, + o = 4, + dx = 0.1, + dy = 0.1, + dz = 0.0, + ) + end + + @testset "invalid dim" begin + @test_throws ArgumentError validateGrid(Dict(:dim => 4); allowPartial = true) + end + + @testset "topology inference" begin + grid = validateGrid(Dict(:x => [0.0, 0.5, 1.0]); allowPartial = true) + @test grid.topology == :nonuniform + + X = zeros(3, 4) + Y = zeros(3, 4) + grid2 = validateGrid(Dict(:m => 2, :n => 3, :nodes => (; X, Y))) + @test grid2.topology == :curvilinear + end + + @testset "boundary condition edge cases" begin + @test_throws ArgumentError makeGrid(m = 4, dx = 0.1, dc = [1.0, 1.0]) + @test_throws ArgumentError makeGrid(m = 4, dx = 0.1, nc = [1.0, 1.0]) + + @test_throws ArgumentError makeGrid( + m = 4, + dx = 0.1, + dc = [1.0, 1.0, 1.0], + nc = [1.0, 1.0, 1.0], + ) + + @test_throws ArgumentError makeGrid( + m = 4, + dx = 0.1, + dc = ones(2, 1), + nc = ones(2, 1), + ) + + grid = makeGrid(m = 4, dx = 0.1, dc = [1.0], nc = [0.0]) + @test grid.bc.dc == [1.0, 1.0] + @test grid.bc.nc == [0.0, 0.0] + @test grid.bc.isPeriodic == [false] + + grid2 = makeGrid( + m = 4, + n = 4, + dx = 0.1, + dy = 0.1, + bc = Dict(:dc => [0.0, 0.0, 1.0, 1.0], :nc => [0.0, 0.0, 1.0, 1.0]), + ) + @test grid2.bc.isPeriodic == [true, false] + end + + @testset "curvilinear validation errors" begin + bad_nodes_missing_y = (; X = zeros(5, 5)) + @test_throws ArgumentError makeGrid( + m = 4, + n = 4, + topology = :curvilinear, + nodes = bad_nodes_missing_y, + ) + + bad_nodes_wrong_x_size = (; X = zeros(4, 5), Y = zeros(5, 5)) + @test_throws ArgumentError makeGrid( + m = 4, + n = 4, + topology = :curvilinear, + nodes = bad_nodes_wrong_x_size, + ) + + bad_nodes_wrong_y_size = (; X = zeros(5, 5), Y = zeros(5, 4)) + @test_throws ArgumentError makeGrid( + m = 4, + n = 4, + topology = :curvilinear, + nodes = bad_nodes_wrong_y_size, + ) + end +end diff --git a/julia/MOLE.jl/test/runtests.jl b/julia/MOLE.jl/test/runtests.jl index 9dbd3da5..3e0d7593 100644 --- a/julia/MOLE.jl/test/runtests.jl +++ b/julia/MOLE.jl/test/runtests.jl @@ -1,6 +1,6 @@ using MOLE using Test, LinearAlgebra, SparseArrays -import MOLE: Operators, BCs +import MOLE: Operators, BCs, Grids @testset "Testing MOLE operators" begin @testset "Testing 1D divergence" begin @@ -26,3 +26,8 @@ end include("BCs/scalarBC.jl") end end + +@testset "Testing Grids" begin + include("Grids/make_grid.jl") + include("Grids/validate_grid.jl") +end