Private types and functions

compute_deviation_from_isotropy(X, Y, Z; radius = 1)

Compute a scalar measure of the deviation from isotropy for a spherical grid (e.g., a conformal cubed-sphere panel), defined by the Cartesian coordinate arrays X, Y, and Z. Each of X, Y, and Z is a 2D array of size (Nx, Ny) holding the Cartesian coordinates of the grid vertices on the sphere, such that the point at (i, j) corresponds to (X[i, j], Y[i, j], Z[i, j]). The grid therefore contains (Nx−1) × (Ny−1) spherical quadrilateral cells.

For each quadrilateral cell, the function computes the great-circle distances of its four edges on the sphere, evaluates the sum of absolute differences between consecutive edge lengths as a measure of cell anisotropy, and then returns the Euclidean norm of these deviations over the entire grid.

Arguments

• X, Y, Z: (Nx, Ny) arrays of Cartesian coordinates of grid vertices on the sphere.

Keyword Arguments

• radius: Sphere radius. Default: 1.

Returns

• A non-negative scalar quantifying the overall deviation from isotropy in the grid. Larger values correspond to more anisotropic grids.

Examples

A regular latitude-longitude grid:

using CubedSphere: compute_deviation_from_isotropy

Nx, Ny = 3, 3
lons = range(-π/4, π/4, length = Nx)
lats = range(-π/6, π/6, length = Ny)

X = [cos(φ) * cos(λ) for λ in lons, φ in lats]
Y = [cos(φ) * sin(λ) for λ in lons, φ in lats]
Z = [sin(φ)          for λ in lons, φ in lats]

compute_deviation_from_isotropy(X, Y, Z)

# output

1.6552138747243959

The vertices of a tetrahedron:

using CubedSphere: compute_deviation_from_isotropy

V = [   0     0     1 ;
     2√2/3    0   -1/3;
     -√2/3  √6/3  -1/3;
     -√2/3 -√6/3  -1/3]

X = reshape(V[:, 1], 2, 2)
Y = reshape(V[:, 2], 2, 2)
Z = reshape(V[:, 3], 2, 2)

isapprox(compute_deviation_from_isotropy(X, Y, Z), 0, atol=1e-14)

# output

true

compute_model_diagnostics(X, Y, Z, minimum_reference_cell_area)

Compute the two model diagnostics for a given mapped grid:

1. Normalized minimum cell width: sqrt(min(cell_area) / minimum_reference_cell_area), using spherical cell areas from compute_cell_areas(X, Y, Z) and a uniform-grid reference area.
2. Deviation from isotropy: the heuristic anisotropy measure from compute_deviation_from_isotropy(X, Y, Z).

Arguments

• X, Y, Z: (Nx, Ny) Cartesian coordinate arrays on the sphere.
• minimum_reference_cell_area: The minimum cell area of a reference (typically uniform) grid.

Returns

• A 2-element vector [normalized_minimum_cell_width, deviation_from_isotropy].

compute_weighted_model_diagnostics(model_diagnostics)

Apply weights from specify_weights_for_model_diagnostics() to the unweighted diagnostics.

Arguments

• model_diagnostics: A 2-element vector [normalized_minimum_cell_width, deviation_from_isotropy].

Returns

• A 2-element vector with weights applied, in the same order as the inputs.

cube_to_sphere(x, y)

Maps the coordinates $(x, y) ∈ [-1, 1] × [-1, 1]$ of the face of a cube to the 3D coordinates $X, Y, Z$ on the sphere via conformal_cubed_sphere_mapping, that is:

X[i, j], Y[i, j], Z[i, j] = conformal_cubed_sphere_mapping(x[i], y[j])

exponential_spacing(N, k₀ByN)

Construct a symmetric set of N face locations on [-1, 1] with exponentially graded spacing away from the domain center. Let k₀ = k₀ByN * N, and define an exponential map on the right half, x(t) = α·exp(t/k₀) + β, anchored so that it passes through (t₀, 0) and (t₁, 1), then mirror about zero. Endpoints are fixed at -1 and +1.

• Odd N (M = (N+1)/2): a face lies at 0 (faces[M] = 0). The right-half faces use t = 1, …, M with x(1) = 0, x(M) = 1, and the left half is the negative mirror.
• Even N (M = N/2): no face at 0. The two central faces straddle zero, with 0 midway between them; the right-half faces use t = 2, …, M+1 anchored by x(1.5) = 0, x(M+1) = 1, and the left half is mirrored.

Arguments

• N: Number of faces (N ≥ 2).
• k₀ByN: Grading parameter scaled by N (the code uses k₀ = k₀ByN * N). Larger k₀ByN yields spacing closer to uniform; smaller k₀ByN increases clustering near the center. Requires k₀ByN > 0.

Returns

• faces: A length-N strictly increasing vector of face coordinates on [-1, 1] with symmetry faces[i] = - faces[N+1-i], and end-points faces[1] = -1, faces[N] = 1.

forward_map(Nx, Ny, spacing, θ)

Evaluate the (weighted) model diagnostics for a non-uniform conformal cubed-sphere panel defined by parameters θ. The steps are:

1. Clamp θ to parameter limits from specify_parameter_limits(spacing).
2. Build a reference uniform grid via conformal_cubed_sphere_coordinates(Nx, Ny; spacing=UniformSpacing()) and compute the minimum_reference_cell_area.
3. Build the non-uniform grid via conformal_cubed_sphere_coordinates(Nx, Ny; spacing, ...), passing the parameters in θ according to spacing type.
4. Compute model diagnostics and then apply weights.

Arguments

• Nx, Ny: Number of grid vertices along the panel coordinates.
• spacing: GeometricSpacing() (interprets θ[1] as ratio^(Nx-1)) or ExponentialSpacing() (interprets θ[1] as k₀/Nx).
• θ: Parameter vector (one element in this implementation).

Returns

• A 2-element vector of weighted model diagnostics.

geometric_spacing(N, ratio_raised_to_N_minus_one)

Construct a symmetric set of N face locations on the interval [-1, 1] with geometrically graded spacing away from the domain center. Let r be the geometric ratio recovered from ratio_raised_to_N_minus_one via r = (ratio_raised_to_N_minus_one)^(1/(N-1)). The gaps between consecutive faces grow by a factor of r as one moves outward from the center, and the layout is mirrored about zero. Endpoints are fixed at -1 and +1.

• Odd N: A face lies at 0. Outward gaps are Δx, Δx*r, Δx*r^2, …, mirrored about 0, and chosen so the rightmost face lands at +1.
• Even N: No face at 0. Central faces are at ±Δx/2; outward gaps are Δx*r, Δx*r^2, …, mirrored, and chosen so the rightmost face lands at +1.

Arguments

• N: Number of faces (N ≥ 2).
• ratio_raised_to_N_minus_one: The value r^(N-1) for some geometric ratio r > 0. Values near 1 yield nearly uniform spacing. (Exactly 1 is not supported by the closed-form formulas used here.)

Returns

• faces: A length-N monotonically increasing vector of face coordinates on [-1, 1] with geometric grading and symmetry: faces[1] = -1, faces[N] = 1, and faces[i] = -faces[N+1-i].

optimize!(Nx, Ny, spacing, θ; N_iterations=10, Δt=1)

Run an Ensemble Kalman Inversion (EKI) to tune the spacing parameter(s) for a non-uniform conformal cubed sphere panel. An ensemble of parameter vectors θ is iteratively updated so that the weighted model diagnostics produced by forward_map match the ideal targets from specify_ideal_weighted_model_diagnostics().

At each iteration:

1. Forward evaluations (parallelized): For each ensemble member θ[n], compute the predicted diagnostics G[n] = forward_map(Nx, Ny, spacing, θ[n]).
2. Ensemble statistics: Form means θ̄ = mean(θ) and G̅ = mean(G), then compute
   • Cross-covariance Cᵘᵖ = cov(θ, G) (shape: nθ × ndata), and
   • Data covariance Cᵖᵖ = cov(G, G) (shape: ndata × ndata),
   using the unbiased (N_ensemble-1) denominator.
3. Perturbed observations: For each member, build y[n] = ideal + Δt * η[n] where η[n] ~ N(0, I). Here Δt sets the perturbation magnitude of the observations.
4. Residuals: r[n] = y[n] - G[n].
5. Implicit update (Kalman-like step): Update each parameter vector via θ[n] ← θ[n] + K * r[n], with K = Cᵘᵖ * (Cᵖᵖ + I/Δt)⁻¹, implemented by solving the linear system with a Cholesky factorization of Cᵖᵖ + I/Δt. The same Δt also acts as an implicit damping/step-size control: smaller Δt ⇒ stronger regularization and smaller updates; larger Δt ⇒ weaker regularization and larger, noisier updates.
6. Monitoring: Report error = ‖mean(r)‖ and store a snapshot of the ensemble.

This function mutates the input ensemble θ in place (it becomes the final ensemble) and records the full ensemble after each iteration.

Arguments

• Nx, Ny: Number of grid vertices along the panel coordinates.
• spacing: GeometricSpacing() or ExponentialSpacing().
• θ: Initial ensemble — a vector of one-element parameter vectors [θ] (length N_ensemble).

Keyword Arguments

• N_iterations: Number of EKI iterations. Default: 10
• Δt: Pseudo-time step that sets both the observation perturbation scale and the implicit damping in the linear solve (Cᵖᵖ + I/Δt).

Returns

• θ_series: A length N_iterations + 1 vector; each entry is a snapshot of the full ensemble (index 1 is the initial ensemble; the last is the final ensemble). The input θ is also mutated to the final state.

specify_ideal_weighted_model_diagnostics()

Return the target (ideal) values for the weighted model diagnostics used by the inversion. By default, the target normalized minimum cell width is 4 and the target deviation from isotropy is 0, then weights are applied in the same order.

Returns

• A 2-element vector of target weighted diagnostics.

specify_parameter_limits(spacing)

Return the lower and upper bounds for the single spacing parameter as a 2-element vector [min, max]. For a consistent interface with possible multi-parameter extensions, this vector is returned wrapped in a one-element array: [[min, max]]. The bounds depend on spacing:

• GeometricSpacing() returns [[min, max]] for ratio^(N-1).
• ExponentialSpacing() returns [[min, max]] for k₀/N.

Returns [θ_limits] where θ_limits == [min, max] for the single parameter in this implementation.

specify_parameters(spacing)

Return a one-element vector [θ] containing the initial guess for the spacing parameters used to build non-uniform conformal cubed sphere panels for different spacings.

• GeometricSpacing(): uses a single parameter interpreted downstream as ratio^(N-1) for geometric grading.
• ExponentialSpacing(): uses a single parameter interpreted downstream as k₀/N for exponential grading.

specify_random_parameters(N_ensemble, spacing)

Draw an ensemble of random parameter vectors within the limits from specify_parameter_limits(spacing). Each ensemble member is sampled uniformly within its parameter bounds.

Arguments

• N_ensemble: Number of ensemble members to generate.
• spacing: GeometricSpacing() or ExponentialSpacing().

Returns

• A vector of length N_ensemble, where each element is a one-element parameter vector [θ] (a single parameter in this implementation).

specify_weights_for_model_diagnostics()

Return the weights applied to the model diagnostics used by the objective function. The two diagnostics are, in order: (1) normalized minimum cell width, (2) deviation from isotropy.

Returns

• A 2-element vector of weights, e.g., [10, 1].