fix and test multimaterial example

Project.toml

@@ -29,12 +29,14 @@ JuliaFormatter = "1"
 NonconvexIpopt = "0.4"
 Reexport = "1"
 StaticArraysCore = "1"
+TopOpt = "0.8"
 UnPack = "1"
 Zygote = "0.6"
 julia = "1"
 
 [extras]
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
+TopOpt = "53a1e1a5-51bb-58a9-8a02-02056cc81109"
 
 [targets]
-test = ["Test"]
+test = ["Test", "TopOpt"]

test/example.jl

@@ -26,7 +26,7 @@ function uncertainComp(x, logEs)
     Es = exp.(logEs)
     # interpolation of properties between materials
     interp = MaterialInterpolation(Es, penalty)
-    MultiMaterialVariables(x, nmats) |> interp |> filter |> comp
+    MultiMaterialVariables(x, nmats) |> tounit |> filter |> interp |> comp
     # return sum(x) + sum(Es)
 end
 # wrap original function in RandomFunction struct
@@ -37,19 +37,29 @@ x0 = fill(M, ncells * (length(logEs) - 1))
 # (Returns propability distribution of the objective for current point)
 d = rf(x0)
 # mass constraint
+constr_penalty = TopOpt.PowerPenalty(1.0)
+constr_interp = MaterialInterpolation(densities, constr_penalty)
 constr = x -> begin
-    ρs = PseudoDensities(MultiMaterialVariables(x, nmats))
-    return sum(element_densities(ρs, densities)) / ncells - 0.3 # unit element volume
+    ρs = constr_interp(MultiMaterialVariables(x, nmats))
+    return sum(ρs.x) / ncells - 0.3 # elements have unit volumes
 end
 function obj(x) # objective for TO problem
     dist = rf(x)
     mean(dist)[1] + 2 * sqrt(cov(dist)[1, 1])
 end
 obj(x0)
-Zygote.gradient(obj, x0)
-FiniteDifferences.grad(central_fdm(5, 1), obj, x0)[1]
+g1 = Zygote.gradient(obj, x0)[1]
+g2 = FiniteDifferences.grad(central_fdm(5, 1), obj, x0)[1]
+norm(g1 - g2)
+
+cg1 = Zygote.gradient(constr, x0)[1]
+cg2 = FiniteDifferences.grad(central_fdm(5, 1), constr, x0)[1]
+norm(cg1 - cg2)
 
 m = Model(obj) # create optimization model
 addvar!(m, zeros(length(x0)), ones(length(x0))) # setup optimization variables
 Nonconvex.add_ineq_constraint!(m, constr) # setup volume inequality constraint
+
 @time r = Nonconvex.optimize(m, TOBSAlg(), x0; options = TOBSOptions())
+# @time r = Nonconvex.optimize(m, IpoptAlg(), x0; options = IpoptOptions())
+# @time r = Nonconvex.optimize(m, MMA87(), x0; options = MMAOptions())

test/runtests.jl

@@ -21,3 +21,7 @@ using ReliabilityOptimization, Test, FiniteDifferences, Zygote
         @test norm(g1 - g2) < 1e-7
     end
 end
+
+@testset "Multi-material TopOpt example" begin
+    include("example.jl")
+end