From fd56e8314c2b3db0596a38d80b17cab232986ef5 Mon Sep 17 00:00:00 2001 From: Charles Kawczynski Date: Thu, 29 Feb 2024 09:05:47 -0500 Subject: [PATCH] Fix some warnings Fix Vararg warning Fix Base.Pairs warning Try to fix MultipleConditions --- src/utilities/convergence_condition.jl | 6 +- test/integrator_utils.jl | 6 +- test/problems.jl | 123 +++++++++++-------------- 3 files changed, 65 insertions(+), 70 deletions(-) diff --git a/src/utilities/convergence_condition.jl b/src/utilities/convergence_condition.jl index 46ed9ff1..a746acb2 100644 --- a/src/utilities/convergence_condition.jl +++ b/src/utilities/convergence_condition.jl @@ -92,8 +92,10 @@ updated_cache((; rate, order)::MinimumRateOfConvergence, cache, val, err, iter) Checks multiple `ConvergenceCondition`s, combining their results with either `all` or `any`. """ -struct MultipleConditions{CC <: Union{typeof(all), typeof(any)}, C <: Tuple{Vararg{<:ConvergenceCondition}}} <: - ConvergenceCondition +struct MultipleConditions{ + CC <: Union{typeof(all), typeof(any)}, + C <: Tuple{Vararg{T} where {T <: ConvergenceCondition}}, +} <: ConvergenceCondition condition_combiner::CC conditions::C end diff --git a/test/integrator_utils.jl b/test/integrator_utils.jl index 7c7646b8..a0becc8f 100644 --- a/test/integrator_utils.jl +++ b/test/integrator_utils.jl @@ -6,7 +6,11 @@ data structure has an instance of any types in `t`. """ function any_reltype(found, obj, name, ets, pc = (); warn = true) for pn in propertynames(obj) - prop = getproperty(obj, pn) + prop = if obj isa Base.Pairs + values(obj) + else + getproperty(obj, pn) + end pc_full = (pc..., ".", pn) pc_string = name * string(join(pc_full)) for et in ets diff --git a/test/problems.jl b/test/problems.jl index 7e42952f..0beccd20 100644 --- a/test/problems.jl +++ b/test/problems.jl @@ -1,16 +1,7 @@ using DiffEqBase, ClimaTimeSteppers, LinearAlgebra, StaticArrays using ClimaCore using ClimaComms -import ClimaCore.Domains as Domains -import ClimaCore.Geometry as Geometry -import ClimaCore.Meshes as Meshes -import ClimaCore.Topologies as Topologies -import ClimaCore.Spaces as Spaces -import ClimaCore.Fields as Fields -import ClimaCore.Operators as Operators - -import Krylov -Krylov.ktypeof(x::Fields.FieldVector) = ClimaComms.array_type(x){eltype(parent(x)), 1} +import ClimaCore: Domains, Geometry, Meshes, Topologies, Spaces, Fields, Operators, Limiters """ Single variable linear ODE @@ -449,6 +440,7 @@ Wfact!(W, Y, p, dtγ, t) = nothing 2D diffusion test problem. See [`2D diffusion problem`](@ref) for more details. """ function climacore_2Dheat_test_cts(::Type{FT}) where {FT} + context = ClimaComms.context() dss_tendency = true n_elem_x = 2 @@ -465,7 +457,7 @@ function climacore_2Dheat_test_cts(::Type{FT}) where {FT} Domains.IntervalDomain(Geometry.YPoint(FT(0)), Geometry.YPoint(FT(1)), periodic = true), ) mesh = Meshes.RectilinearMesh(domain, n_elem_x, n_elem_y) - topology = Topologies.Topology2D(mesh) + topology = Topologies.Topology2D(context, mesh) quadrature = Spaces.Quadratures.GLL{n_poly + 1}() space = Spaces.SpectralElementSpace2D(topology, quadrature) (; x, y) = Fields.coordinate_field(space) @@ -564,6 +556,7 @@ end # Implemented in flux form, with an optional limiter and hyperdiffusion. # TODO: Use this as an integration test. function deformational_flow_test(::Type{FT}; use_limiter = true, use_hyperdiffusion = true) where {FT} + context = ClimaComms.context() # Table III # Note: the paper uses "a" in place of "R" R = FT(6371220) # radius of Earth [m] @@ -593,30 +586,28 @@ function deformational_flow_test(::Type{FT}; use_limiter = true, use_hyperdiffus # hyperviscosity coefficient [m^4] (specified in the limiter paper) D₄ = FT(6.6e14) - centers = ClimaCore.Geometry.LatLongZPoint.(rad2deg(φ_c), rad2deg.((λ_c1, λ_c2)), FT(0)) + centers = Geometry.LatLongZPoint.(rad2deg(φ_c), rad2deg.((λ_c1, λ_c2)), FT(0)) # custom discretization (paper's discretization results in a very slow test) vert_nelems = 10 horz_nelems = 4 horz_npoly = 3 - vert_domain = ClimaCore.Domains.IntervalDomain( - ClimaCore.Geometry.ZPoint(FT(0)), - ClimaCore.Geometry.ZPoint(z_top); - boundary_names = (:bottom, :top), - ) - vert_mesh = ClimaCore.Meshes.IntervalMesh(vert_domain, nelems = vert_nelems) - vert_cent_space = ClimaCore.Spaces.CenterFiniteDifferenceSpace(vert_mesh) + vert_domain = + Domains.IntervalDomain(Geometry.ZPoint(FT(0)), Geometry.ZPoint(z_top); boundary_names = (:bottom, :top)) + vert_mesh = Meshes.IntervalMesh(vert_domain, nelems = vert_nelems) + z_topology = Topologies.IntervalTopology(context, vert_mesh) + vert_cent_space = Spaces.CenterFiniteDifferenceSpace(z_topology) - horz_domain = ClimaCore.Domains.SphereDomain(R) - horz_mesh = ClimaCore.Meshes.EquiangularCubedSphere(horz_domain, horz_nelems) - horz_topology = ClimaCore.Topologies.Topology2D(horz_mesh) - horz_quad = ClimaCore.Spaces.Quadratures.GLL{horz_npoly + 1}() - horz_space = ClimaCore.Spaces.SpectralElementSpace2D(horz_topology, horz_quad) + horz_domain = Domains.SphereDomain(R) + horz_mesh = Meshes.EquiangularCubedSphere(horz_domain, horz_nelems) + horz_topology = Topologies.Topology2D(context, horz_mesh) + horz_quad = Spaces.Quadratures.GLL{horz_npoly + 1}() + horz_space = Spaces.SpectralElementSpace2D(horz_topology, horz_quad) - cent_space = ClimaCore.Spaces.ExtrudedFiniteDifferenceSpace(horz_space, vert_cent_space) - cent_coords = ClimaCore.Fields.coordinate_field(cent_space) - face_space = ClimaCore.Spaces.FaceExtrudedFiniteDifferenceSpace(cent_space) + cent_space = Spaces.ExtrudedFiniteDifferenceSpace(horz_space, vert_cent_space) + cent_coords = Fields.coordinate_field(cent_space) + face_space = Spaces.FaceExtrudedFiniteDifferenceSpace(cent_space) # initial density (Equation 8) cent_ρ = @. p_0 / (R_d * T_0) * exp(-cent_coords.z / H) @@ -627,7 +618,7 @@ function deformational_flow_test(::Type{FT}; use_limiter = true, use_hyperdiffus φ = deg2rad(coord.lat) ds = map(centers) do center - r = ClimaCore.Geometry.great_circle_distance(coord, center, horz_space.global_geometry) + r = Geometry.great_circle_distance(coord, center, horz_space.global_geometry) return min(1, (r / R_t)^2 + ((z - z_c) / Z_t)^2) end in_slot = z > z_c && φ_c - FT(0.125) < φ < φ_c + FT(0.125) @@ -640,11 +631,11 @@ function deformational_flow_test(::Type{FT}; use_limiter = true, use_hyperdiffus return (; q1, q2, q3, q4, q5) end - init_state = ClimaCore.Fields.FieldVector(; cent_ρ, cent_ρq = cent_ρ .* cent_q) + init_state = Fields.FieldVector(; cent_ρ, cent_ρq = cent_ρ .* cent_q) # current wind vector (Equations 15--26) - current_cent_wind_vector = ClimaCore.Fields.Field(ClimaCore.Geometry.UVWVector{FT}, cent_space) - current_face_wind_vector = ClimaCore.Fields.Field(ClimaCore.Geometry.UVWVector{FT}, face_space) + current_cent_wind_vector = Fields.Field(Geometry.UVWVector{FT}, cent_space) + current_face_wind_vector = Fields.Field(Geometry.UVWVector{FT}, face_space) function wind_vector(coord, ρ, t) z = coord.z φ = deg2rad(coord.lat) @@ -667,12 +658,12 @@ function deformational_flow_test(::Type{FT}; use_limiter = true, use_hyperdiffus ω = ω_0 * sin(λ′) * cos(φ) * cos(2 * FT(π) * t / τ) * s w = -ω / (g * ρ) - return ClimaCore.Geometry.UVWVector(u, v, w) + return Geometry.UVWVector(u, v, w) end - horz_div = ClimaCore.Operators.Divergence() - horz_wdiv = ClimaCore.Operators.WeakDivergence() - horz_grad = ClimaCore.Operators.Gradient() + horz_div = Operators.Divergence() + horz_wdiv = Operators.WeakDivergence() + horz_grad = Operators.Gradient() cent_χ = similar(cent_q) function T_lim!(tendency, state, _, t) @. current_cent_wind_vector = wind_vector(cent_coords, state.cent_ρ, t) @@ -680,24 +671,21 @@ function deformational_flow_test(::Type{FT}; use_limiter = true, use_hyperdiffus @. tendency.cent_ρq = -horz_div(state.cent_ρq * current_cent_wind_vector) use_hyperdiffusion || return nothing @. cent_χ = horz_wdiv(horz_grad(state.cent_ρq / state.cent_ρ)) - ClimaCore.Spaces.weighted_dss!(cent_χ) + Spaces.weighted_dss!(cent_χ) @. tendency.cent_ρq += -D₄ * horz_wdiv(state.cent_ρ * horz_grad(cent_χ)) return nothing end - limiter = ClimaCore.Limiters.QuasiMonotoneLimiter(cent_q; rtol = FT(0)) + limiter = Limiters.QuasiMonotoneLimiter(cent_q; rtol = FT(0)) function lim!(state, _, t, ref_state) use_limiter || return nothing - ClimaCore.Limiters.compute_bounds!(limiter, ref_state.cent_ρq, ref_state.cent_ρ) - ClimaCore.Limiters.apply_limiter!(state.cent_ρq, state.cent_ρ, limiter) + Limiters.compute_bounds!(limiter, ref_state.cent_ρq, ref_state.cent_ρ) + Limiters.apply_limiter!(state.cent_ρq, state.cent_ρ, limiter) return nothing end - vert_div = ClimaCore.Operators.DivergenceF2C() - vert_interp = ClimaCore.Operators.InterpolateC2F( - top = ClimaCore.Operators.Extrapolate(), - bottom = ClimaCore.Operators.Extrapolate(), - ) + vert_div = Operators.DivergenceF2C() + vert_interp = Operators.InterpolateC2F(top = Operators.Extrapolate(), bottom = Operators.Extrapolate()) function T_exp!(tendency, state, _, t) @. current_face_wind_vector = wind_vector(face_coords, vert_interp(state.cent_ρ), t) @. tendency.cent_ρ = -vert_div(vert_interp(state.cent_ρ) * current_face_wind_vector) @@ -705,7 +693,7 @@ function deformational_flow_test(::Type{FT}; use_limiter = true, use_hyperdiffus end function dss!(state, _, t) - ClimaCore.Spaces.weighted_dss!(state.q) + Spaces.weighted_dss!(state.q) end function analytic_sol(t) @@ -731,6 +719,7 @@ end # (https://gmd.copernicus.org/articles/5/887/2012/gmd-5-887-2012.pdf) # Implemented in flux form, with an optional limiter and hyperdiffusion. function horizontal_deformational_flow_test(::Type{FT}; use_limiter = true, use_hyperdiffusion = true) where {FT} + context = ClimaComms.context() # constants (using the same notation as deformational_flow_test) R = FT(6371220) # radius of Earth [m] τ = 60 * 60 * 24 * FT(12) # period of motion (12 days) [s] @@ -740,18 +729,18 @@ function horizontal_deformational_flow_test(::Type{FT}; use_limiter = true, use_ R_t = R / 2 # horizontal half-width of tracers [m] D₄ = FT(6.6e14) # hyperviscosity coefficient [m^4] (specified in the limiter paper) - centers = ClimaCore.Geometry.LatLongPoint.(rad2deg(φ_c), rad2deg.((λ_c1, λ_c2))) + centers = Geometry.LatLongPoint.(rad2deg(φ_c), rad2deg.((λ_c1, λ_c2))) # 1.5° resolution on the equator: 360° / (4 * nelems * npoly) = 1.5° nelems = 20 npoly = 3 - domain = ClimaCore.Domains.SphereDomain(R) - mesh = ClimaCore.Meshes.EquiangularCubedSphere(domain, nelems) - topology = ClimaCore.Topologies.Topology2D(mesh) - quad = ClimaCore.Spaces.Quadratures.GLL{npoly + 1}() - space = ClimaCore.Spaces.SpectralElementSpace2D(topology, quad) - coords = ClimaCore.Fields.coordinate_field(space) + domain = Domains.SphereDomain(R) + mesh = Meshes.EquiangularCubedSphere(domain, nelems) + topology = Topologies.Topology2D(context, mesh) + quad = Spaces.Quadratures.GLL{npoly + 1}() + space = Spaces.SpectralElementSpace2D(topology, quad) + coords = Fields.coordinate_field(space) # initial conditions (Section 2.2) ρ = ones(space) @@ -760,15 +749,15 @@ function horizontal_deformational_flow_test(::Type{FT}; use_limiter = true, use_ λ = deg2rad(coord.long) hs = map(centers) do center - center′ = ClimaCore.Geometry.CartesianPoint(center, space.global_geometry) - coord′ = ClimaCore.Geometry.CartesianPoint(coord, space.global_geometry) + center′ = Geometry.CartesianPoint(center, space.global_geometry) + coord′ = Geometry.CartesianPoint(coord, space.global_geometry) dist_squared = (coord′.x1 - center′.x1)^2 + (coord′.x2 - center′.x2)^2 + (coord′.x3 - center′.x3)^2 # Note: the paper doesn't divide by R^2, which only works if R = 1 return FT(0.95) * exp(-5 * dist_squared / R^2) end gaussian_hills = hs[1] + hs[2] rs = map(centers) do center - return ClimaCore.Geometry.great_circle_distance(coord, center, space.global_geometry) + return Geometry.great_circle_distance(coord, center, space.global_geometry) end cosine_bells = if rs[1] < R_t FT(0.1) + FT(0.9) * (1 + cos(FT(π) * rs[1] / R_t)) / 2 @@ -791,10 +780,10 @@ function horizontal_deformational_flow_test(::Type{FT}; use_limiter = true, use_ return (; gaussian_hills, cosine_bells, slotted_cylinders) end - init_state = ClimaCore.Fields.FieldVector(; ρ, ρq = ρ .* q) + init_state = Fields.FieldVector(; ρ, ρq = ρ .* q) # current wind vector (Section 2.3) - current_wind_vector = ClimaCore.Fields.Field(ClimaCore.Geometry.UVVector{FT}, space) + current_wind_vector = Fields.Field(Geometry.UVVector{FT}, space) function wind_vector(coord, t) φ = deg2rad(coord.lat) λ = deg2rad(coord.long) @@ -805,12 +794,12 @@ function horizontal_deformational_flow_test(::Type{FT}; use_limiter = true, use_ u = k * sin(λ′)^2 * sin(2 * φ) * cos(FT(π) * t / τ) + 2 * FT(π) * R / τ * cos(φ) v = k * sin(2 * λ′) * cos(φ) * cos(FT(π) * t / τ) - return ClimaCore.Geometry.UVVector(u, v) + return Geometry.UVVector(u, v) end - div = ClimaCore.Operators.Divergence() - wdiv = ClimaCore.Operators.WeakDivergence() - grad = ClimaCore.Operators.Gradient() + div = Operators.Divergence() + wdiv = Operators.WeakDivergence() + grad = Operators.Gradient() χ = similar(q) function T_lim!(tendency, state, _, t) @. current_wind_vector = wind_vector(coords, t) @@ -818,22 +807,22 @@ function horizontal_deformational_flow_test(::Type{FT}; use_limiter = true, use_ @. tendency.ρq = -div(state.ρq * current_wind_vector) use_hyperdiffusion || return nothing @. χ = wdiv(grad(state.ρq / state.ρ)) - ClimaCore.Spaces.weighted_dss!(χ) + Spaces.weighted_dss!(χ) @. tendency.ρq += -D₄ * wdiv(state.ρ * grad(χ)) return nothing end - limiter = ClimaCore.Limiters.QuasiMonotoneLimiter(q; rtol = FT(0)) + limiter = Limiters.QuasiMonotoneLimiter(q; rtol = FT(0)) function lim!(state, _, t, ref_state) use_limiter || return nothing - ClimaCore.Limiters.compute_bounds!(limiter, ref_state.ρq, ref_state.ρ) - ClimaCore.Limiters.apply_limiter!(state.ρq, state.ρ, limiter) + Limiters.compute_bounds!(limiter, ref_state.ρq, ref_state.ρ) + Limiters.apply_limiter!(state.ρq, state.ρ, limiter) return nothing end function dss!(state, _, t) - ClimaCore.Spaces.weighted_dss!(state.ρ) - ClimaCore.Spaces.weighted_dss!(state.ρq) + Spaces.weighted_dss!(state.ρ) + Spaces.weighted_dss!(state.ρq) end function analytic_sol(t)