Update NOMAD interface to NOMAD 4.5.1 (#83)

Author: salomonl

Project.toml

@@ -1,7 +1,7 @@
 name = "NOMAD"
 uuid = "02130f1c-4665-5b79-af82-ff1385104aa0"
 authors = ["Alexis Montoison (@amontoison), Pierrick Pascal (@ppascal97), Ludovic Salomon (@salomonl)"]
-version = "2.4.2"
+version = "3.0.0"
 
 [deps]
 Libdl = "8f399da3-3557-5675-b5ff-fb832c97cbdb"
@@ -12,6 +12,6 @@ Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 Tulip = "6dd1b50a-3aae-11e9-10b5-ef983d2400fa"
 
 [compat]
-NOMAD_jll = "= 4.3.1"
+NOMAD_jll = "= 4.5.1"
 Tulip = "0.9"
 julia = "1.10"

docs/Project.toml

@@ -1,6 +1,8 @@
 [deps]
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
+Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 
 [compat]
-Documenter = "0.27"
-julia = "1.6"
+Documenter = "1.0"
+julia = "1.10"
+Plots = "1.41.1"

docs/make.jl

@@ -10,8 +10,6 @@ makedocs(
   pages = ["Home" => "index.md",
            "Parameters" => "nomadProblem.md",
            "Run Optimization" => "run_nomad.md",
-           #"Results" => "nomadResults.md",
-           #"Surrogates" => "surrogates.md",
            "Tutorial" => "tutorial.md"
           ]
 )

docs/src/nomadParameters.md

@@ -69,7 +69,9 @@ p.options.display_stats = ["BBE", "EVAL", "SOL", "OBJ", "CONS_H"] # some display
 
 # Solution
 result = solve(p, [0.0;2.0])
-println("Solution: ", result.x_best_feas)
+println("Optimization status: ", result.status)
+println("Solution: ", result.x_sol)
+println("is feasible: ", result.feasible)
 ```
 
 ## Linear equality constraints
@@ -131,9 +133,10 @@ x0 = [0.57186958424864897665429452899843;
 
 # Solution
 result = solve(p, x0)
-println("Solution: ", result.x_best_feas)
-println("Satisfy Ax = b: ", A * result.x_best_feas ≈ b)
-println("And inside bound constraints: ", all(-10.0 .<= result.x_best_feas .<= 10.0))
+println("status: ", result.status)
+println("Solution: ", result.x_sol)
+println("Satisfy Ax = b: ", A * result.x_sol ≈ b)
+println("And inside bound constraints: ", all(-10.0 .<= result.x_sol .<= 10.0))
 ```
 
 The reader can take a look at the [`test`](https://github.com/bbopt/NOMAD.jl/tree/master/test) folder for more complex examples.
@@ -203,3 +206,151 @@ For more details about the parameters used in this section, we refer the reader
 
 [C. Audet, A. Ianni, S. Le Digabel and C. Tribes, Reducing the number of function evaluations in mesh
 adaptive direct search algorithms, *SIAM Journal on Optimization*, 24(2), 621-642, 2014.](https://doi.org/10.1137/120895056)
+
+## Multiobjective optimization
+
+`NOMAD.jl` can solve multiobjective optimization problems of the form:
+```math
+\begin{array}{rl}
+  (MBB) \ \ \
+  \displaystyle \min_{x \in \mathbb{R}^n} & f(x) = \left(f_1(x), f_2(x), \ldots, f_m(x)\right)^\top\\
+  s.t.
+  & c(x) \leq 0 \\
+  & lb \leq x \leq ub \\
+\end{array}
+```
+where $m \geq 2$ is the number of objectives.
+
+Unlike single-objective optimization, the set of solutions to such a problem is generally not a singleton.
+It is composed of several solutions whose objective values represent the best trade-offs that the
+algorithm may achieve at the end of the optimization.
+
+`NOMAD.jl` switches to multiobjective mode when a `NomadProblem` is given more than one objective `OBJ` as
+parameter. However, `NOMAD.jl` only supports up to `5` objectives. Keep in mind that the larger the number of
+objectives, the more difficult the problem is to solve.
+
+Finally, the multiobjective mode of `NOMAD.jl` does not support certain options available in
+single-objective optimization. Specifically, in multiobjective optimization, `NOMAD.jl` deactivates
+all search strategies except the speculative search, the `NM` search and the `Quadratic model`
+search. The default `eval_sort_type` and `direction_type` options are switched to `"DIR_LAST_SUCCESS"`
+and `"ORTHO N+1 QUAD"`, respectively.
+
+We consider the following example, the optimal design of a welded beam.
+
+```math
+\begin{array}{rl}
+  (WB) \ \ \
+  \displaystyle \min_{x = (h, l, t, b) \in \mathbb{R}^4} & \left(1.10471 h^2 l + 0.04811 t b (14.0 + l), \dfrac{2.1952}{t^3 b}\right)^\top \\
+  s.t.
+  & c_1(x) = \tau(x) - 13600 \leq 0  \\
+  & c_2(x) = \sigma(x) - 30000 \leq 0\\
+  & c_3(x) = h - b \leq 0\\
+  & c_4(x) = 6000 - P_c(x) \leq 0\\
+  & h, b \in [0.125, 5], l, t \in [0.1, 10]\\
+\end{array}
+```
+where $h, l, t, b$ are the four design parameters to optimize.
+
+The following equations describe the physics constraints of the system (stress and buckling terms):
+
+```math
+\begin{array}{rl}
+  \tau(x) & = \sqrt{(\tau')^2 + (\tau'')^2 + \dfrac{l \tau' \tau''}{\sqrt{0.25 (l^2 + (h + t)^2))}}}, \\ 
+  \tau' & = \dfrac{6000}{\sqrt{2} h l},  \\
+  \tau'' & = \dfrac{6000 (14 + 0.5 l) \sqrt{0.25 (l^2 + (h + t)^2))}}{2 \left(0.707 h l (l^2/12 + 0.25 (h + t)^2)\right)},\\
+  \sigma(x) & = \dfrac{504000}{t^2 b},\\
+  P_c(x) & = 64746.022 (1 - 0.0282346 t) t b^3.
+\end{array}
+```
+
+$f_1$ corresponds to the construction costs and $f_2$ the end deflection (which corresponds to rigidity) of the beam. Both must be minimized.
+
+```@example welded_beam_design
+using NOMAD
+
+# blackbox
+function welded_beam(x)
+    # Variables
+    h = x[1]
+    l = x[2]
+    t = x[3]
+    b = x[4]
+
+    # Objectives
+    f1 = 1.10471 * h^2 * l + 0.04811 * t * b * (14.0 + l)
+    f2 = 2.1952 / (t^3 * b)
+
+    # Physics equations
+    tau_p = 6000.0 / (sqrt(2) * h * l)
+    tau_pp = 6000 * (14 + 0.5 * l) * sqrt(0.25 * (l^2 + (h + t)^2))
+    tau_pp /= (2 * (0.707 * h * l * (l^2 / 12.0 + 0.25 * (h + t)^2)))
+    tau = sqrt(tau_p^2 + tau_pp^2 + l * tau_p * tau_pp / sqrt(0.25 * (l^2 + (h + t)^2)))
+    sigma = 504000 / (t^2 * b)
+    Pc = 64746.022 * (1 - 0.0282346 * t) * t * b^3
+
+    # Constraints
+    c1 = tau - 13600
+    c2 = sigma - 30000
+    c3 = h - b
+    c4 = 6000 - Pc
+
+    bb_outputs = [f1, f2, c1, c2, c3, c4]
+    success = true
+    count_eval = true
+    return (success, count_eval, bb_outputs)
+end
+
+# Define the problem
+lb = [0.125, 0.1, 0.1, 0.125]
+ub = [5.0, 10.0, 10.0, 5.0]
+pb = NomadProblem(4, 6, ["OBJ", "OBJ", "PB", "PB", "PB", "PB"],
+                  welded_beam,
+                  lower_bound=[0.125, 0.1, 0.1, 0.125],
+                  upper_bound=[5.0, 10.0, 10.0, 5.0])
+
+# Set some options
+pb.options.display_degree = 2
+pb.options.max_bb_eval = 1500
+
+# As for the single-objective case, you could deactivate this
+# option to go faster, but the performance may be worse
+# pb.options.quad_model_search = false
+
+# Start from a middle box point
+result = solve(pb, (lb + ub) / 2)
+
+# The solution set is returned as a matrix of dimensions n x nb_solutions,
+# where n is the dimension of the problem.
+println("Optimization status: ", result.status)
+println("The algorithm has found ", size(result.x_sol, 2), " solutions:")
+for (ind, (x, bbo)) in enumerate(zip(eachcol(result.x_sol), eachcol(result.bbo_sol)))
+  println("Solution ", ind, ": ", x,
+          "; f(x) = ", bbo[1:2],
+          "; c(x) = ", bbo[3:end])
+end
+println("They are ", result.feasible ? "feasible" : "infeasible")
+
+# In multiobjective optimization, it can be interesting to see the
+# set of trade-offs in the objective space
+using Plots
+fig = scatter(result.bbo_sol[1, :], result.bbo_sol[2, :],
+              xlabel="Cost", ylabel="Deflection",
+              title="Pareto front approximation")
+```
+
+For more information about the multiobjective algorithm (DMulti-MADS), please refer to the following articles:
+
+[J. Bigeon, S. Le Digabel, & L. Salomon, DMulti-MADS: Mesh adaptive direct multisearch for bound-constrained
+ blackbox multiobjective optimization, *Computational Optimization and Applications*, 79(2), 301-338, 2021.](https://doi.org/10.1007/s10589-021-00272-9)
+
+[J. Bigeon, S. Le Digabel, & L. Salomon, Handling of constraints in multiobjective blackbox optimization,
+ *Computational Optimization and Applications*, 89(1), 69-113, 2024.](https://doi.org/10.1007/s10589-024-00588-2)
+
+[S. Le Digabel, A. Lesage-Landry, L. Salomon, & C. Tribes (2025),
+ Efficient search strategies for constrained multiobjective blackbox optimization,
+ arXiv preprint arXiv:2504.02986., 2025.](https://doi.org/10.48550/arXiv.2504.02986)
+
+The problem is taken from:
+
+[K. Deb, Evolutionary algorithms for multi-criterion optimization in engineering design,
+ *Evolutionary algorithms in engineering and computer science*, 2, 135-161, (1999)]()