NSGA-II From Scratch: Selection, Crowding Distance and the Pareto Front

Why multi-objective optimization?

Real engineering never optimizes just one thing. A lighter airplane wing is probably weaker. A more aerodynamically efficient airfoil may be harder to manufacture. You need to find the set of designs where you cannot improve one objective without worsening another — the Pareto front. NSGA-II (Non-dominated Sorting Genetic Algorithm II), published by Deb et al. in 2002, is the most widely used algorithm for this task with over 40,000 academic citations.

Unlike classical optimization algorithms that convert multiple objectives into a single one using weights (weighted sum), NSGA-II keeps the objectives separate and finds the entire Pareto front in a single run. This is a fundamental difference: a weighted sum gives you one solution for one choice of weights, while NSGA-II gives you the whole trade-off surface. The engineer can then browse the solutions and select the one that best balances the competing objectives based on application-specific priorities — a process that often reveals unexpected design possibilities that a single-objective formulation would never surface.

The concept of Pareto optimality was formalized by Vilfredo Pareto in the 19th century in the context of economics, but its application to engineering design optimization has only become practical in the last two decades with the rise of evolutionary algorithms and cheap parallel computing. Today, NSGA-II is the default choice for multi-objective optimization in CFD-driven design, structural optimization, and generative design workflows.

Concept 1: Pareto dominance

Solution x₁ dominates x₂ if x₁ is no worse in all objectives AND strictly better in at least one. Pareto front = set of non-dominated solutions. Example: (C_d = 0.010, C_l = 0.50) dominates (C_d = 0.012, C_l = 0.50). But (C_d = 0.010, C_l = 0.50) and (C_d = 0.008, C_l = 0.45) are both Pareto-optimal — neither dominates the other.

Mathematically (assuming minimization): f_i(x₁) ≤ f_i(x₂) for all i, and there exists at least one j such that f_j(x₁) < f_j(x₂). This is a partial order — not all solutions are comparable. Two solutions where each is better in one objective and worse in another are mutually non-dominating. This is exactly what makes multi-objective optimization fundamentally different from single-objective optimization: there is no total ordering of solutions.

The Pareto front is the set of all non-dominated solutions in the objective space. In a two-objective problem, it appears as a curve (typically convex toward the ideal point). In a three-objective problem, it is a surface. The shape of the Pareto front reveals the nature of the trade-offs between objectives — a steep section of the front means a small improvement in one objective costs a large deterioration in another, while a flat section means the opposite. This information is itself valuable to the designer beyond just selecting a final solution.

Concept 2: Non-dominated sorting

Sort into fronts. Front 1 = non-dominated, remove, Front 2 = next non-dominated. O(M*N²) algorithm. This is the first fundamental mechanism of NSGA-II. It classifies the population into fronts of increasing dominance depth:

Front 1: all non-dominated solutions in the population — the current best approximation of the Pareto front.
Remove Front 1. Front 2: non-dominated solutions among the remaining ones.
Repeat until all solutions have been assigned to a front.

The algorithm is O(M*N²): for each solution, we compute (a) domination counter n_p (how many solutions dominate it), (b) set S_p (solutions it dominates). Solutions with n_p = 0 belong to Front 1. These are removed and the domination counters of the solutions they dominated are decremented. Any whose counter reaches zero are placed in Front 2, and the process continues. This method is both elegant and efficient, avoiding the naive O(M*N³) approach of repeatedly scanning all pairs.

The key insight of Deb et al. was that this sorting can be done incrementally, exploiting the fact that removing Front 1 only affects the domination status of solutions dominated by members of Front 1 — not all solutions. In practice, this makes the non-dominated sorting step negligible in cost compared to the objective function evaluations (typically CFD or FEA simulations), which dominate the overall computational budget.

Concept 3: Crowding distance

Within same front, prefer solutions in sparse regions. For each objective: sort, compute normalized distance between neighbors, sum across objectives. Boundary solutions get infinite crowding distance. Promotes diversity — prevents the algorithm from converging to a single cluster of similar solutions.

Non-dominated sorting tells us which front each solution belongs to, but not how to compare solutions WITHIN the same front. Crowding distance is a measure of how "isolated" a solution is in the objective space. For each objective, solutions are sorted by that objective. For each solution, the normalized distance between its two neighbors is summed. This is repeated for all objectives and accumulated. Boundary solutions receive infinite crowding distance to preserve the extremes of the Pareto front.

This mechanism is essential for two reasons. First, it promotes diversity: solutions in sparsely populated regions are preferred, ensuring the algorithm explores the full extent of the Pareto front rather than collapsing to a narrow cluster. Second, it enables the elite preservation step — when a front is partially accepted to fill the next generation, crowding distance provides the tie-breaking criterion to decide which solutions from that front survive. Without crowding distance, NSGA-II would degenerate into a random selection within the last accepted front, producing irregular and unpredictable coverage of the Pareto front.

The algorithm step by step

Initialize random population P₀ of size N.
Evaluate objectives for all individuals.
Non-dominated sort P_t into fronts.
Binary tournament selection: pick 2 individuals randomly, winner has better (lower-index) front. If same front, winner has higher crowding distance. This forms the mating pool.
Simulated Binary Crossover (SBX) and Polynomial Mutation to create offspring Q_t of size N.
Combine parents and offspring: R_t = P_t U Q_t (size 2N).
Non-dominated sort R_t. Fill new population P_t+1 by adding complete fronts in order until reaching N. If a front enters partially, sort by crowding distance and take the best.
Repeat from step 3 for G generations.

This mechanism is called elitism: the best parents compete directly with their offspring for a place in the next generation, guaranteeing that the quality of the Pareto front never deteriorates. Without elitism (as in the original NSGA), good solutions found in one generation can be lost in the next due to the stochastic nature of crossover and mutation. NSGA-II's elitism ensures monotonic improvement of the approximated Pareto front.

The SBX crossover operator creates offspring near parents with a controllable spread (governed by the distribution index eta_c). Higher eta_c values produce offspring closer to the parents, resulting in more local search. Polynomial mutation applies a similar principle to individual variables, with distribution index eta_m. These operators are specifically designed for real-coded variables (as opposed to binary-coded GAs) and produce bounded offspring that respect variable constraints without requiring repair mechanisms.

Python implementation

Below is a complete, minimal implementation of the three core NSGA-II components: dominance check, non-dominated sorting, and crowding distance computation. These functions can be integrated into any genetic algorithm framework to convert it from single-objective to multi-objective optimization.

import numpy as np

def dominates(a, b):
    """Check whether solution a dominates b (minimization)."""
    return all(a <= b) and any(a < b)

def non_dominated_sort(objectives):
    """Return list of fronts, each front is a list of indices."""
    n = len(objectives)
    domination_count = np.zeros(n, dtype=int)
    dominated_set = [[] for _ in range(n)]
    fronts = [[]]
    
    for i in range(n):
        for j in range(n):
            if i == j:
                continue
            if dominates(objectives[i], objectives[j]):
                dominated_set[i].append(j)
            elif dominates(objectives[j], objectives[i]):
                domination_count[i] += 1
        if domination_count[i] == 0:
            fronts[0].append(i)
    
    i = 0
    while len(fronts[i]) > 0:
        next_front = []
        for p in fronts[i]:
            for q in dominated_set[p]:
                domination_count[q] -= 1
                if domination_count[q] == 0:
                    next_front.append(q)
        i += 1
        fronts.append(next_front)
    
    return fronts[:-1]  # Remove last empty front

def crowding_distance(objectives, front):
    """Compute crowding distance for solutions in a front."""
    n_obj = objectives.shape[1]
    distances = np.zeros(len(front))
    for m in range(n_obj):
        idx = np.array(front)
        sorted_idx = idx[np.argsort(objectives[idx, m])]
        dist_pos = {v: k for k, v in enumerate(front)}
        distances[dist_pos[sorted_idx[0]]] = np.inf
        distances[dist_pos[sorted_idx[-1]]] = np.inf
        f_range = objectives[sorted_idx[-1], m] - objectives[sorted_idx[0], m]
        if f_range == 0:
            continue
        for k in range(1, len(front) - 1):
            pos = dist_pos[sorted_idx[k]]
            distances[pos] += (objectives[sorted_idx[k+1], m] -
                              objectives[sorted_idx[k-1], m]) / f_range
    return distances

The dominates function is the atomic operation of the entire algorithm — it is called O(N²) times during sorting. Its implementation as a vectorized NumPy comparison makes it fast enough for population sizes up to several hundred. Note that the function assumes minimization for all objectives; for maximization, simply negate the objective values before passing them to NSGA-II.

The non_dominated_sort function implements the O(M*N²) algorithm from the original NSGA-II paper. It first builds the domination graph by comparing all pairs of solutions, then iteratively extracts fronts by collecting solutions whose domination counter reaches zero. The returned list of fronts is ordered from best (Front 1) to worst.

The crowding_distance function computes diversity scores for solutions within a single front. For each objective, solutions are sorted and the normalized distance between neighbors is accumulated. Boundary solutions receive infinite distance to ensure they are always preferred, preserving the spread of the Pareto front approximation. The normalization by f_range makes the metric scale-invariant across objectives with different units and magnitudes.

Pareto front visualization

For problems with 2 objectives, the Pareto front is a curve in the objective space. For 3 objectives, it is a surface. For higher dimensions, visualization becomes challenging — parallel coordinate plots or dimensionality reduction (PCA) are typically used.

The evolution of the Pareto front across generations is one of the most informative visualizations. Generation 1 shows a scattered cloud of points; by Generation 50, a well-defined front emerges approaching the true Pareto-optimal surface. The hypervolume indicator — the volume of objective space dominated by the approximated front relative to a reference point — provides a scalar metric of convergence and diversity. A stabilizing hypervolume curve indicates that the algorithm has converged and further generations are unlikely to produce significant improvement.

For three-objective problems, 3D scatter plots with a reference plane or convex hull visualization effectively communicate the trade-off surface. Interactive visualization tools that allow rotating, zooming, and selecting individual designs from the front are particularly powerful for design reviews, enabling stakeholders to explore the trade-off space directly rather than having the engineer pre-select a handful of candidate designs.

NSGA-II in practice: parameters

Population size: 50-200 for 2-3 objectives, up to 500 for 5+ objectives. Larger populations maintain diversity but increase computational cost linearly.
Generations: 50-300. Convergence is typically reached when the hypervolume indicator stabilizes — plotting hypervolume vs generation is the standard diagnostic.
Crossover probability: 0.9. Mutation probability: 1/n_vars (one variable mutated per offspring on average).
Distribution indices: eta_c = 20 (crossover), eta_m = 20 (mutation). Higher values produce offspring closer to parents, favoring local search over exploration.
Computational cost: dominated by objective evaluation (e.g., running CFD). NSGA-II overhead is negligible compared to simulation cost — a microsecond of sorting versus minutes or hours of CFD.

The population size choice involves a trade-off: too small and the algorithm may miss large portions of the Pareto front due to insufficient diversity; too large and convergence slows because more generations are needed to propagate good genes through the population. A rule of thumb is 50-100 for two objectives and 100-200 for three objectives. For problems with expensive evaluations (CFD, FEA), smaller populations with more generations are preferred to maximize the information gained per simulation.

Limitations and extensions

NSGA-II struggles beyond 3-4 objectives. The reasons are twofold: first, as the number of objectives increases, almost all solutions become non-dominated (the dominance relation loses discriminatory power), making front-based selection nearly random. Second, crowding distance loses effectiveness in high dimensions because the volume grows exponentially and solutions naturally spread out, making density estimation unreliable.

NSGA-III (Deb & Jain, 2014) addresses this with reference-point-based selection. Instead of crowding distance, it projects solutions onto a set of uniformly distributed reference directions on a normalized hyperplane and selects solutions that are closest to each reference direction. This ensures coverage of the entire Pareto front even in high-dimensional objective spaces. NSGA-III has been successfully applied to problems with up to 15 objectives.

For expensive objectives (CFD, FEA), surrogate-assisted NSGA-II is preferred: train a Gaussian Process model on already-evaluated designs and use it to pre-filter candidates before running the real simulation. The GP provides both a prediction and an uncertainty estimate, enabling principled decisions about which designs are worth evaluating. This approach can reduce the number of required CFD evaluations by an order of magnitude, making multi-objective optimization feasible for problems where each evaluation takes hours.

NSGA-II from scratch: the algorithm that optimizes without choosing