using TopOpt, LinearAlgebra
using StatsFuns: logsumexpGlobal Stress: Relaxation and Aggregation for Stress Constraints
Description
This tutorial solves the classical stress-constrained topology optimization problem: minimize volume subject to a von Mises stress limit,
\[\min_{\rho}\; V(\rho) \quad \text{s.t.} \quad \sigma^{vm}_e \le \sigma_{\lim} \;\; \forall e\]
Stress constraints are fundamentally harder than compliance constraints for two independent reasons (Verbart et al., 2017; Le et al., 2010):
- Singular optima. The microscopic stress \(\sigma = C_0 : \varepsilon\) stays finite as an element’s density \(\rho \to 0\), so the stress constraint blocks any gradient-driven path that removes a load-bearing member: the true optima lie in degenerate subspaces that gradient methods cannot reach (Sved & Ginos, 1968; Cheng & Guo, 1997).
- Locality. Stress is an element-wise state variable, giving thousands of local constraints — as many as there are design variables — which is expensive to handle directly (see the
local_stresstutorial).
The standard toolkit, demonstrated here, combines:
- Relaxation to remove the singularity: the qp-approach uses a relaxed stress \(\tilde\sigma_e = \rho_e^{q}\,\sigma_e\) with \(q < p\) (Bruggi, 2008; Le et al., 2010 use \(q = 1/2\) with SIMP \(p = 3\)), while ε-relaxation replaces the constraint by \(\rho_e(\sigma_e/\sigma_{\lim} - 1) \le \epsilon\) (Cheng & Guo, 1997; Duysinx & Bendsøe, 1998). Both make the constraint automatically satisfied at sufficiently low density.
- Aggregation to collapse the local constraints into one global smooth measure: the p-norm (Duysinx & Sigmund, 1998) or the Kreisselmeier–Steinhauser (KS) function (Yang & Chen, 1996), with adaptive normalization to control the true maximum (Le et al., 2010).
We compare three formulations on the classical L-bracket benchmark and show that all produce the signature feature of stress-constrained design: the rounded re-entrant corner, which a compliance-driven design keeps sharp.
Setup
The benchmark problem
The L-bracket: fixed top edge, downward load at the end of the horizontal arm. The re-entrant corner is a stress concentrator, which is what makes this the standard stress benchmark.
A point force on a continuum produces a singular stress field that the mesh regularizes to ~f / element_size; the resulting “peak stress” is a mesh-dependent load-introduction artifact, not a property of the topology. We therefore distribute the load over 5 nodes with load_width, which removes the singularity. (Alternatives used in the literature: exclude a fixed region around the load from the aggregation, or use regional stress measures.)
E, ν, f = 1.0, 0.3, 1.0
problem = LBeam(
Val{:Linear};
length=50, height=50, upperslab=25, lowerslab=25, E=E, ν=ν, force=f,
load_width=5, # distribute the load over 5 nodes
)
xmin = 1e-3
# Fixed SIMP exponent p = 3, as standard for stress problems (Le et al., 2010)
solver = FEASolver(DirectSolver, problem; xmin=xmin, penalty=PowerPenaltyFun(3.0))
filter = DensityFilterFun(solver; rmin=2.0)
volfrac = VolumeFun(solver)
N = length(solver.vars)Microscopic vs. relaxed stress
von_mises_stress_function computes the microscopic stress \(\sigma = C_0 : \varepsilon\) with the base Young’s modulus (the physically consistent choice for porous materials; Duysinx & Bendsøe, 1998). Because it is finite at zero density, it exhibits the singularity problem. Passing stress_exponent = q > 0 returns instead the relaxed stress \(\tilde\sigma_e = \rho_e^q \sigma_e\), which vanishes at zero density and thereby makes singular optima reachable. The relaxation acts on the physical density (not the penalized stiffness density), so it is independent of the SIMP exponent p used for the stiffness.
σ_micro = von_mises_stress_function(solver)
σ_rel = von_mises_stress_function(solver; stress_exponent=0.5)
x_solid = ones(N)
σ_solid = σ_rel(filter(PseudoDensities(x_solid)))
println("peak stress of the solid L-bracket: ", maximum(σ_solid))peak stress of the solid L-bracket: 0.5984185129578113The stress limit
The stress limit is a material property, not a multiple of the solid design’s average stress: fully-stressed low-volume designs run hotter than the solid block. Here we set σ_lim to 1.75× the solid block’s peak, which places the optimum at roughly 40–60% volume — the interesting regime for this benchmark.
σlim = 1.75 * maximum(σ_solid)
obj = x -> volfrac(filter(PseudoDensities(x)))
alg = MMA87(; dualoptimizer=ConjugateGradient())
function run_chunks(
constr, x0, chunks, maxiter; update! = (x, j) -> nothing, label = "", keep_best = true
)
model = Model(obj)
addvar!(model, zeros(N), ones(N))
add_ineq_constraint!(model, constr)
x = copy(x0)
for j in 1:chunks
update!(x, j)
tol = j == chunks ? Nonconvex.Tolerance(; kkt=1e-4) : Nonconvex.Tolerance(; kkt=1e-3)
options = MMAOptions(; maxiter, tol, keep_best)
r = optimize(model, alg, x; options)
x = r.minimizer
println("$label chunk $j: V=$(round(obj(x), digits=3)), constr=$(round(constr(x), digits=4))")
end
return x
endApproach 1: relaxed stress + normalized p-norm + adaptive normalization
The normalized p-norm \(\sigma_{PN} = \big(\frac{1}{N}\sum_e \tilde\sigma_e^P\big)^{1/P}\) is a lower bound on the maximum stress (equality only for a uniform stress state), so it cannot be constrained directly. Following Le et al. (2010), we rescale it by an adaptive factor \(c \approx \max(\tilde\sigma)/\sigma_{PN}\) — updated between optimization chunks — so that \(c\,\sigma_{PN}\) tracks the true maximum. The p-norm’s gradient is spread over all moderately stressed elements, which makes this the most robust of the three formulations in practice.
P = 8
c = Ref(maximum(σ_solid) / (norm(σ_solid, P) / N^(1 / P)))
constr_pnorm = x -> begin
s = σ_rel(filter(PseudoDensities(x)))
c[] * (norm(s, P) / N^(1 / P)) / σlim - 1
end
# Adaptive normalization: c tracks max(σ̃)/σ_PN between chunks
update_pnorm! = x -> begin
s = σ_rel(filter(PseudoDensities(x)))
c[] = maximum(s) / (norm(s, P) / N^(1 / P))
end
x1 = run_chunks(constr_pnorm, x_solid, 4, 70; update! = (x, j) -> update_pnorm!(x), label = "p-norm")[ Info: iter obj Δobj violation kkt_residual [ Info: 0 1.0e+00 Inf 0.0e+00 4.3e-01 [ Info: 1 6.6e-01 3.4e-01 1.9e-01 1.3e-02 [ Info: 2 4.8e-01 1.8e-01 3.1e+00 4.1e-01 [ Info: 3 1.7e-01 3.1e-01 3.1e+01 3.5e-01 [ Info: 4 6.0e-02 1.1e-01 7.0e+01 1.1e-01 [ Info: 5 2.7e-01 2.1e-01 4.3e+01 2.6e+38 [ Info: 6 4.6e-01 1.9e-01 2.7e+01 8.3e+36 [ Info: 7 2.1e-01 2.5e-01 5.1e+01 3.4e-01 [ Info: 8 8.1e-02 1.3e-01 7.2e+01 3.5e-02 [ Info: 9 3.6e-02 4.5e-02 9.4e+01 7.7e-02 [ Info: 10 3.0e-01 2.7e-01 2.4e+01 4.2e+37 [ Info: 11 1.8e-01 1.3e-01 3.4e+01 5.3e-01 [ Info: 12 2.3e-01 5.5e-02 2.9e+01 1.9e+00 [ Info: 13 1.8e-01 5.5e-02 2.9e+01 1.3e-01 [ Info: 14 3.3e-01 1.5e-01 1.5e+01 5.9e+37 [ Info: 15 3.7e-01 4.4e-02 1.2e+01 7.6e+36 [ Info: 16 4.1e-01 3.6e-02 1.4e+01 4.4e+36 [ Info: 17 4.5e-01 4.1e-02 4.7e+00 5.7e+37 [ Info: 18 5.0e-01 4.9e-02 3.7e+00 4.0e+36 [ Info: 19 5.0e-01 7.7e-03 2.6e+00 1.7e+36 [ Info: 20 5.2e-01 1.9e-02 1.9e+00 1.2e+36 [ Info: 21 5.6e-01 3.4e-02 1.3e+00 1.0e+36 [ Info: 22 5.7e-01 1.3e-02 9.8e-01 5.5e+35 [ Info: 23 5.9e-01 2.2e-02 8.1e-01 4.8e+35 [ Info: 24 6.1e-01 1.4e-02 6.2e-01 3.4e+35 [ Info: 25 6.1e-01 8.8e-03 6.4e-01 3.9e+35 [ Info: 26 6.2e-01 9.7e-03 5.4e-01 2.5e+35 [ Info: 27 6.3e-01 4.1e-03 5.1e-01 1.9e+35 [ Info: 28 6.4e-01 1.1e-02 3.8e-01 1.5e+35 [ Info: 29 6.5e-01 5.7e-03 3.1e-01 1.1e+35 [ Info: 30 6.5e-01 8.9e-03 2.5e-01 7.9e+34 [ Info: 31 6.7e-01 1.2e-02 2.0e-01 5.4e+34 [ Info: 32 6.8e-01 1.2e-02 1.6e-01 3.5e+34 [ Info: 33 6.9e-01 7.9e-03 1.2e-01 1.9e+34 [ Info: 34 7.0e-01 1.3e-02 7.1e-02 4.5e+35 [ Info: 35 7.1e-01 1.1e-02 2.9e-02 3.7e+34 [ Info: 36 7.2e-01 1.1e-02 0.0e+00 2.5e+35 [ Info: 37 6.8e-01 4.2e-02 0.0e+00 3.7e-03 [ Info: 38 6.6e-01 2.3e-02 0.0e+00 4.6e-03 [ Info: 39 6.3e-01 2.1e-02 0.0e+00 4.5e-03 [ Info: 40 6.2e-01 1.7e-02 5.2e-03 3.2e-03 [ Info: 41 6.0e-01 1.7e-02 2.6e-01 1.4e-01 [ Info: 42 4.8e-01 1.2e-01 1.3e+00 3.8e-02 [ Info: 43 3.7e-01 1.1e-01 2.7e+00 4.7e-02 [ Info: 44 3.0e-01 7.0e-02 3.6e+00 6.9e-02 [ Info: 45 4.5e-01 1.5e-01 1.2e+00 6.7e+00 [ Info: 46 4.8e-01 3.5e-02 9.2e-01 1.7e+00 [ Info: 47 5.8e-01 9.2e-02 3.6e-01 4.6e+35 [ Info: 48 6.3e-01 5.6e-02 1.6e-01 6.1e+35 [ Info: 49 6.7e-01 4.2e-02 4.9e-02 7.5e+34 [ Info: 50 7.0e-01 2.3e-02 0.0e+00 1.6e+35 [ Info: 51 6.4e-01 6.1e-02 0.0e+00 5.8e-03 [ Info: 52 6.1e-01 2.6e-02 0.0e+00 5.8e-03 [ Info: 53 5.9e-01 2.5e-02 0.0e+00 5.7e-03 [ Info: 54 5.7e-01 1.8e-02 0.0e+00 5.3e-03 [ Info: 55 5.5e-01 1.4e-02 0.0e+00 5.2e-03 [ Info: 56 5.4e-01 1.2e-02 0.0e+00 4.2e-03 [ Info: 57 5.2e-01 1.7e-02 2.6e-02 1.0e-02 [ Info: 58 5.0e-01 2.2e-02 4.1e-01 1.2e-01 [ Info: 59 3.8e-01 1.2e-01 3.3e+00 6.9e-02 [ Info: 60 2.7e-01 1.2e-01 6.4e+00 1.8e-02 [ Info: 61 1.8e-01 8.1e-02 1.3e+01 2.7e-02 [ Info: 62 1.5e-01 3.6e-02 3.9e+01 2.3e-01 [ Info: 63 1.2e-01 2.4e-02 3.5e+01 2.4e-01 [ Info: 64 9.7e-02 2.7e-02 3.8e+01 1.1e-01 [ Info: 65 3.3e-01 2.3e-01 1.6e+01 1.0e+37 [ Info: 66 2.6e-01 7.6e-02 2.3e+01 4.6e-01 [ Info: 67 3.5e-01 9.2e-02 2.0e+01 4.5e+37 [ Info: 68 3.8e-01 3.5e-02 1.1e+01 1.3e+37 [ Info: 69 2.7e-01 1.1e-01 1.7e+01 7.9e-01 [ Info: 70 2.1e-01 5.7e-02 2.2e+01 2.7e-01 p-norm chunk 1: V=0.679, constr=-0.0054 [ Info: iter obj Δobj violation kkt_residual [ Info: 0 6.8e-01 Inf 0.0e+00 3.5e-01 [ Info: 1 4.9e-01 1.9e-01 3.5e-01 3.8e-02 [ Info: 2 2.2e-01 2.7e-01 7.6e+00 1.3e-01 [ Info: 3 1.3e-01 9.3e-02 2.9e+01 2.7e-01 [ Info: 4 3.7e-01 2.4e-01 1.6e+01 1.3e+38 [ Info: 5 4.5e-01 7.7e-02 1.1e+01 1.2e+37 [ Info: 6 4.7e-01 2.4e-02 1.4e+01 7.8e+37 [ Info: 7 3.2e-01 1.6e-01 2.3e+01 8.5e+00 [ Info: 8 1.6e-01 1.6e-01 2.3e+01 3.7e-02 [ Info: 9 9.2e-02 6.4e-02 4.4e+01 3.1e-01 [ Info: 10 2.1e-01 1.2e-01 4.2e+01 1.5e+37 [ Info: 11 1.1e-01 1.0e-01 3.8e+01 1.9e-01 [ Info: 12 1.5e-01 4.4e-02 4.6e+01 3.8e+00 [ Info: 13 7.1e-02 8.3e-02 4.2e+01 3.9e-02 [ Info: 14 2.0e-01 1.3e-01 2.4e+01 7.4e+37 [ Info: 15 3.3e-01 1.2e-01 1.7e+01 4.4e+37 [ Info: 16 3.8e-01 4.9e-02 6.9e+00 8.8e+36 [ Info: 17 4.9e-01 1.1e-01 3.3e+00 5.9e+37 [ Info: 18 5.4e-01 5.7e-02 1.6e+00 1.6e+36 [ Info: 19 6.2e-01 7.1e-02 6.8e-01 7.8e+36 [ Info: 20 6.7e-01 5.7e-02 2.7e-01 2.1e+37 [ Info: 21 7.2e-01 5.0e-02 0.0e+00 6.6e+35 [ Info: 22 6.1e-01 1.2e-01 3.2e-01 6.7e-02 [ Info: 23 5.0e-01 1.1e-01 9.1e-01 6.6e-02 [ Info: 24 4.7e-01 2.8e-02 1.2e+00 1.9e-01 [ Info: 25 3.9e-01 8.1e-02 2.1e+00 1.7e-01 [ Info: 26 3.3e-01 6.3e-02 2.7e+00 9.4e-02 [ Info: 27 4.3e-01 1.0e-01 2.4e+00 7.2e+00 [ Info: 28 4.1e-01 1.6e-02 1.9e+00 6.8e-01 [ Info: 29 4.1e-01 3.5e-03 1.3e+00 4.8e-01 [ Info: 30 4.9e-01 8.2e-02 8.6e-01 9.3e+36 [ Info: 31 5.3e-01 4.0e-02 3.8e-01 4.1e+35 [ Info: 32 5.9e-01 5.9e-02 1.4e-01 4.6e+35 [ Info: 33 6.3e-01 3.9e-02 0.0e+00 1.6e+36 [ Info: 34 5.9e-01 3.7e-02 7.5e-04 4.7e-03 [ Info: 35 5.7e-01 2.3e-02 0.0e+00 5.0e-03 [ Info: 36 5.5e-01 1.9e-02 0.0e+00 3.1e-03 [ Info: 37 5.3e-01 1.8e-02 0.0e+00 3.2e-03 [ Info: 38 5.2e-01 1.4e-02 0.0e+00 2.6e-03 [ Info: 39 5.1e-01 1.4e-02 0.0e+00 2.6e-03 [ Info: 40 4.9e-01 1.3e-02 0.0e+00 2.1e-03 [ Info: 41 4.8e-01 1.3e-02 0.0e+00 2.2e-03 [ Info: 42 4.7e-01 1.3e-02 6.1e-03 1.5e-03 [ Info: 43 4.6e-01 1.0e-02 1.6e-02 9.3e-03 [ Info: 44 4.4e-01 2.2e-02 7.1e-02 1.2e-02 [ Info: 45 4.3e-01 1.4e-03 2.8e-02 7.5e-03 [ Info: 46 4.3e-01 1.7e-03 3.7e-03 1.6e-03 [ Info: 47 4.3e-01 3.8e-03 0.0e+00 1.5e-03 [ Info: 48 4.2e-01 4.0e-03 0.0e+00 1.2e-03 [ Info: 49 4.2e-01 3.7e-03 0.0e+00 1.2e-03 [ Info: 50 4.2e-01 3.4e-03 0.0e+00 1.1e-03 [ Info: 51 4.1e-01 3.0e-03 0.0e+00 1.0e-03 [ Info: 52 4.1e-01 2.6e-03 0.0e+00 9.5e-04 p-norm chunk 2: V=0.412, constr=-0.0021 [ Info: iter obj Δobj violation kkt_residual [ Info: 0 4.1e-01 Inf 1.2e-02 1.2e-02 [ Info: 1 4.1e-01 1.5e-03 0.0e+00 1.1e-03 [ Info: 2 4.0e-01 7.6e-03 9.2e-03 2.7e-03 [ Info: 3 3.9e-01 1.4e-02 6.3e-02 1.3e-02 [ Info: 4 3.7e-01 2.2e-02 6.8e-01 1.5e-01 [ Info: 5 2.1e-01 1.5e-01 4.2e+00 7.2e-02 [ Info: 6 1.8e-01 3.4e-02 1.2e+01 2.9e-01 [ Info: 7 8.3e-02 9.6e-02 5.7e+01 1.1e-01 [ Info: 8 2.8e-01 1.9e-01 4.0e+01 6.3e+38 [ Info: 9 5.2e-01 2.4e-01 6.1e+00 1.4e+37 [ Info: 10 3.1e-01 2.1e-01 2.0e+01 3.0e-01 [ Info: 11 1.2e-01 1.9e-01 4.8e+01 4.3e-02 [ Info: 12 6.7e-02 5.4e-02 4.3e+01 7.4e-02 [ Info: 13 2.2e-01 1.5e-01 3.3e+01 1.4e+38 [ Info: 14 3.8e-01 1.6e-01 2.5e+01 7.5e+38 [ Info: 15 2.5e-01 1.2e-01 2.5e+01 2.4e+00 [ Info: 16 1.2e-01 1.4e-01 5.1e+01 8.0e-02 [ Info: 17 1.7e-01 5.2e-02 3.2e+01 4.2e+00 [ Info: 18 1.8e-01 1.2e-02 3.0e+01 8.5e-01 [ Info: 19 1.1e-01 7.2e-02 4.2e+01 4.1e-01 [ Info: 20 2.3e-01 1.2e-01 2.4e+01 5.5e+37 [ Info: 21 3.1e-01 8.2e-02 1.6e+01 1.0e+37 [ Info: 22 4.2e-01 1.0e-01 1.1e+01 2.3e+37 [ Info: 23 4.6e-01 4.0e-02 7.5e+00 2.4e+36 [ Info: 24 4.9e-01 3.1e-02 5.9e+00 3.2e+36 [ Info: 25 5.2e-01 3.8e-02 3.9e+00 3.3e+36 [ Info: 26 5.4e-01 1.1e-02 2.8e+00 2.7e+36 [ Info: 27 5.7e-01 3.8e-02 2.1e+00 1.9e+36 [ Info: 28 6.0e-01 2.2e-02 1.6e+00 1.6e+36 [ Info: 29 6.1e-01 1.5e-02 1.2e+00 9.8e+35 [ Info: 30 6.3e-01 1.9e-02 8.5e-01 5.9e+35 [ Info: 31 6.4e-01 1.2e-02 6.6e-01 3.2e+35 [ Info: 32 6.6e-01 1.8e-02 4.7e-01 1.7e+35 [ Info: 33 6.7e-01 1.4e-02 3.2e-01 9.8e+34 [ Info: 34 6.9e-01 1.7e-02 1.9e-01 3.6e+34 [ Info: 35 7.0e-01 1.4e-02 9.5e-02 4.8e+35 [ Info: 36 7.2e-01 1.2e-02 8.0e-03 5.6e+35 [ Info: 37 6.6e-01 5.9e-02 1.5e-03 9.3e-03 [ Info: 38 6.3e-01 3.1e-02 0.0e+00 6.9e-03 [ Info: 39 6.0e-01 3.1e-02 0.0e+00 5.1e-03 [ Info: 40 5.7e-01 2.2e-02 0.0e+00 3.6e-03 [ Info: 41 5.5e-01 2.0e-02 0.0e+00 2.5e-03 [ Info: 42 5.4e-01 1.8e-02 3.9e-03 2.0e-03 [ Info: 43 5.2e-01 1.7e-02 6.4e-02 2.1e-02 [ Info: 44 4.7e-01 4.9e-02 5.7e-01 1.0e-01 [ Info: 45 3.9e-01 7.5e-02 1.2e+00 7.7e-02 [ Info: 46 3.2e-01 7.8e-02 2.2e+00 4.0e-02 [ Info: 47 3.7e-01 5.3e-02 1.1e+00 6.3e-01 [ Info: 48 3.9e-01 1.8e-02 7.6e-01 4.8e-01 [ Info: 49 4.7e-01 8.7e-02 4.1e-01 5.4e+35 [ Info: 50 5.4e-01 6.2e-02 2.3e-01 5.3e+37 [ Info: 51 5.3e-01 1.1e-02 5.3e-02 3.5e-01 [ Info: 52 5.1e-01 2.0e-02 9.7e-03 1.5e-02 [ Info: 53 4.9e-01 1.2e-02 0.0e+00 6.7e-03 [ Info: 54 4.8e-01 1.3e-02 0.0e+00 3.6e-03 [ Info: 55 4.7e-01 9.9e-03 0.0e+00 2.7e-03 [ Info: 56 4.6e-01 7.0e-03 0.0e+00 2.4e-03 [ Info: 57 4.6e-01 5.2e-03 0.0e+00 2.3e-03 [ Info: 58 4.5e-01 5.0e-03 0.0e+00 2.4e-03 [ Info: 59 4.5e-01 4.9e-03 0.0e+00 2.5e-03 [ Info: 60 4.4e-01 4.9e-03 0.0e+00 2.4e-03 [ Info: 61 4.4e-01 5.1e-03 0.0e+00 2.5e-03 [ Info: 62 4.3e-01 5.5e-03 0.0e+00 1.9e-03 [ Info: 63 4.3e-01 5.3e-03 5.7e-04 2.1e-03 [ Info: 64 4.2e-01 5.5e-03 2.0e-03 1.4e-03 [ Info: 65 4.2e-01 4.1e-03 1.4e-03 1.4e-03 [ Info: 66 4.2e-01 3.9e-03 2.3e-03 1.1e-03 [ Info: 67 4.1e-01 3.4e-03 7.7e-03 5.6e-03 [ Info: 68 4.0e-01 1.0e-02 4.1e-02 9.1e-03 [ Info: 69 4.0e-01 1.6e-03 2.6e-02 7.1e-03 [ Info: 70 4.0e-01 3.0e-03 5.5e-03 2.6e-03 p-norm chunk 3: V=0.419, constr=0.0014 [ Info: iter obj Δobj violation kkt_residual [ Info: 0 4.2e-01 Inf 0.0e+00 4.8e-02 [ Info: 1 4.0e-01 2.0e-02 2.7e-02 5.5e-03 [ Info: 2 3.7e-01 2.9e-02 5.8e-01 8.1e-02 [ Info: 3 2.5e-01 1.2e-01 4.1e+00 1.4e-01 [ Info: 4 1.2e-01 1.3e-01 3.5e+01 3.9e-01 [ Info: 5 3.4e-02 8.9e-02 4.7e+01 1.5e-02 [ Info: 6 1.9e-01 1.6e-01 3.4e+01 4.4e+38 [ Info: 7 2.4e-01 4.8e-02 4.1e+01 2.4e+38 [ Info: 8 3.3e-01 9.3e-02 2.6e+01 9.2e+37 [ Info: 9 2.1e-01 1.3e-01 3.4e+01 3.1e-01 [ Info: 10 6.0e-02 1.5e-01 5.6e+01 3.1e-02 [ Info: 11 8.6e-02 2.6e-02 4.4e+01 4.5e-01 [ Info: 12 5.3e-02 3.3e-02 4.1e+01 5.2e-02 [ Info: 13 2.4e-01 1.9e-01 2.9e+01 1.9e+38 [ Info: 14 2.9e-01 5.0e-02 1.7e+01 6.3e+37 [ Info: 15 4.3e-01 1.4e-01 6.8e+00 2.2e+36 [ Info: 16 5.0e-01 6.9e-02 4.4e+00 3.7e+37 [ Info: 17 5.3e-01 3.4e-02 2.1e+00 1.5e+37 [ Info: 18 5.3e-01 2.8e-03 2.7e+00 2.6e+02 [ Info: 19 5.8e-01 4.6e-02 1.1e+00 1.1e+37 [ Info: 20 6.0e-01 1.9e-02 5.0e-01 3.6e+36 [ Info: 21 6.3e-01 3.6e-02 2.1e-01 4.5e+35 [ Info: 22 6.6e-01 2.8e-02 2.5e-02 1.9e+35 [ Info: 23 5.9e-01 7.2e-02 1.1e-01 5.8e-02 [ Info: 24 5.0e-01 8.5e-02 4.3e-01 5.6e-02 [ Info: 25 4.2e-01 8.5e-02 1.5e+00 1.0e-01 [ Info: 26 2.7e-01 1.5e-01 9.3e+00 1.1e-01 [ Info: 27 1.7e-01 9.8e-02 3.7e+01 1.6e-01 [ Info: 28 9.5e-02 7.7e-02 4.0e+01 3.0e-02 [ Info: 29 1.6e-01 6.4e-02 2.7e+01 9.2e+00 [ Info: 30 8.0e-02 8.0e-02 4.4e+01 1.1e-01 [ Info: 31 1.6e-01 8.0e-02 4.3e+01 1.2e+38 [ Info: 32 3.2e-01 1.6e-01 6.9e+00 1.8e+37 [ Info: 33 4.2e-01 9.6e-02 1.6e+01 7.1e+36 [ Info: 34 2.4e-01 1.7e-01 1.3e+01 3.5e-01 [ Info: 35 3.5e-01 1.0e-01 1.2e+01 1.6e+01 [ Info: 36 2.2e-01 1.2e-01 8.8e+00 1.3e-01 [ Info: 37 3.3e-01 1.1e-01 4.6e+00 2.8e+36 [ Info: 38 3.8e-01 4.6e-02 3.9e+00 4.9e+36 [ Info: 39 4.2e-01 4.2e-02 2.5e+00 1.8e+36 [ Info: 40 4.6e-01 3.6e-02 1.5e+00 8.1e+35 [ Info: 41 5.0e-01 4.3e-02 1.0e+00 5.1e+35 [ Info: 42 5.3e-01 3.2e-02 6.2e-01 2.6e+35 [ Info: 43 5.7e-01 3.2e-02 3.7e-01 1.3e+35 [ Info: 44 6.0e-01 3.2e-02 2.2e-01 5.5e+34 [ Info: 45 6.3e-01 3.3e-02 1.1e-01 8.9e+35 [ Info: 46 6.6e-01 2.7e-02 2.4e-02 5.3e+34 [ Info: 47 6.3e-01 3.2e-02 0.0e+00 2.3e-02 [ Info: 48 5.9e-01 3.2e-02 0.0e+00 7.1e-03 [ Info: 49 5.7e-01 2.0e-02 0.0e+00 5.0e-03 [ Info: 50 5.6e-01 1.6e-02 0.0e+00 4.1e-03 [ Info: 51 5.5e-01 1.2e-02 0.0e+00 3.3e-03 [ Info: 52 5.3e-01 1.1e-02 0.0e+00 2.9e-03 [ Info: 53 5.2e-01 1.2e-02 0.0e+00 3.1e-03 [ Info: 54 5.1e-01 1.1e-02 0.0e+00 3.5e-03 [ Info: 55 5.0e-01 9.5e-03 0.0e+00 3.6e-03 [ Info: 56 4.9e-01 9.5e-03 5.4e-03 3.4e-03 [ Info: 57 4.8e-01 8.2e-03 9.6e-02 4.0e-02 [ Info: 58 3.8e-01 1.1e-01 1.8e+00 1.0e-01 [ Info: 59 3.1e-01 6.9e-02 2.7e+00 5.6e-02 [ Info: 60 2.7e-01 4.3e-02 5.6e+00 1.4e-01 [ Info: 61 1.7e-01 1.0e-01 2.9e+01 6.3e-02 [ Info: 62 1.2e-01 4.9e-02 4.0e+01 7.2e-02 [ Info: 63 9.4e-02 2.2e-02 3.6e+01 1.4e-01 [ Info: 64 2.0e-01 1.0e-01 3.1e+01 1.2e+37 [ Info: 65 3.3e-01 1.3e-01 2.8e+01 5.2e+37 [ Info: 66 1.9e-01 1.4e-01 2.8e+01 9.3e-01 [ Info: 67 3.3e-01 1.4e-01 1.1e+01 1.3e+37 [ Info: 68 3.5e-01 1.5e-02 1.3e+01 1.5e+37 [ Info: 69 4.3e-01 7.8e-02 7.1e+00 5.1e+36 [ Info: 70 4.5e-01 2.3e-02 3.7e+00 3.3e+36 p-norm chunk 4: V=0.534, constr=-0.0061
1875-element Vector{Float64}:
0.0
0.0
0.0
0.07623110389593508
0.2232857692913714
0.21707349475438906
0.3540673405955637
0.4968341098489188
0.6491224826967136
0.7559387738006419
⋮
0.09769509451099298
0.12298587679853494
0.23244167994006232
0.642117046522064
0.5949334360297507
0.5827938557009922
0.8961734661157196
0.9840552560534952
1.0Approach 2: relaxed stress + KS aggregation
The KS function \(\sigma_{KS} = \max(s) + \log\sum_e e^{\gamma(s_e - \max(s))}/\gamma\) is an upper bound on the maximum, with slack at most \(\log(N)/\gamma\) — so no normalization factor is needed, and \(\sigma_{KS} \le 1\) certifies \(\max(\tilde\sigma) \le \sigma_{\lim}\). The price is conservatism (the design below keeps extra material) and sharper gradients: large γ tightens the bound but concentrates the gradient on the few peak elements, which destabilizes MMA. We use a moderate γ = 20 and warm-start from the p-norm design.
γ = 20.0
constr_ks = x -> begin
s = σ_rel(filter(PseudoDensities(x))) ./ σlim
logsumexp(γ .* s) / γ - 1
end
# keep_best=false: the warm start is infeasible under the conservative KS
# measure, and we want MMA to move (it converges to feasibility). A single
# chunk: restarting MMA between chunks on the sharp KS constraint makes the
# trajectory oscillate.
x2 = run_chunks(constr_ks, x1, 1, 80; label = "KS", keep_best = false)[ Info: iter obj Δobj violation kkt_residual [ Info: 0 5.3e-01 Inf 1.9e-01 1.9e-01 [ Info: 1 4.1e-01 1.2e-01 1.0e+01 4.1e+00 [ Info: 2 1.7e-01 2.4e-01 4.4e+01 2.3e-02 [ Info: 3 3.0e-01 1.3e-01 6.9e+01 6.9e+37 [ Info: 4 1.2e-01 1.8e-01 5.4e+01 8.6e-02 [ Info: 5 2.4e-01 1.2e-01 5.0e+01 1.2e+39 [ Info: 6 1.0e-01 1.4e-01 6.4e+01 1.7e-02 [ Info: 7 2.5e-01 1.5e-01 3.8e+01 2.4e+38 [ Info: 8 7.9e-02 1.7e-01 5.7e+01 1.3e-02 [ Info: 9 2.5e-01 1.7e-01 2.0e+01 6.5e+37 [ Info: 10 1.5e-01 9.3e-02 3.7e+01 4.9e-02 [ Info: 11 8.7e-02 6.6e-02 7.8e+01 1.4e-01 [ Info: 12 2.4e-01 1.6e-01 2.9e+01 7.1e+37 [ Info: 13 1.1e-01 1.3e-01 4.4e+01 3.0e-02 [ Info: 14 1.7e-01 6.4e-02 4.3e+01 5.7e+38 [ Info: 15 2.6e-01 8.0e-02 2.1e+01 2.1e+38 [ Info: 16 2.3e-01 2.5e-02 2.8e+01 1.6e+38 [ Info: 17 2.0e-01 3.0e-02 2.3e+01 4.2e+38 [ Info: 18 2.4e-01 4.2e-02 1.9e+01 9.1e+36 [ Info: 19 2.3e-01 6.9e-03 2.1e+01 1.7e+38 [ Info: 20 2.5e-01 2.0e-02 1.4e+01 1.3e+38 [ Info: 21 2.6e-01 3.3e-03 1.2e+01 7.9e+37 [ Info: 22 2.7e-01 1.1e-02 2.1e+01 1.2e+38 [ Info: 23 2.7e-01 2.7e-03 9.9e+00 3.4e+37 [ Info: 24 2.9e-01 2.0e-02 8.9e+00 1.2e+37 [ Info: 25 3.1e-01 1.6e-02 8.8e+00 9.8e+36 [ Info: 26 3.2e-01 1.5e-02 8.7e+00 3.3e+37 [ Info: 27 3.3e-01 1.0e-02 8.7e+00 1.6e+37 [ Info: 28 3.4e-01 9.0e-03 6.9e+00 2.6e+37 [ Info: 29 3.5e-01 4.3e-03 6.7e+00 2.5e+37 [ Info: 30 3.3e-01 1.4e-02 7.7e+00 5.2e+37 [ Info: 31 3.5e-01 1.7e-02 7.0e+00 2.9e+37 [ Info: 32 3.5e-01 2.9e-03 6.8e+00 2.5e+37 [ Info: 33 3.6e-01 1.2e-02 6.9e+00 3.1e+37 [ Info: 34 3.8e-01 1.8e-02 6.3e+00 2.8e+37 [ Info: 35 3.9e-01 1.2e-02 8.8e+00 3.3e+37 [ Info: 36 3.9e-01 6.1e-04 5.2e+00 2.3e+37 [ Info: 37 4.1e-01 1.7e-02 5.6e+00 1.6e+37 [ Info: 38 4.0e-01 6.2e-03 4.6e+00 3.2e+37 [ Info: 39 4.1e-01 6.0e-03 4.2e+00 2.4e+37 [ Info: 40 4.3e-01 1.7e-02 3.7e+00 2.2e+37 [ Info: 41 4.1e-01 1.7e-02 3.7e+00 3.0e+37 [ Info: 42 4.3e-01 1.4e-02 3.3e+00 7.8e+36 [ Info: 43 4.3e-01 5.3e-03 3.1e+00 1.2e+37 [ Info: 44 4.4e-01 5.7e-03 3.0e+00 1.4e+37 [ Info: 45 4.4e-01 6.1e-03 2.9e+00 8.4e+36 [ Info: 46 4.5e-01 7.2e-03 2.7e+00 1.0e+37 [ Info: 47 4.6e-01 8.1e-03 2.4e+00 1.1e+37 [ Info: 48 4.7e-01 9.9e-03 2.5e+00 8.1e+36 [ Info: 49 4.7e-01 5.6e-03 2.4e+00 7.8e+36 [ Info: 50 4.8e-01 3.2e-03 2.1e+00 7.2e+36 [ Info: 51 4.9e-01 1.2e-02 1.9e+00 5.0e+36 [ Info: 52 5.0e-01 1.0e-02 1.7e+00 4.7e+36 [ Info: 53 5.1e-01 1.0e-02 1.6e+00 4.0e+36 [ Info: 54 5.2e-01 1.0e-02 1.4e+00 3.2e+36 [ Info: 55 5.3e-01 1.1e-02 1.3e+00 2.3e+36 [ Info: 56 5.4e-01 1.1e-02 1.2e+00 2.1e+36 [ Info: 57 5.5e-01 1.1e-02 1.1e+00 1.8e+36 [ Info: 58 5.6e-01 9.8e-03 1.0e+00 1.6e+36 [ Info: 59 5.7e-01 7.7e-03 8.5e-01 1.3e+36 [ Info: 60 5.8e-01 1.5e-02 7.2e-01 8.6e+35 [ Info: 61 6.0e-01 1.3e-02 6.2e-01 7.1e+35 [ Info: 62 6.1e-01 1.3e-02 5.4e-01 5.5e+35 [ Info: 63 6.2e-01 1.1e-02 4.6e-01 4.6e+35 [ Info: 64 6.2e-01 4.6e-03 4.5e-01 4.0e+35 [ Info: 65 6.3e-01 4.4e-03 3.7e-01 3.0e+35 [ Info: 66 6.4e-01 1.1e-02 3.1e-01 2.0e+35 [ Info: 67 6.5e-01 6.9e-03 2.6e-01 1.8e+35 [ Info: 68 6.5e-01 5.8e-03 2.2e-01 1.1e+35 [ Info: 69 6.6e-01 8.4e-03 1.7e-01 7.0e+34 [ Info: 70 6.7e-01 8.1e-03 1.2e-01 3.6e+34 [ Info: 71 6.8e-01 9.4e-03 8.1e-02 1.6e+34 [ Info: 72 6.9e-01 8.4e-03 4.6e-02 4.5e+35 [ Info: 73 7.0e-01 7.6e-03 1.5e-02 2.5e+36 [ Info: 74 6.9e-01 8.6e-03 4.0e-03 4.8e-02 [ Info: 75 6.7e-01 1.3e-02 0.0e+00 9.9e-03 [ Info: 76 6.6e-01 1.5e-02 0.0e+00 2.6e-02 [ Info: 77 6.4e-01 1.7e-02 0.0e+00 2.8e-02 [ Info: 78 6.2e-01 1.7e-02 0.0e+00 3.2e-02 [ Info: 79 6.1e-01 1.7e-02 0.0e+00 1.8e-02 [ Info: 80 5.9e-01 1.7e-02 0.0e+00 3.5e-02 KS chunk 1: V=0.59, constr=-0.0038
1875-element Vector{Float64}:
0.09459443733681551
0.1325106848119762
0.21137418615177767
0.0
0.0
0.24084202805839341
0.3637660253818712
0.3771145481584302
0.49139757093442543
0.5845368271802315
⋮
0.9673471412905384
0.9932125667087237
1.0
1.0
0.9712072788947137
0.8846871881585193
0.5124608054226102
0.43935379275067865
0.4203507563773462Approach 3: ε-relaxation of the microscopic stress
ε-relaxation acts on the constraint rather than the stress measure: \(g_e = \rho_e(\sigma_e/\sigma_{\lim} - 1) - \epsilon \le 0\), with the true microscopic stress. At zero density the constraint reads \(-\epsilon \le 0\) — always satisfied — and as \(\epsilon \to 0\) the original local constraints are recovered. We pass the raw filtered densities (no xmin floor): g is linear in ρ, which keeps it AD-safe and floor-independent. The relaxed values are signed, so we aggregate them with the KS function (a p-norm is inappropriate for signed values); dividing by ε first keeps the KS slack \(\log(N)/\gamma\) small relative to ε. We run a single fixed ε = 0.25. Driving ε → 0 tightens the true stress bound \(\sigma_{\lim}(1 + \epsilon)\) on solid elements, but needs a slow continuation with many iterations — the global optimum of the relaxed problem can jump discontinuously as ε → 0 (Stolpe & Svanberg, 2001) — which is why production codes pair ε-relaxation with augmented-Lagrangian solvers (Fancello & Pereira, 2006) rather than plain MMA.
ε = 0.25
constr_eps = x -> begin
xfil = filter(PseudoDensities(x))
g = epsilon_relaxed(σ_micro(xfil), xfil.x, σlim, ε) ./ ε
logsumexp(γ .* g) / γ
end
x3 = run_chunks(constr_eps, x_solid, 1, 120; label = "ε-relaxed")[ Info: iter obj Δobj violation kkt_residual [ Info: 0 1.0e+00 Inf 0.0e+00 2.7e+00 [ Info: 1 6.0e-01 4.0e-01 9.0e-01 1.4e-02 [ Info: 2 1.9e-01 4.1e-01 5.0e+01 1.8e-01 [ Info: 3 9.3e-03 1.8e-01 1.3e+02 5.1e-02 [ Info: 4 8.8e-02 7.9e-02 8.3e+01 2.2e+39 [ Info: 5 1.5e-01 5.9e-02 1.3e+02 3.9e+39 [ Info: 6 2.3e-01 7.9e-02 7.4e+01 3.5e+39 [ Info: 7 2.5e-01 2.4e-02 9.9e+01 1.6e+39 [ Info: 8 9.0e-02 1.6e-01 8.0e+01 1.9e-01 [ Info: 9 2.8e-01 1.9e-01 3.6e+01 1.7e+38 [ Info: 10 1.2e-01 1.6e-01 6.6e+01 7.0e-02 [ Info: 11 6.4e-02 5.7e-02 8.2e+01 2.2e-01 [ Info: 12 2.5e-01 1.9e-01 4.8e+01 6.6e+38 [ Info: 13 3.8e-01 1.3e-01 3.7e+01 9.1e+37 [ Info: 14 3.5e-01 3.1e-02 2.8e+01 6.4e+37 [ Info: 15 4.1e-01 5.8e-02 2.7e+01 5.0e+38 [ Info: 16 3.4e-01 6.6e-02 3.1e+01 5.2e+38 [ Info: 17 2.8e-01 5.8e-02 2.8e+01 2.1e+00 [ Info: 18 1.9e-01 9.5e-02 4.4e+01 8.3e-02 [ Info: 19 1.3e-01 6.2e-02 8.8e+01 1.2e+00 [ Info: 20 9.0e-02 3.5e-02 1.1e+02 1.2e-01 [ Info: 21 5.9e-02 3.1e-02 9.8e+01 2.0e-02 [ Info: 22 2.0e-01 1.4e-01 5.6e+01 8.9e+38 [ Info: 23 2.1e-01 1.5e-02 7.9e+01 6.6e+38 [ Info: 24 2.8e-01 6.2e-02 6.9e+01 1.7e+38 [ Info: 25 3.2e-01 4.8e-02 3.3e+01 6.2e+38 [ Info: 26 3.4e-01 1.6e-02 3.0e+01 1.4e+38 [ Info: 27 3.2e-01 1.7e-02 2.3e+01 4.9e+38 [ Info: 28 3.4e-01 1.5e-02 3.0e+01 1.5e+38 [ Info: 29 3.5e-01 9.2e-03 2.2e+01 5.2e+38 [ Info: 30 3.2e-01 2.4e-02 2.6e+01 4.1e+38 [ Info: 31 3.5e-01 2.8e-02 2.3e+01 3.0e+37 [ Info: 32 3.9e-01 4.1e-02 2.7e+01 1.3e+38 [ Info: 33 4.3e-01 3.8e-02 2.3e+01 1.7e+38 [ Info: 34 4.4e-01 1.5e-02 1.9e+01 4.8e+38 [ Info: 35 4.7e-01 2.9e-02 1.5e+01 1.8e+38 [ Info: 36 4.8e-01 4.9e-03 1.3e+01 1.3e+38 [ Info: 37 5.0e-01 2.2e-02 1.1e+01 1.1e+38 [ Info: 38 4.7e-01 3.2e-02 1.4e+01 1.8e+38 [ Info: 39 4.9e-01 2.1e-02 1.4e+01 5.7e+37 [ Info: 40 5.0e-01 1.5e-02 1.1e+01 8.7e+37 [ Info: 41 4.6e-01 4.0e-02 9.1e+00 1.3e+38 [ Info: 42 5.1e-01 4.3e-02 9.4e+00 2.3e+37 [ Info: 43 5.3e-01 1.8e-02 7.7e+00 5.8e+37 [ Info: 44 5.0e-01 2.9e-02 8.0e+00 9.3e+37 [ Info: 45 5.3e-01 3.3e-02 7.1e+00 2.8e+37 [ Info: 46 5.5e-01 1.8e-02 7.0e+00 3.1e+37 [ Info: 47 5.5e-01 4.7e-03 6.4e+00 4.7e+37 [ Info: 48 5.6e-01 1.0e-02 5.5e+00 2.8e+37 [ Info: 49 5.6e-01 9.0e-04 5.0e+00 3.1e+37 [ Info: 50 5.7e-01 2.5e-03 4.6e+00 2.7e+37 [ Info: 51 5.8e-01 1.1e-02 5.0e+00 2.3e+37 [ Info: 52 5.5e-01 2.5e-02 5.7e+00 4.3e+37 [ Info: 53 5.6e-01 7.2e-03 4.9e+00 1.8e+37 [ Info: 54 5.8e-01 1.8e-02 3.8e+00 1.2e+37 [ Info: 55 5.8e-01 3.8e-03 3.7e+00 1.6e+37 [ Info: 56 5.9e-01 8.9e-03 3.3e+00 1.3e+37 [ Info: 57 6.0e-01 8.0e-03 3.4e+00 1.6e+37 [ Info: 58 5.9e-01 3.6e-03 2.8e+00 2.0e+37 [ Info: 59 6.1e-01 1.3e-02 2.5e+00 8.8e+36 [ Info: 60 6.1e-01 1.9e-03 2.2e+00 9.8e+36 [ Info: 61 6.2e-01 1.2e-02 1.9e+00 3.8e+36 [ Info: 62 6.3e-01 4.9e-03 1.7e+00 4.2e+36 [ Info: 63 6.3e-01 7.8e-03 1.4e+00 3.0e+36 [ Info: 64 6.4e-01 1.1e-02 1.1e+00 1.8e+36 [ Info: 65 6.5e-01 9.4e-03 9.0e-01 1.3e+36 [ Info: 66 6.6e-01 7.0e-03 7.6e-01 2.1e+36 [ Info: 67 6.7e-01 7.9e-03 5.3e-01 4.8e+35 [ Info: 68 6.8e-01 1.0e-02 3.5e-01 2.1e+35 [ Info: 69 6.9e-01 8.6e-03 1.9e-01 8.9e+34 [ Info: 70 7.0e-01 9.3e-03 4.9e-02 1.8e+36 [ Info: 71 6.9e-01 6.2e-03 0.0e+00 3.1e-02 [ Info: 72 6.8e-01 1.2e-02 0.0e+00 1.0e-02 [ Info: 73 6.7e-01 1.3e-02 0.0e+00 2.7e-02 [ Info: 74 6.5e-01 1.7e-02 1.5e-03 8.1e-03 [ Info: 75 6.3e-01 1.7e-02 1.4e-02 3.0e-02 [ Info: 76 6.1e-01 1.9e-02 1.1e-01 7.1e-02 [ Info: 77 5.6e-01 4.9e-02 6.8e-01 2.0e-02 [ Info: 78 5.6e-01 6.0e-03 6.0e-01 3.2e-01 [ Info: 79 5.2e-01 3.5e-02 2.0e+00 2.8e-01 [ Info: 80 5.1e-01 1.3e-02 1.7e+00 8.0e+35 [ Info: 81 4.9e-01 2.3e-02 2.7e+00 1.5e+00 [ Info: 82 4.9e-01 3.2e-03 2.2e+00 8.5e+35 [ Info: 83 4.9e-01 1.2e-03 3.0e+00 2.0e+36 [ Info: 84 5.1e-01 1.7e-02 4.0e+00 2.5e+36 [ Info: 85 4.8e-01 2.4e-02 3.2e+00 8.3e+36 [ Info: 86 5.2e-01 3.7e-02 1.7e+00 9.0e+35 [ Info: 87 5.5e-01 3.2e-02 1.1e+00 1.1e+36 [ Info: 88 5.8e-01 2.8e-02 6.5e-01 2.3e+35 [ Info: 89 5.9e-01 1.1e-02 7.1e-01 2.9e+35 [ Info: 90 6.1e-01 1.8e-02 6.5e-01 5.9e+36 [ Info: 91 6.3e-01 2.0e-02 2.8e-01 8.4e+34 [ Info: 92 6.5e-01 1.4e-02 7.1e-02 5.8e+35 [ Info: 93 6.4e-01 9.7e-03 6.7e-03 4.9e-02 [ Info: 94 6.1e-01 2.1e-02 4.1e-02 8.0e-02 [ Info: 95 6.0e-01 2.0e-02 6.4e-02 1.3e-02 [ Info: 96 5.9e-01 7.7e-03 2.3e-02 1.9e-02 [ Info: 97 5.8e-01 1.2e-02 1.7e-02 2.0e-02 [ Info: 98 5.6e-01 1.2e-02 0.0e+00 1.1e-02 [ Info: 99 5.5e-01 1.3e-02 3.0e-03 2.9e-02 [ Info: 100 5.4e-01 1.1e-02 0.0e+00 1.1e-02 [ Info: 101 5.3e-01 1.1e-02 0.0e+00 2.1e-02 [ Info: 102 5.2e-01 1.1e-02 7.7e-04 1.5e-02 [ Info: 103 5.1e-01 1.1e-02 9.7e-03 2.8e-02 [ Info: 104 5.0e-01 9.1e-03 3.2e-02 7.3e-02 [ Info: 105 4.6e-01 3.5e-02 1.4e+00 8.5e-02 [ Info: 106 4.1e-01 4.9e-02 5.2e+00 2.8e-02 [ Info: 107 3.7e-01 4.1e-02 5.3e+00 7.7e-03 [ Info: 108 3.5e-01 2.0e-02 9.0e+00 4.9e-01 [ Info: 109 3.1e-01 3.9e-02 2.2e+01 1.5e-01 [ Info: 110 2.9e-01 2.4e-02 1.8e+01 5.7e-02 [ Info: 111 2.7e-01 2.1e-02 2.1e+01 1.0e-01 [ Info: 112 2.4e-01 3.4e-02 3.1e+01 2.7e-02 [ Info: 113 2.2e-01 1.6e-02 3.6e+01 1.9e-01 [ Info: 114 2.5e-01 2.9e-02 2.6e+01 2.9e+38 [ Info: 115 2.6e-01 9.6e-03 3.2e+01 1.6e+38 [ Info: 116 2.2e-01 3.8e-02 3.4e+01 2.5e-01 [ Info: 117 2.9e-01 6.9e-02 1.8e+01 7.8e+37 [ Info: 118 2.1e-01 8.4e-02 3.0e+01 3.1e-01 [ Info: 119 2.6e-01 5.2e-02 2.9e+01 2.1e+38 [ Info: 120 2.8e-01 2.3e-02 2.3e+01 7.3e+37 ε-relaxed chunk 1: V=0.649, constr=0.0015
1875-element Vector{Float64}:
0.5777483778167181
0.8041315810927907
0.8049768123539591
0.845521174881268
0.8533638103278612
0.8533638103278612
0.8533924849142122
0.8536382770911134
0.8337017642065413
0.8604940953473451
⋮
0.7654606320255918
0.8399562894641959
0.8817906829720652
0.23442108141669013
0.4129535712706737
0.7958733902582941
0.9246440934172928
0.9678842390734571
0.9588035113176117Results and verification
All three formulations converge to the same topology family — a diagonal strut feeding the load into a smoothly rounded re-entrant corner:
function report(label, x)
topo = filter(PseudoDensities(x)).x
srel = σ_rel(filter(PseudoDensities(x)))
smic = σ_micro(filter(PseudoDensities(x)))
solid = topo .> 0.5
println(
"$label: V = $(round(obj(x), digits=3)), max relaxed σ = $(round(maximum(srel), digits=3)), " *
"max microscopic σ (solid elements) = $(round(maximum(smic[solid]), digits=3)), " *
"σlim = $(round(σlim, digits=3))",
)
return topo
end
topo1 = report("p-norm ", x1)
topo2 = report("KS ", x2)
topo3 = report("ε-relaxed", x3)p-norm : V = 0.534, max relaxed σ = 1.174, max microscopic σ (solid elements) = 1.651, σlim = 1.047
KS : V = 0.59, max relaxed σ = 0.908, max microscopic σ (solid elements) = 1.105, σlim = 1.047
ε-relaxed: V = 0.649, max relaxed σ = 1.121, max microscopic σ (solid elements) = 1.405, σlim = 1.0471875-element Vector{Float64}:
0.715404853430725
0.7677760568164756
0.8134921131287172
0.8374696547119405
0.846494521821564
0.8494951399902215
0.8502778238500724
0.8511435168520584
0.8539552276331573
0.8561368024505925
⋮
0.7643367900003815
0.8031607522638377
0.7577529154143697
0.649439704802439
0.6434593488558055
0.7641764782220632
0.8883116905072358
0.9440101829267941
0.9553957482894964The microscopic stress is reported over solid elements only: with a density filter, the gray boundary layer is penalized material whose true stress exceeds the relaxed measure by \(\rho^{-q}\), a well-documented artifact of filtered SIMP stress evaluation (da Silva et al., 2019).
using CairoMakie
fig = visualize(problem; topology=topo1, default_exagg_scale=0.0)Precompiling packages... 6611.3 ms ✓ QuartoNotebookWorkerMakieExt (serial) 1 dependency successfully precompiled in 7 seconds Precompiling packages... 5470.3 ms ✓ QuartoNotebookWorkerCairoMakieExt (serial) 1 dependency successfully precompiled in 6 seconds
fig = visualize(problem; topology=topo2, default_exagg_scale=0.0)
fig = visualize(problem; topology=topo3, default_exagg_scale=0.0)
Discussion
- Relaxation is not optional. Without
stress_exponent(orepsilon_relaxed), the microscopic stress is finite in vanishing material and the stress-constrained problem has singular optima that MMA cannot reach — the optimization stalls or collapses. - Aggregation needs normalization. The raw p-norm
norm(σ, P)scales withN^(1/P)and mixes bulk and peak stress; constrain the normalizedN^(-1/P) * norm(σ, P)plus an adaptive factor, or the KS function. - The KS slack
log(N)/γis conservative by design (approach 2 uses more material); the p-norm with adaptivecis tighter but needs the normalization update. - Point loads and supports are singular stress sources; distribute them (
load_widthhere) or exclude their vicinity from the aggregation. - Deeper discussion of the theory and the failure modes of the naive formulation: see
STRESS_CONSTRAINED_TO.mdin the repository root.
References
- Le, Norato, Bruns, Ha & Tortorelli (2010). Stress-based topology optimization for continua. Struct. Multidiscip. Optim. 41:605–620.
- Duysinx & Bendsøe (1998). Topology optimization of continuum structures with local stress constraints. Int. J. Numer. Meth. Eng. 43:1453–1478.
- Cheng & Guo (1997). ε-relaxed approach in structural topology optimization. Struct. Optim. 13:258–266.
- Bruggi (2008). On an alternative approach to stress constraints relaxation in topology optimization. Struct. Multidiscip. Optim. 36:125–141.
- Duysinx & Sigmund (1998). New developments in handling stress constraints in optimal material distribution. 7th AIAA/USAF/NASA/ISSMo Symp.
- Yang & Chen (1996). Stress-based topology optimization. Struct. Optim. 12:98–105.
- Verbart, Langelaar & van Keulen (2017). A unified aggregation and relaxation approach for stress-constrained topology optimization. Struct. Multidiscip. Optim. 55:663–679.
- Stolpe & Svanberg (2001). On the trajectories of the epsilon-relaxation approach for stress-constrained truss topology optimization. Struct. Multidiscip. Optim. 21:140–151.
- da Silva, Beck & Sigmund (2019). Stress-constrained topology optimization considering uniform manufacturing uncertainties. Comput. Methods Appl. Mech. Eng. 344:512–537.
- Sved & Ginos (1968). Structural optimization under multiple loading. Int. J. Mech. Sci. 10:803–805.