Skip to content

Commit

Permalink
formatting
Browse files Browse the repository at this point in the history
  • Loading branch information
mohamed82008 committed Nov 15, 2023
1 parent a8ecbc6 commit 842becd
Show file tree
Hide file tree
Showing 2 changed files with 75 additions and 56 deletions.
51 changes: 35 additions & 16 deletions src/NonconvexTOBS.jl
Original file line number Diff line number Diff line change
Expand Up @@ -15,17 +15,26 @@ struct TOBSOptions
nt::NamedTuple # one of the fields of the nt is milp_options
end

function TOBSOptions(
;movelimit::Real = 0.1, # move limit parameter
function TOBSOptions(;
movelimit::Real = 0.1, # move limit parameter
pastN::Int = 20, # number of past iterations for moving average calculation
convParam::Real = 0.001, # convergence parameter (upper bound of error)
constrRelax::Real = 0.1, # constraint relaxation parameter
timeLimit::Real = 1.0,
optimizer = Cbc.Optimizer,
maxiter::Int = 200,
timeStable::Bool = true
timeStable::Bool = true,
)
return TOBSOptions((; movelimit, pastN, convParam, constrRelax, timeLimit, optimizer, maxiter, timeStable))
return TOBSOptions((;
movelimit,
pastN,
convParam,
constrRelax,
timeLimit,
optimizer,
maxiter,
timeStable,
))
end

@params mutable struct TOBSWorkspace <: Workspace
Expand All @@ -34,20 +43,29 @@ end
options::TOBSOptions
end
function TOBSWorkspace(
model::VecModel, x0::AbstractVector = NonconvexCore.getinit(model);
options = TOBSOptions(), kwargs...,
model::VecModel,
x0::AbstractVector = NonconvexCore.getinit(model);
options = TOBSOptions(),
kwargs...,
)
return TOBSWorkspace(model, copy(x0), options)
end
@params struct TOBSResult <: AbstractResult
minimizer
minimum
error
minimizer::Any
minimum::Any
error::Any
end

function optimize!(workspace::TOBSWorkspace)
@unpack model, x0, options = workspace
@unpack movelimit, pastN, convParam, constrRelax, timeLimit, optimizer, maxiter, timeStable = options.nt
@unpack movelimit,
pastN,
convParam,
constrRelax,
timeLimit,
optimizer,
maxiter,
timeStable = options.nt
milp_solver = JuMP.optimizer_with_attributes(optimizer)
numVars = length(NonconvexCore.getinit(model))
count = 1 # iteration counter
Expand All @@ -70,7 +88,8 @@ function optimize!(workspace::TOBSWorkspace)
JuMP.set_optimizer_attribute(m, "seconds", timeLimit)
if !skip_step || count == 1
if count > 1
currentConstr, jacConstr = NonconvexCore.value_jacobian(model.ineq_constraints, x)
currentConstr, jacConstr =
NonconvexCore.value_jacobian(model.ineq_constraints, x)
end
violation = norm(currentConstr)
if (violation <= best_sol[4] - 1e-8 || violation < 1e-8 && objval < best_sol[2])
Expand All @@ -88,7 +107,7 @@ function optimize!(workspace::TOBSWorkspace)
JuMP.@constraint(m, deltaX .<= absdeltaX)
JuMP.@constraint(m, .-deltaX .<= absdeltaX)
# Constrain amount of change per iteration
JuMP.@constraint(m, sum(absdeltaX) <= movelimit*numVars)
JuMP.@constraint(m, sum(absdeltaX) <= movelimit * numVars)
# Constraint relaxation
Δ = map(currentConstr) do c
abs(c) < constrRelax ? -c : -constrRelax * c
Expand Down Expand Up @@ -122,12 +141,12 @@ function optimize!(workspace::TOBSWorkspace)
objval, objgrad = NonconvexCore.value_gradient(getobjective(model), x)
# Apply "time stabilization"
if timeStable
objgrad = (objgrad + pastGrad)/2
objgrad = (objgrad + pastGrad) / 2
pastGrad = copy(objgrad)
end
objHist[end] = objval
if count > pastN
er = abs(sum([objHist[i] - objHist[i - 1] for i in 2:pastN])) / sum(objHist)
er = abs(sum([objHist[i] - objHist[i-1] for i = 2:pastN])) / sum(objHist)
end
@info "iter = $count, obj = $(round.(objHist[end]; digits=3)), constr_vio_norm = $(round(norm(currentConstr), digits=3)), er = $(round(er, digits=3))"
end
Expand All @@ -136,8 +155,8 @@ function optimize!(workspace::TOBSWorkspace)
return TOBSResult(best_sol[1], best_sol[2], er)
end

function Workspace(model::VecModel, optimizer::TOBSAlg, x0::AbstractVector; kwargs...,)
function Workspace(model::VecModel, optimizer::TOBSAlg, x0::AbstractVector; kwargs...)
return TOBSWorkspace(model, x0; kwargs...)
end

end
end
80 changes: 40 additions & 40 deletions test/runtests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2,44 +2,44 @@ using NonconvexTOBS, TopOpt, Test

@testset "Example" begin

# Recreating cantilever problem from original paper

E = 1.0 # Young’s modulus
v = 0.3 # Poisson’s ratio
f = 1.0 # downward force
rmin = 6.0 # filter radius
xmin = 0.001 # minimum density
V = 0.5 # maximum volume fraction
p = 3.0 # topological optimization penalty

# Define FEA problem
problem_size = (160, 100) # size of rectangular mesh
x0 = fill(1.0, prod(problem_size)) # initial design
problem = PointLoadCantilever(Val{:Linear}, problem_size, (1.0, 1.0), E, v, f)

# define FEA solver and auxiliary functions
solver = FEASolver(Direct, problem; xmin=xmin)
cheqfilter = DensityFilter(solver; rmin=rmin) # filter function
comp = TopOpt.Compliance(problem, solver) # compliance function

obj(x) = comp(cheqfilter(x)) # compliance objective
constr(x) = sum(cheqfilter(x)) / length(x) - V # volume fraction constraint

# Optimization setup
using NonconvexUtils
tobj = TraceFunction(obj, on_grad = true)
m = Model(tobj) # create optimization model
addvar!(m, zeros(length(x0)), ones(length(x0))) # setup optimization variables
Nonconvex.add_ineq_constraint!(m, constr) # setup volume inequality constraint
options = TOBSOptions() # optimization options with default values
TopOpt.setpenalty!(solver, p)

# Perform TOBS optimization
@time r = Nonconvex.optimize(m, TOBSAlg(), x0; options=options)

# Results
@show obj(r.minimizer)
@show constr(r.minimizer)
topology = r.minimizer
# Recreating cantilever problem from original paper

E = 1.0 # Young’s modulus
v = 0.3 # Poisson’s ratio
f = 1.0 # downward force
rmin = 6.0 # filter radius
xmin = 0.001 # minimum density
V = 0.5 # maximum volume fraction
p = 3.0 # topological optimization penalty

# Define FEA problem
problem_size = (160, 100) # size of rectangular mesh
x0 = fill(1.0, prod(problem_size)) # initial design
problem = PointLoadCantilever(Val{:Linear}, problem_size, (1.0, 1.0), E, v, f)

# define FEA solver and auxiliary functions
solver = FEASolver(Direct, problem; xmin = xmin)
cheqfilter = DensityFilter(solver; rmin = rmin) # filter function
comp = TopOpt.Compliance(problem, solver) # compliance function

obj(x) = comp(cheqfilter(x)) # compliance objective
constr(x) = sum(cheqfilter(x)) / length(x) - V # volume fraction constraint

# Optimization setup
using NonconvexUtils
tobj = TraceFunction(obj, on_grad = true)
m = Model(tobj) # create optimization model
addvar!(m, zeros(length(x0)), ones(length(x0))) # setup optimization variables
Nonconvex.add_ineq_constraint!(m, constr) # setup volume inequality constraint
options = TOBSOptions() # optimization options with default values
TopOpt.setpenalty!(solver, p)

# Perform TOBS optimization
@time r = Nonconvex.optimize(m, TOBSAlg(), x0; options = options)

# Results
@show obj(r.minimizer)
@show constr(r.minimizer)
topology = r.minimizer

end

0 comments on commit 842becd

Please sign in to comment.