using TopOpt, LinearAlgebra
using WGLMakie
WGLMakie.activate!(; resize_to=:parent)
using Bonito
if haskey(ENV, "QUARTO_PROJECT_DIR")
Bonito.Page(exportable=true, offline=true)
else
Bonito.browser_display()
end
display_app(app) = display(app)Fixed Element Regions with FixedElementProjectorFun
Description
Many practical designs reserve parts of the domain: a solid region where a load is applied or a bolt hole is drilled, a void region that must stay empty for access, or a keep-out zone near another component. FixedElementProjectorFun maps free design variables to a full density vector in which selected elements are pinned to solid (black, density 1) or void (white, density 0), so the optimizer only varies the remaining elements.
This tutorial demonstrates the workflow on a cantilever, fixing the support column and the load-application region as solid and carving a keep-out void, then running a SIMP optimization over the free elements.
Setup
WGLMakie.activate!(; resize_to=:parent) selects the browser renderer and fills the Quarto output column. Bonito.Page(exportable=true, offline=true) embeds the assets needed by visualize(...; static=true) in the Quarto output, so the visualization does not require a running Julia process.
Define the problem and the fixed regions
E, ν, f = 1.0, 0.3, 1.0
nels = (40, 20)
problem = PointLoadCantilever(nels, (1.0, 1.0), E, ν, f)
nx, ny = nels
nel = getncells(problem)
# Cells are ordered row-major (x fastest): index = (row-1)*nx + col.
cell_block(xr, yr) = [r * nx + c for r in (yr[1]-1):(yr[2]-1) for c in xr]
# Solid (black): the left support column and the right load column.
black_cells = cell_block(1:1, 1:ny)
# Void (white): a rectangular keep-out hole in the upper-right quadrant.
white_cells = cell_block(30:40, 13:20)
projector = get_fixed_element_projector(problem, black_cells, white_cells)
nfree = get_free_variable_count(projector)
println("fixed solid: $(count(projector.black)), fixed void: $(count(projector.white)), free: $nfree")The projector throws an ArgumentError if the same element is listed as both black and white, or if the free-variable vector has the wrong length — a fail-fast guard against silent mis-projection.
Objective and constraint over free variables
solver = FEASolver(DirectSolver, problem; xmin=0.001, penalty=PowerPenaltyFun(3.0))
comp = ComplianceFun(solver)
volfrac = VolumeFun(solver)
V = 0.4
obj = x_free -> comp(PseudoDensities(projector(x_free)))
constr = x_free -> volfrac(PseudoDensities(projector(x_free))) - VNote that the volume constraint counts the fixed solid elements too: the fixed regions already consume part of the volume budget, so the optimizer spreads the remainder over the free elements.
Check the initial layout
x0_free = fill(0.5, nfree)
ρ0 = projector(x0_free)
println("fixed black respected: $(all(ρ0[projector.black] .== 1.0))")
println("fixed white respected: $(all(ρ0[projector.white] .== 0.0))")
println("free elements at 0.5: $(all(ρ0[projector.free] .== 0.5))")fixed black respected: true
fixed white respected: true
free elements at 0.5: truefig0 = visualize(problem; static=true, topology=ρ0)
display_app(fig0)undeformed mesh
load arrows
support arrows
Run SIMP optimization over the free elements
model = Model(obj)
addvar!(model, zeros(nfree), ones(nfree))
add_ineq_constraint!(model, constr)
r = optimize(model, MMA87(), x0_free;
options=MMAOptions(; maxiter=150, tol=Tolerance(; kkt=1e-4, f=1e-4)))
x_opt = projector(r.minimizer)
println("final compliance = $(round(obj(r.minimizer), digits=3)), volume = $(round(constr(r.minimizer) + V, digits=3))")
println("fixed regions still respected: $(all(x_opt[projector.black] .== 1.0) && all(x_opt[projector.white] .== 0.0))")[ Info: iter obj Δobj violation kkt_residual [ Info: 0 3.2e+02 Inf 8.7e-02 3.4e+01 [ Info: 1 4.2e+02 9.5e+01 0.0e+00 1.6e+02 [ Info: 2 2.5e+02 1.7e+02 0.0e+00 4.0e+01 [ Info: 3 1.8e+02 6.5e+01 0.0e+00 1.4e+01 [ Info: 4 1.6e+02 2.7e+01 0.0e+00 5.5e+00 [ Info: 5 1.4e+02 1.5e+01 0.0e+00 3.4e+00 [ Info: 6 1.3e+02 1.2e+01 0.0e+00 3.1e+00 [ Info: 7 1.2e+02 1.2e+01 0.0e+00 2.9e+00 [ Info: 8 1.1e+02 1.1e+01 0.0e+00 2.6e+00 [ Info: 9 9.8e+01 8.0e+00 0.0e+00 1.9e+00 [ Info: 10 9.2e+01 5.5e+00 0.0e+00 1.2e+00 [ Info: 11 8.8e+01 4.2e+00 0.0e+00 8.7e-01 [ Info: 12 8.5e+01 3.1e+00 0.0e+00 7.3e-01 [ Info: 13 8.2e+01 2.5e+00 0.0e+00 5.5e-01 [ Info: 14 8.0e+01 1.9e+00 0.0e+00 3.9e-01 [ Info: 15 7.9e+01 1.1e+00 0.0e+00 2.3e-01 [ Info: 16 7.9e+01 6.6e-01 0.0e+00 1.5e-01 [ Info: 17 7.8e+01 2.9e-01 0.0e+00 1.4e-01 [ Info: 18 7.8e+01 1.6e-01 0.0e+00 1.4e-01 [ Info: 19 7.8e+01 8.1e-02 0.0e+00 1.7e-01 [ Info: 20 7.8e+01 3.3e-02 0.0e+00 1.3e-01 [ Info: 21 7.8e+01 1.4e-02 0.0e+00 8.1e-02 [ Info: 22 7.8e+01 8.0e-03 0.0e+00 4.5e-02 [ Info: 23 7.8e+01 4.3e-03 0.0e+00 2.3e-02 [ Info: 24 7.8e+01 3.4e-03 0.0e+00 3.2e-02 [ Info: 25 7.8e+01 1.3e-02 0.0e+00 9.2e-02 [ Info: 26 7.8e+01 1.5e-02 0.0e+00 1.1e-02 [ Info: 27 7.8e+01 8.5e-04 0.0e+00 7.7e-03 [ Info: 28 7.8e+01 1.2e-04 0.0e+00 4.8e-03 [ Info: 29 7.8e+01 5.1e-05 0.0e+00 3.1e-03 [ Info: 30 7.8e+01 3.0e-05 0.0e+00 2.0e-03 [ Info: 31 7.8e+01 1.6e-05 0.0e+00 1.2e-03 [ Info: 32 7.8e+01 8.2e-06 0.0e+00 6.9e-04 [ Info: 33 7.8e+01 3.4e-06 0.0e+00 3.0e-04 [ Info: 34 7.8e+01 1.1e-06 0.0e+00 1.6e-04 [ Info: 35 7.8e+01 8.9e-08 0.0e+00 2.4e-04 [ Info: 36 7.8e+01 5.1e-07 0.0e+00 4.3e-04 [ Info: 37 7.8e+01 4.0e-07 0.0e+00 4.7e-04 [ Info: 38 7.8e+01 5.9e-07 0.0e+00 2.2e-04 [ Info: 39 7.8e+01 4.8e-07 0.0e+00 1.4e-05 final compliance = 78.078, volume = 0.4 fixed regions still respected: true
Visualize the result
fig = visualize(problem; static=true, topology=x_opt)
display_app(fig)undeformed mesh
load arrows
support arrows
The design loads the fixed solid columns at support and tip and routes material around the keep-out hole, leaving it empty. Gradients flow only through the free elements: the projector’s rrule returns zero cotangent on black/white elements, so the fixed regions never change.