---

name: "pymoo" description: "Python framework for single- and multi-objective optimization with evolutionary algorithms. Define vectorized objectives and constraints; solve with NSGA-II, NSGA-III, MOEA/D, GAs, or differential evolution. Analyze Pareto fronts, visualize trade-offs, customize operators and callbacks. For engineering design, hyperparameter search, and conflicting objectives. Alternatives: scipy.optimize (single-objective, gradient), platypus, jMetalPy (Java)." license: "Apache-2.0"

pymoo

Overview

pymoo provides a unified API for multi-objective optimization via population-based evolutionary algorithms. Users define a problem by subclassing Problem or ElementwiseProblem, specifying objectives (n_obj), decision variables (n_var), and optional constraints (n_ieq_constr). Algorithms like NSGA-II and NSGA-III return a Result object containing the Pareto-optimal population, objective values, and decision variable values. pymoo separates problem definition, algorithm configuration, operator selection, and analysis — each component is independently replaceable.

When to Use

• Optimizing a design with two or more conflicting objectives (e.g., minimizing cost while maximizing performance)
• Running evolutionary algorithms (GA, DE, PSO) as black-box optimizers when gradients are unavailable
• Performing multi-objective hyperparameter search for ML models where accuracy and inference time trade off
• Computing Pareto fronts for portfolio optimization or multi-criteria decision analysis
• Customizing crossover/mutation operators for domain-specific solution encodings (binary, permutation, real-valued)
• Benchmarking optimization algorithms on standard test problems (ZDT, DTLZ, CTP)
• Use scipy.optimize instead for single-objective, gradient-available, smooth optimization

Prerequisites

• Python packages: pymoo, numpy, matplotlib
• Data requirements: objective function(s) and optional constraint functions; variable bounds
• Environment: CPU sufficient for most problems; GPU not used by pymoo core

pip install pymoo numpy matplotlib

Quick Start

import numpy as np
from pymoo.core.problem import Problem
from pymoo.algorithms.moo.nsga2 import NSGA2
from pymoo.optimize import minimize

class SimpleBiObjective(Problem):
    def __init__(self):
        super().__init__(n_var=2, n_obj=2, xl=np.array([-2, -2]), xu=np.array([2, 2]))

    def _evaluate(self, X, out, *args, **kwargs):
        f1 = X[:, 0] ** 2 + X[:, 1] ** 2
        f2 = (X[:, 0] - 1) ** 2 + X[:, 1] ** 2
        out["F"] = np.column_stack([f1, f2])

algorithm = NSGA2(pop_size=100)
res = minimize(SimpleBiObjective(), algorithm, ("n_gen", 200), seed=1, verbose=False)
print(f"Pareto front size: {len(res.F)}")
print(f"Objective range: F1=[{res.F[:,0].min():.3f}, {res.F[:,0].max():.3f}]")

Core API

Module 1: Problem Definition

Define optimization problems via subclassing. Use Problem for vectorized evaluation (faster), ElementwiseProblem for scalar evaluation (simpler to write).

import numpy as np
from pymoo.core.problem import Problem, ElementwiseProblem

# Vectorized problem (preferred for performance)
class ZDT1(Problem):
    """ZDT1 benchmark: 30 variables, 2 objectives, known Pareto front."""
    def __init__(self):
        super().__init__(n_var=30, n_obj=2, xl=0.0, xu=1.0)

    def _evaluate(self, X, out, *args, **kwargs):
        f1 = X[:, 0]
        g = 1 + 9 * X[:, 1:].mean(axis=1)
        f2 = g * (1 - np.sqrt(f1 / g))
        out["F"] = np.column_stack([f1, f2])

# Elementwise problem with inequality constraints
class ConstrainedProblem(ElementwiseProblem):
    def __init__(self):
        super().__init__(n_var=2, n_obj=1, n_ieq_constr=2,
                         xl=np.array([-5, -5]), xu=np.array([5, 5]))

    def _evaluate(self, x, out, *args, **kwargs):
        out["F"] = (x[0] - 1) ** 2 + (x[1] - 2) ** 2  # objective
        out["G"] = np.array([
            x[0] + x[1] - 2,   # g1 <= 0
            x[0] ** 2 - x[1],  # g2 <= 0
        ])

print(f"ZDT1: {ZDT1().n_var} vars, {ZDT1().n_obj} objectives")

# Mixed-variable problem: some integer, some real
from pymoo.core.variable import Real, Integer, Choice

class MixedProblem(ElementwiseProblem):
    def __init__(self):
        vars = {
            "x": Real(bounds=(-2, 2)),
            "n": Integer(bounds=(1, 10)),
        }
        super().__init__(vars=vars, n_obj=1)

    def _evaluate(self, X, out, *args, **kwargs):
        x, n = X["x"], X["n"]
        out["F"] = (x - n) ** 2

Module 2: Algorithm Selection

pymoo provides 20+ algorithms. Key choices by problem type:

from pymoo.algorithms.moo.nsga2 import NSGA2
from pymoo.algorithms.moo.nsga3 import NSGA3
from pymoo.algorithms.moo.moead import MOEAD
from pymoo.algorithms.soo.nonconvex.ga import GA
from pymoo.algorithms.soo.nonconvex.de import DE
from pymoo.util.ref_dirs import get_reference_directions

# NSGA-II: best for 2-3 objectives, most widely used
nsga2 = NSGA2(pop_size=100)

# NSGA-III: designed for 3+ objectives; needs reference directions
ref_dirs = get_reference_directions("das-dennis", 3, n_partitions=12)  # ~91 dirs
nsga3 = NSGA3(pop_size=len(ref_dirs), ref_dirs=ref_dirs)

# MOEA/D: decomposition-based, good for many objectives
moead = MOEAD(ref_dirs=ref_dirs, n_neighbors=15, prob_neighbor_mating=0.7)

# GA: single-objective genetic algorithm
ga = GA(pop_size=100)

# DE: Differential Evolution, good for continuous problems
de = DE(pop_size=100, variant="DE/rand/1/bin", CR=0.9, F=0.8)

print("Algorithms initialized")

Module 3: Operators (Crossover & Mutation)

Operators define how solutions evolve. Replace defaults to match variable type.

from pymoo.operators.crossover.sbx import SBX
from pymoo.operators.mutation.pm import PM
from pymoo.operators.crossover.pntx import TwoPointCrossover
from pymoo.operators.mutation.bitflip import BitflipMutation
from pymoo.operators.sampling.rnd import FloatRandomSampling, BinaryRandomSampling

# Real-valued: Simulated Binary Crossover + Polynomial Mutation (defaults for NSGA-II)
alg_real = NSGA2(
    pop_size=100,
    sampling=FloatRandomSampling(),
    crossover=SBX(prob=0.9, eta=15),     # eta: distribution index (higher = closer to parents)
    mutation=PM(eta=20),                  # eta: higher = smaller perturbation
    eliminate_duplicates=True
)

# Binary encoding
alg_bin = GA(
    pop_size=50,
    sampling=BinaryRandomSampling(),
    crossover=TwoPointCrossover(),
    mutation=BitflipMutation(prob=0.02),
)

print("Custom operators configured")

Module 4: Termination Criteria

Control when the algorithm stops.

from pymoo.termination.default import DefaultMultiObjectiveTermination
from pymoo.termination import get_termination

# Simple: fixed number of generations or evaluations
term_gen = get_termination("n_gen", 500)      # stop after 500 generations
term_eval = get_termination("n_eval", 10000)  # stop after 10,000 function evaluations

# Convergence-based (recommended for multi-objective)
term_conv = DefaultMultiObjectiveTermination(
    xtol=1e-8,      # design space tolerance
    cvtol=1e-6,     # constraint violation tolerance
    ftol=0.0025,    # objective space tolerance
    period=30,      # check every 30 generations
    n_max_gen=500,  # hard limit
    n_max_evals=100_000,
)

print("Termination criteria set")

Module 5: Result Analysis and Pareto Front

from pymoo.optimize import minimize
import numpy as np

problem = ZDT1()
algorithm = NSGA2(pop_size=100)
res = minimize(problem, algorithm, ("n_gen", 200), seed=42, verbose=False)

# Access results
print(f"Pareto front solutions: {len(res.F)}")
print(f"Objective values (first 3):\n{res.F[:3]}")
print(f"Decision variables (first 3):\n{res.X[:3]}")
print(f"Algorithm generations: {res.algorithm.n_gen}")

# Filter for feasibility (if constraints exist)
if res.G is not None:
    feasible = (res.G <= 0).all(axis=1)
    print(f"Feasible solutions: {feasible.sum()}/{len(feasible)}")

# Performance indicators
from pymoo.indicators.hv import HV
from pymoo.indicators.igd import IGD

ref_point = np.array([1.1, 1.1])  # reference point for HV (must dominate all solutions)
hv = HV(ref_point=ref_point)
print(f"Hypervolume indicator: {hv(res.F):.4f}")

Module 6: Visualization

import matplotlib.pyplot as plt
from pymoo.visualization.scatter import Scatter

# Scatter plot for 2D/3D Pareto fronts
plot = Scatter(title="ZDT1 Pareto Front")
plot.add(res.F, color="blue", label="NSGA-II result")
plot.show()

# Manual matplotlib plot
fig, ax = plt.subplots(figsize=(6, 5))
ax.scatter(res.F[:, 0], res.F[:, 1], s=10, color="steelblue", alpha=0.8)
ax.set_xlabel("Objective 1 (f₁)")
ax.set_ylabel("Objective 2 (f₂)")
ax.set_title("Pareto Front — ZDT1")
plt.tight_layout()
plt.savefig("pareto_front.pdf", bbox_inches="tight")
print("Saved pareto_front.pdf")

# Parallel Coordinate Plot for 3+ objectives
from pymoo.visualization.pcp import PCP

# Generate 3-objective result for visualization
from pymoo.problems import get_problem
dtlz2 = get_problem("dtlz2")
ref_dirs = get_reference_directions("das-dennis", 3, n_partitions=12)
res3 = minimize(dtlz2, NSGA3(pop_size=len(ref_dirs), ref_dirs=ref_dirs),
                ("n_gen", 200), seed=1)

pcp = PCP(title="DTLZ2 — 3 Objectives", labels=["f1", "f2", "f3"])
pcp.add(res3.F)
pcp.show()

Key Concepts

Pareto Dominance

Solution a dominates b if a is no worse than b on all objectives and strictly better on at least one. The Pareto front is the set of non-dominated solutions — there is no single "best" solution, only trade-offs. NSGA-II uses non-dominated sorting + crowding distance to maintain a diverse Pareto approximation.

Constraint Handling

pymoo uses the constraint violation approach: infeasible solutions are penalized but kept in the population. A solution with constraint violation G[i] > 0 is dominated by any feasible solution regardless of objective values. This means the algorithm first drives the population toward feasibility, then optimizes objectives.

Common Workflows

Workflow 1: Two-Objective Engineering Design

import numpy as np
from pymoo.core.problem import Problem
from pymoo.algorithms.moo.nsga2 import NSGA2
from pymoo.optimize import minimize
import matplotlib.pyplot as plt

# Beam design: minimize weight and minimize deflection
class BeamDesign(Problem):
    """
    Variables: x[0] = width (0.1–5 cm), x[1] = height (0.5–10 cm)
    Obj 1: minimize cross-sectional area (weight proxy)
    Obj 2: minimize deflection (1/I, where I = bh³/12)
    """
    def __init__(self):
        super().__init__(n_var=2, n_obj=2,
                         xl=np.array([0.1, 0.5]),
                         xu=np.array([5.0, 10.0]))

    def _evaluate(self, X, out, *args, **kwargs):
        b, h = X[:, 0], X[:, 1]
        area = b * h                      # objective 1: area (minimize)
        I = b * h**3 / 12
        deflection = 1 / I               # objective 2: deflection (minimize)
        out["F"] = np.column_stack([area, deflection])

res = minimize(BeamDesign(), NSGA2(pop_size=100), ("n_gen", 300), seed=1)

fig, axes = plt.subplots(1, 2, figsize=(12, 5))
axes[0].scatter(res.F[:, 0], res.F[:, 1], s=15, c="steelblue")
axes[0].set_xlabel("Cross-sectional area (weight)")
axes[0].set_ylabel("Deflection (1/I)")
axes[0].set_title("Pareto Front")

axes[1].scatter(res.X[:, 0], res.X[:, 1], s=15, c="coral")
axes[1].set_xlabel("Width b (cm)")
axes[1].set_ylabel("Height h (cm)")
axes[1].set_title("Design Space")

plt.tight_layout()
plt.savefig("beam_design.pdf", bbox_inches="tight")
print(f"Pareto solutions: {len(res.F)}")

Workflow 2: Algorithm Comparison with Callback

import numpy as np
from pymoo.core.problem import Problem
from pymoo.algorithms.moo.nsga2 import NSGA2
from pymoo.algorithm