04 Optimization

Discrete choices, nonlinear physics, several objectives.

Molecular design mixes integer decisions (which core, which anion, how many groups) with continuous variables (compositions, temperatures, flows) and nonlinear, often nonconvex, property models. This chapter covers the methods for such problems, deterministic and stochastic, single- and multi-objective, with and without uncertainty.

Computed for this page The Branin test function on a log scale (low values bright) with two Bayesian-optimization runs computed here: five random starts (gray), then a Gaussian-process surrogate and expected improvement on a 40 × 40 grid (teal; Eqs. 4.4, 8.2–8.3). White ring: best point so far. Crosses: the three global minima, f* = 0.398.

You will see

  • MINLP and why nonconvexity matters
  • Branch-and-bound, outer approximation, generalized Benders decomposition
  • Genetic algorithms and Pareto ranking (NSGA-II)
  • Bayesian optimization and expected improvement
  • Chance constraints and robust design

Labs in this chapter

4.1 Problem classes

minx, y f(x, y)   s.t.   g(x, y) ≤ 0,   h(x, y) = 0,   x ∈ X ⊂ ℝn,   y ∈ {0,1}m
(4.1)

With only continuous variables and linear functions, Eq. 4.1 is an LP; with nonlinear functions, an NLP; with integers and linear functions, an MILP; with both integers and nonlinear functions, an MINLP. Convexity matters: in a convex problem every local optimum is global. Activity-coefficient models make molecular design nonconvex, so local solvers can stop at poor designs, and certifying a global optimum requires bounding techniques.

4.2 Branch-and-bound

Branch-and-bound relaxes the integrality of y, solves the relaxation for a lower bound (for minimization), and branches on a fractional variable. A node is pruned when its bound is worse than the best feasible solution found. For nonconvex problems, spatial branch-and-bound also branches on continuous variables and builds convex underestimators; branch-and-reduce adds range reduction, the approach behind the BARON solver (Tawarmalani & Sahinidis, 2005); range reduction is used by most deterministic global solvers. The gap between the best bound and the best solution certifies how far from optimal the answer can be.

4.3 Decomposition: outer approximation and generalized Benders

Outer approximation (Duran & Grossmann, 1986) alternates two problems. Fixing yk gives an NLP whose solution xk is feasible and yields an upper bound. Linearizing at that point gives a mixed-integer master problem whose solution is a lower bound for convex MINLPs:

min η   s.t.   f(xk, yk) + ∇fT[x − xk; y − yk] ≤ η,   g(xk, yk) + ∇gT[x − xk; y − yk] ≤ 0,   ∀k
(4.2)

Generalized Benders decomposition (Geoffrion, 1972) projects onto the integer space and builds the master from Lagrangian cuts instead of full linearizations: fewer constraints, usually weaker bounds. The decomposition used in the 2013 study is a different, problem-specific kind: constraints are applied in stages, structural first, then pure-component, then mixture, so most candidates are eliminated before the expensive models run.

4.4 Genetic algorithms

A genetic algorithm evolves a population of encoded candidates: selection favors fitter parents, crossover recombines them, mutation adds random change, and elitism carries the best forward unchanged. It needs no derivatives and tolerates discontinuous or noisy objectives, but it certifies nothing about optimality. The 2013 study’s implementation (MATLAB, most parameters at default values) used a population of 20, a crossover fraction of 0.8, a mutation probability of 0.2 and two elite candidates per generation, with fitness equal to the objective.

Lab 4.1 Computed for this page

The 2013 genetic algorithm, reimplemented on the Chapter 2 design space

Objective: maximize conductivity subject to η < 65 cP (the paper’s melting-point constraint, Eqs. 35–36, is not implemented here). Population (20), crossover fraction (0.8), mutation probability (0.2) and elite count (2) as published; selection, encoding and the infeasibility penalty are this page’s choices. The dashed line is the optimum found by enumerating all candidates.

Chromosome: cation core, two chain lengths and anion. Crossover is uniform over genes; mutation resamples a gene. The chart opens on one completed 30-generation run from a random start; New random start returns to generation 0. Results vary with the start.

4.5 Multiple objectives and the Pareto front

Design objectives conflict: capacity against viscosity, performance against toxicity, cost against life-cycle impact. A design x dominates x′ if it is no worse in every objective and strictly better in at least one; the non-dominated set is the Pareto front.

x ≺ x′  ⇔  fk(x) ≤ fk(x′) ∀k   and   ∃k: fk(x) < fk(x′)
(4.3)

Weighted sums can miss parts of a nonconvex front; the ε-constraint method (optimize one objective, bound the others) traces it fully. NSGA-II (Deb et al., 2002) ranks a population by non-domination and preserves spread with a crowding distance. The Pareto chart on the Science page is a multi-objective replot of the founder’s 2016 thermal-storage results.

Figure 4.1 Computed for this page Published parameters

When two objectives agree, the front collapses

Every structure of the Chapter 2 design space with conductivity parameters, scored for viscosity (minimize) and conductivity (maximize). A point is Pareto-optimal if no other structure is at least as good on both objectives and better on one. Here a single liquid dominates all the others, because in these two models lower viscosity goes with higher conductivity.

    Group-contribution predictions at 298.15 K (2013 paper, Tables 4–5), not measurements. Real conflicts come from objectives with different physics, such as performance against toxicity, cost or life-cycle impact; the front then becomes a curve or surface, which NSGA-II or the ε-constraint method trace.

    4.6 Bayesian optimization

    When each evaluation is expensive, such as a synthesis and measurement, the number of evaluations matters more than computing time. Bayesian optimization fits a Gaussian-process surrogate (Chapter 8) and chooses the next point by an acquisition function. Expected improvement (Mockus, Tiešis & Žilinskas, 1978; used in the EGO algorithm of Jones, Schonlau & Welch, 1998), for minimization:

    EI(x) = (fmin − μ(x)) Φ(Z) + σ(x) φ(Z),    Z = fmin − μ(x)σ(x)
    (4.4)

    μ, σ: surrogate mean and standard deviation; Φ, φ: standard normal CDF and density; fmin: best value so far. The first term rewards predicted improvement; the second rewards uncertainty. With constraints, EI is multiplied by the predicted probability that each constraint holds.

    Lab 4.2 Computed for this page

    Bayesian optimization on the Branin function

    The Branin function is a standard test problem with known global minimum f* = 0.397887 at three points. A Gaussian process with a squared-exponential kernel (length scale chosen from nine candidates by marginal likelihood at each step) is refitted after every evaluation; the next point maximizes the chosen acquisition function on a 60 × 60 grid.

    Acquisition functions, minimizationEI = (fmin − μ)Φ(Z) + σφ(Z), Z = (fmin − μ)/σ · PI = Φ((fmin − ξ − μ)/σ), ξ = 0.01 × standard deviation of observed f · LCB: choose the minimum of μ − κσ
    f(x1, x2) = (x2 − 5.1x12/(4π2) + 5x1/π − 6)2 + 10(1 − 1/(8π))cos x1 + 10, on x1 ∈ [−5, 10], x2 ∈ [0, 15]; minima at (−π, 12.275), (π, 2.275) and (3π, 2.475). Source: Surjanovic & Bingham, Virtual Library of Simulation Experiments. Five random starting points, then the chosen acquisition. EI: Mockus, Tiešis & Žilinskas (1978), as used by Jones, Schonlau & Welch (1998); PI: Kushner (1964); confidence-bound selection: Srinivas, Krause, Kakade & Seeger (2010), used here with a constant κ rather than their increasing schedule.

    4.7 Optimization under uncertainty

    Every property in a design problem is a prediction with error. Three standard formulations carry that error into the decision:

    stochastic:   minx Eξ[f(x, ξ)];    chance-constrained:   P[g(x, ξ) ≤ 0] ≥ 1 − ε;    robust:   minx maxξ∈U f(x, ξ)
    (4.5)

    In molecular design, a chance constraint reads naturally: the predicted chronic EC50 must exceed a threshold with 95% probability, not merely on average. That requires calibrated predictive distributions, which is why uncertainty quantification (Chapter 8) is part of the optimization layer, not an add-on. Chapter 5 writes the full design problem this way. Planned

    References

    1. Duran, M. A.; Grossmann, I. E. An outer-approximation algorithm for a class of mixed-integer nonlinear programs. Math. Program. 1986, 36, 307–339. doi:10.1007/BF02592064
    2. Geoffrion, A. M. Generalized Benders decomposition. J. Optim. Theory Appl. 1972, 10, 237–260. doi:10.1007/BF00934810
    3. Tawarmalani, M.; Sahinidis, N. V. A polyhedral branch-and-cut approach to global optimization. Math. Program. 2005, 103, 225–249. doi:10.1007/s10107-005-0581-8
    4. Deb, K.; Pratap, A.; Agarwal, S.; Meyarivan, T. A fast and elitist multiobjective genetic algorithm: NSGA-II. IEEE Trans. Evol. Comput. 2002, 6, 182–197. doi:10.1109/4235.996017
    5. Jones, D. R.; Schonlau, M.; Welch, W. J. Efficient global optimization of expensive black-box functions. J. Global Optim. 1998, 13, 455–492. doi:10.1023/A:1008306431147
    6. Mockus, J.; Tiešis, V.; Žilinskas, A. The application of Bayesian methods for seeking the extremum. In Towards Global Optimization 2; Dixon, L. C. W., Szegő, G. P., Eds.; North-Holland, 1978; pp 117–129.
    7. Kushner, H. J. A new method of locating the maximum point of an arbitrary multipeak curve in the presence of noise. J. Basic Eng. 1964, 86, 97–106. doi:10.1115/1.3653121
    8. Srinivas, N.; Krause, A.; Kakade, S. M.; Seeger, M. Gaussian process optimization in the bandit setting: no regret and experimental design. Proc. 27th Int. Conf. Mach. Learn. (ICML), 2010. arXiv:0912.3995
    9. Surjanovic, S.; Bingham, D. Virtual Library of Simulation Experiments: Branin function. Simon Fraser University. sfu.ca
    10. Karunanithi, A. T.; Mehrkesh, A. AIChE J. 2013, 59, 4627–4640 (genetic-algorithm settings). doi:10.1002/aic.14228
    11. Gardner, J. R.; Kusner, M. J.; Xu, Z.; Weinberger, K. Q.; Cunningham, J. P. Bayesian optimization with inequality constraints. Proc. 31st Int. Conf. Mach. Learn., PMLR 32, 2014, 937–945. proceedings.mlr.press

    Privacy

    Terms