Functions

ComplianceFun(solver::AbstractFEASolver)

Differentiable structural compliance objective J = Fᵀ U = Σ ρ_e u_eᵀ K_e u_e.

Construct with ComplianceFun(solver). Call as comp(PseudoDensities(x)) where x is the filtered, optionally projected design. The closed-form gradient dJ/dx_e = -u_eᵀ K_e u_e · dρ_e/dx_e is propagated via a ChainRulesCore.rrule.

Only valid for StiffnessTopOptProblem with homogeneous Dirichlet BCs.

See [1] §2.1 for compliance minimization and [4] for the SIMP interpolation used in the stiffness assembly.

VolumeFun(solver::AbstractFEASolver; fraction=true)

Differentiable volume (or volume fraction) function. When fraction=true (default) returns V/total_volume; otherwise returns the absolute volume V = Σ x_e * cellvolume_e.

Construct with VolumeFun(solver). Call as vol(PseudoDensities(x)).

The volume constraint is the standard resource constraint in SIMP topology optimization; see [1] §2.1.

ThermalComplianceFun{T, TS<:AbstractFEASolver, TC<:AbstractVector{T}, TG<:AbstractVector{T}}

Thermal compliance objective function for heat transfer topology optimization.

Mathematical Formulation

For steady-state heat conduction with Dirichlet temperature v on ΓD and heat flux q on ΓN: -∇·(k(ρ)∇T) = 0 in Ω T = v on ΓD k∇T·n = q on ΓN

where k(ρ) = k_min + ρ^p (k_0 - k_min) (SIMP interpolation) and the heat source q is independent of ρ.

Thermal compliance is J = Q^T T where Q is the assembled load vector (the heat flux contribution on the free degrees of freedom) and T is the solved temperature. With homogeneous Dirichlet BCs (v = 0), J = T^T K(ρ) T; with inhomogeneous Dirichlet BCs the two differ and only Q^T T is the correct compliance (the T^T K T form leaks the prescribed-temperature energy).

Gradient (adjoint method): solve the condensed system K_cond λ = -Q_cond where Q_cond is Q zeroed on prescribed DOFs, then dJ/dxe = (λe^T Ke Te) · dρe/dx_e

Q is independent of ρ, so the standard adjoint applies; the only caveat is that the adjoint state λ is not -T when the Dirichlet values are nonzero (because the Dirichlet lift K(ρ) v enters the forward residual).

See [1] §1.3 and §4.1 for thermal topology optimization, and [2] for SIMP-based heat conduction optimization.

Usage

heatflux = Dict{String,Float64}("top" => 100.0)  # heat flux on boundary (W/m²)problem = HeatConductionProblem(Val{:Linear}, nels, sizes, k; Tleft=0.0, Tright=0.0, heatflux=heatflux)solver = FEASolver(DirectSolver, problem; xmin=0.001)comp = ThermalComplianceFun(solver)val = comp(PseudoDensities(ones(length(solver.vars))))

MeanComplianceFun(solver, scenarios)

Mean compliance over multiple load cases. Returns the average of the per-scenario compliance values. scenarios is typically produced by generate_scenarios.

BlockComplianceFun(solver, scenarios)

Per-load-case compliance vector for multi-load problems. Returns a vector of compliance values, one per load case.

DisplacementFun(solver::AbstractFEASolver)

Differentiable nodal displacement function. Solves the FEA system and returns the displacement vector u. Useful for displacement-constrained optimization and as input to stress functions.

Construct with DisplacementFun(solver). Call as u = disp(PseudoDensities(x)).

The adjoint-based gradient solves K λ = Δ (reusing the factorization) and computes du/dx_e = -dρ_e/dx_e · u_eᵀ K_e λ. See [1] §2.1 for the adjoint method in topology optimization.

StressTensorFun(solver::AbstractFEASolver)

Element-wise microscopic stress tensor. Computes the symmetric stress tensor for each element from the nodal displacements using the base Young's modulus. Call as σ = σf(u) where u is the displacement vector (e.g. from DisplacementFun). Returns a vector of symmetric matrices, one per element.

For 2D problems the stress is computed with the plane-strain constitutive law (σzz = λ·(εxx + εyy) ≠ 0), consistent with the stiffness matrix assembly. The returned tensor is 3×3 for 2D problems and 3×3 for 3D problems. See [5] §4.2 for the plane-strain formulation and [6] for stress-constrained topology optimization.

ElementStressTensorFun(solver::AbstractFEASolver)

Element stress tensor operator that also stores per-element metadata for use in stress-constrained optimization and ML applications.

von_mises_stress_function(solver::AbstractFEASolver; stress_exponent=0)

Return a function that computes the element-wise von Mises stress from the design. Applies the penalty and interpolation, solves the FEA system, and computes the von Mises stress for each element using the base Young's modulus.

Call as σv = σvf(PseudoDensities(x)). Returns a vector of von Mises stress values, one per element.

With stress_exponent = 0 (default) the returned stress is the microscopic stress σ = C_0 : ε of [6], which is finite at zero density and therefore exhibits the stress-singularity phenomenon: optimal designs with vanishing members are unreachable for gradient-based optimizers. With stress_exponent = q > 0 the returned stress is the relaxed (penalized) stress σ̃_e = ρ_e^q σ_e, where ρ_e is the physical density (xmin + (1 - xmin) x_e). The relaxation acts on the physical density, not the penalized stiffness density, so it is independent of the SIMP penalty and of the PENALTY_BEFORE_INTERPOLATION preference; the xmin floor keeps ρ_e^q finite and differentiable at zero density. Because σ̃_e → 0 as ρ_e → 0, low-density elements automatically satisfy any stress bound, which removes the singular optima. q = 0.5 with SIMP stiffness exponent p = 3 is the choice of [7]; the equivalent constraint-side formulation is the qp-approach of [8] with relaxation ε_qp = p - q. See also [9] for the alternative ε-relaxation (available here as epsilon_relaxed) and [10] for a review of both.

epsilon_relaxed(σv, ρ, σlim, ε)

Element-wise ε-relaxed stress constraint values (each value must be ≤ 0):

g_e = ρ_e (σ_e / σ_lim - 1) - ε

with element von Mises stresses σv, densities ρ (pass the raw filtered design variables, without the xmin floor: g is linear in ρ, so at ρ = 0 the constraint equals -ε exactly and the relaxation is independent of the void floor; the floor still enters through the computed stress field), stress limit σlim, and relaxation parameter ε > 0. Multiplying the constraint residual by ρ_e widens the degenerate subspaces of the feasible domain: for sufficiently small ρ_e the constraint is satisfied regardless of the stress, which makes the singular optima of stress-constrained topology optimization [11] reachable by gradient-based optimizers [9] and [6]. As ε → 0 the original constraints σ_e ≤ σ_lim on the material domain are recovered; a continuation strategy decreasing ε between optimization runs is common, though the global optimum of the relaxed problem can jump discontinuously with ε [12]. The qp-approach [8] is the alternative relaxation, available through the stress_exponent keyword of von_mises_stress_function. The returned values are signed, so aggregate them with a signed-safe smooth-max such as the Kreisselmeier– Steinhauser function (logsumexp(γ .* g) / γ), not a p-norm.

TrussStressFun(problem, solver)

Element-wise macroscopic axial stress for truss problems. Computes the axial stress in each truss member from nodal displacements and the penalized design. Call as σ = σf(u, ρ) where u is the displacement vector and ρ is the penalized/interpolated design.

The axial stress is computed as σ_e = -(R_e · Ke_e · u_e)[1] / A_e, where R_e is the local-to-global transformation matrix, Ke_e is the element stiffness matrix, and A_e is the cross-sectional area. Compressive stress is negative, tensile stress is positive. See [13] for the truss finite element formulation.

ElementKFun(solver::AbstractFEASolver)

Element-wise stiffness matrices from the input design variables. Applies the penalty and interpolation, then computes each element stiffness matrix via quadrature. Useful in buckling-constrained optimization.

Call as Kes = kesf(PseudoDensities(x)). Returns a vector of symmetric matrices, one per element.

AssembleKFun(problem)

Assemble element-wise matrices into a global sparse matrix. Useful in buckling-constrained optimization.

Call as K = assemble(Kes) where Kes is a vector of element matrices. Returns a SparseMatrixCSC.

apply_boundary_with_zerodiag!(Kσ, ch)

Apply boundary condition to a matrix. Zero-out the corresponding [i,:] and [:,j] with i, j ∈ ch.prescribed_dofs.

This function is typically used with the stress stiffness matrix Kσ. More info about this can be found at: https://github.com/JuliaTopOpt/TopOpt.jl/wiki/Applying-boundary-conditions-to-the-stress-stiffness-matrix

apply_boundary_with_meandiag!(K, ch)

Apply boundary condition to a matrix. Zero-out the corresponding [i,:] and [:,j] with i, j ∈ ch.prescribeddofs, then fill in K[i,i] for i ∈ ch.prescribeddofs with the mean diagonal of the original matrix.

TrussElementKσFun(problem, solver)

Element-wise stress/geometric stiffness matrices for truss domains. Used in buckling-constrained truss optimization.

Call as Kσs = Kσsf(u, ρ) where u is the displacement vector and ρ is the penalized/interpolated design. Returns a vector of symmetric matrices.

NeuralNetworkFun(nn, problem)

Re-parametrizes the design in terms of a neural network's weights and biases. nn is a Flux.jl model whose first layer takes 2 (or 3) coordinates for 2D (or 3D) and whose last layer returns a scalar in [0, 1]. In prediction mode the network is called on each element's centroid to produce that element's design variable.

TrainFunctionFun(nn_model)

The training function used in the re-parameterized topology optimization formulation. Takes the vector of neural-network weights/biases p and returns the vector of element-wise design variables x.

PredictFunctionFun(nn_model)

Prediction function that applies the neural network to each element centroid to produce the design variable. Used for evaluating a trained model.

MaterialInterpolationFun(values, penalty)

Maps a softmax over per-material decision variables to a physical material property (e.g. Young's modulus or density). values is a vector of material property values (length nmats, including void as the first entry). penalty is applied to the softmax output.

MultiMaterialVariablesFun(y, nmats)

Wraps the raw per-cell, per-material decision variables y (length ncells * (nmats - 1)) for use with MaterialInterpolationFun.

element_densities(mv::MultiMaterialVariablesFun)

Extract the per-element density vector from multi-material variables.

tounit(mv::MultiMaterialVariablesFun)

Convert MultiMaterialVariablesFun to unit-sum densities via softmax, so the per-cell material fractions sum to 1.

FixedElementProjectorFun

Callable struct that maps free design variables to a full density vector. Implements the projection with ChainRulesCore.rrule for automatic differentiation.

Fields

• nel::Int: Total number of elements
• black::BitVector: Fixed solid elements (density = 1)
• white::BitVector: Fixed void elements (density = 0)
• free::BitVector: Free design elements

get_fixed_element_projector(problem, black_cells, white_cells) -> FixedElementProjectorFun

Create a projector that maps free design variables to a full density vector.

Arguments

• problem: A TopOpt problem (e.g., PointLoadCantilever, HalfMBB, etc.)
• black_cells: Indices of elements fixed to solid (density = 1)
• white_cells: Indices of elements fixed to void (density = 0)

Returns

• FixedElementProjectorFun: A callable struct that maps x_free -> ρ_full

Example

problem = PointLoadCantilever((60, 20), (1.0, 1.0), 1.0, 0.3, 1.0)nel = prod(problem.nels)# Fix first row solid and last row voidblack = 1:problem.nels[1]  # First rowwhite = (nel - problem.nels[1] + 1):nel  # Last row# Create projectorprojector = get_fixed_element_projector(problem, black, white)# Use in optimizationx_free = fill(0.5, count(projector.free))  # Initialize free variablesρ = projector(x_free)  # Full density vector

get_fixed_element_projector(nel::Int, black_cells, white_cells) -> FixedElementProjectorFun

Create a projector given the total number of elements.

Arguments

• nel::Int: Total number of elements
• black_cells: Indices of elements fixed to solid (density = 1)
• white_cells: Indices of elements fixed to void (density = 0)

Returns

• FixedElementProjectorFun: A callable struct that maps x_free -> ρ_full

get_free_variables(p::FixedElementProjectorFun) -> BitVector

Return a BitVector indicating which elements are free design variables.

get_free_variable_count(p::FixedElementProjectorFun) -> Int

Return the number of free design variables.