Optimization Algorithms in Density-Based Topology Optimization

Density-based topology optimization relies on continuous design variables to represent material distribution, and various optimization algorithms have been developed to update these densities effectively.

One of the most widely used methods is the Optimality Criteria (OC) method, which updates the density based on a closed-form expression derived from the Karush-Kuhn-Tucker (KKT) conditions. The OC method is simple to implement and provides fast convergence in many structural problems. It typically uses a multiplicative update rule that enforces a volume constraint while improving the compliance.

Another approach is the Modified Optimality Criteria (MOC) method, which extends the OC framework by incorporating additional strategies for handling constraints, improving stability, or enhancing convergence. For example, some MOC variants integrate projection and filtering directly into the update step or apply move limits and continuation schemes to control intermediate densities more robustly.

In addition to these, gradient-based methods using standard optimization libraries (e.g., MMA: Method of Moving Asymptotes) are also popular in density-based formulations, particularly when dealing with multiple constraints or noncompliance objectives.

Each of these methods offers different trade-offs in terms of implementation complexity, convergence speed, and robustness. The choice of optimization algorithm can significantly affect the quality and performance of the final design.

Optimality Criteria

Optimality Criteria Method (OC)

The Optimality Criteria (OC) method is a widely used algorithm in density-based topology optimization due to its simplicity and effectiveness. It provides a closed-form update rule for the design variables based on the optimality conditions derived from the Lagrangian of the problem.

Assume the goal is to minimize compliance \(C(\rho)\) subject to a volume constraint:

\[\begin{split}\min_{\rho} \quad C(\rho) = \mathbf{f}^T \mathbf{u}(\rho) \\ \text{subject to} \quad \sum_e \rho_e v_e \leq V_{\text{max}}, \quad \rho_{\text{min}} \leq \rho_e \leq 1\end{split}\]

The KKT conditions yield an optimality criterion involving the derivative of the objective function and a Lagrange multiplier \(\lambda\). Based on this, the OC method defines the update rule for each element \(e\) as:

\[\rho_e^{(t+1)} = \text{clip}\left( \rho_e^{(t)} \cdot \left( \frac{-\partial C / \partial \rho_e}{\lambda} \right)^{\eta},\ \rho_{\text{min}},\ 1 \right)\]

where:

  • \(\partial C / \partial \rho_e\) is the sensitivity of compliance with respect to the design variable,

  • \(\lambda\) is chosen to satisfy the volume constraint (e.g., by bisection),

  • \(\eta\) is a numerical damping or scaling factor (commonly \(\eta = 1\)),

  • clip ensures the updated value stays within allowable bounds.

The OC update rule is applied iteratively, and move limits may also be enforced to restrict the change in \(\rho_e\) between iterations.

Advantages

  • Simplicity: Easy to implement and understand.

  • Efficiency: Fast convergence in practice for compliance minimization.

  • No gradient solver needed: Only objective sensitivities are required, not the full gradient descent machinery.

Disadvantages

  • Not general-purpose: Hard to extend to problems with multiple constraints or noncompliance objectives.

  • Heuristic parameters: Requires tuning of parameters like damping factor, move limits.

  • Volume constraint enforcement: Inexact unless a Lagrange multiplier search is properly implemented.

  • Limited to multiplicative updates: Cannot easily incorporate modern optimization strategies (e.g., trust regions or line search).

Despite its limitations, the OC method remains popular for compliance-based problems, especially when a quick and robust heuristic method is desired for early prototyping or academic demonstration.

Modified Optimality Criteria (MOC) Variants

In density-based topology optimization, the Modified Optimality Criteria (MOC) method can be implemented in several ways. I implemented two variants: a log-space update and a direct additive update. These differ in how they incorporate sensitivity information and handle volume constraints.

1. Log-space Update Method

This method modifies the OC update by applying it in log-space. Instead of directly updating the density \(\rho\), it updates \(\log \rho\), promoting better numerical stability and capturing multiplicative behavior.

The Lagrangian gradient \(dL\) is computed as the sum of the compliance sensitivity \(dC\) and a volume penalty term \(\lambda_v\). The update is then applied in log-space:

\[\log \rho^{(t+1)} = \log \rho^{(t)} + \eta \cdot \log \left( \frac{-dC}{\lambda_v} \right)\]

The updated density is then recovered via exponentiation:

\[\rho^{(t+1)} = \rho^{(t)} \cdot \left( \frac{-dC}{\lambda_v} \right)^{\eta}\]

Here:

  • \(\rho^{(t)}\) is the current density field,

  • \(dC\) is the sensitivity of the compliance with respect to density,

  • \(\lambda_v\) is the derivative of the volume penalty,

  • \(\eta\) is a step size or learning rate.

This approach improves numerical stability and ensures that the density remains positive throughout the optimization process. Importantly, this means the optimization behaves multiplicatively rather than additively.

In this multiplicative regime, using a simple linear average to smooth the dual variable \(\lambda_v\) can lead to instability or a misrepresentation of the underlying dynamics. Instead, the Exponential Moving Average (EMA) is more appropriate because:

  • EMA applies a multiplicative-style memory, giving more weight to recent changes while still preserving the trend of past updates.

  • It provides a smooth yet responsive estimate of the dual variable, crucial for stable convergence.

  • It avoids sharp oscillations that may occur when \(\lambda_v\) changes abruptly, which could cause large exponential shifts in the density update.

The EMA of the dual variable is computed as:

\[\lambda_v^{(t)} = \lambda_\text{decay} \cdot \lambda_v^{(t-1)} + (1 - \lambda_\text{decay}) \mu_{p} \cdot \hat{\lambda}_v^{(t)}\]

where \(\hat{\lambda}_v^{(t)}\) is the current value (e.g., based on volume constraint violation) and \(\lambda_\text{decay} \in (0, 1]\) is the smoothing factor.

Since MOC operates in a logarithmic or ratio-based update regime, EMA naturally complements this behavior by ensuring the dual variable evolves in a similarly smooth and proportional manner.

The objective function and the original constraint remain unchanged in the formulation:

However, in the implementation of MOC with EMA, the control of the dual variable \(\lambda_v^{(t)}\) is based on a relative (ratio-based) constraint violation, rather than an absolute difference.

Advantages:

  • Efficient in-place update suitable for large-scale problems.

  • Ensures positivity of design variables automatically.

  • Straightforward to implement using standard vector operations.

Disadvantages:

  • Requires careful tuning of \(\eta\) and \(\lambda_v\).

  • Does not enforce volume constraints exactly—relies on penalty balancing.

  • Convergence behavior may vary depending on the problem and filter.

2. Linear-Space Update Method

This variant formulates the update as an explicit increment \(\Delta \rho\) added to the current density. The update is defined as:

\[\rho^{(t+1)} = \rho^{(t)} + \Delta \rho\]

where:

\[\Delta \rho = -\eta \cdot (dC + \lambda_v)\]

This form is structurally simpler and better suited for integration with additional constraint handling mechanisms. The update increment \(\Delta \rho\) may be further clipped to enforce move limits or bound constraints.

Advantages:

  • Direct control over the update magnitude.

  • Easier to incorporate projection, filtering, or move limits.

  • Simple to interpret and modify in algorithmic experiments.

Disadvantages:

  • Can violate positivity if not carefully bounded.

  • Volume constraint is only approximately satisfied unless post-processing is added.

  • Requires stabilization strategies (e.g., clipping or damping) for robust performance.

Summary

Both update strategies aim to descend along the objective gradient while respecting volume constraints and maintaining stability. The choice depends on implementation goals:

  • Use the log-space update when prioritizing positivity and multiplicative structure.

  • Use the additive update when flexibility, constraint control, or custom damping is desired.

Log-space Lagrangian Method

Log-space Lagrangian Method

This method is a variant of the Modified Optimality Criteria (MOC) approach that performs density updates in log-space rather than directly in the physical domain. The goal is to improve numerical stability and ensure that the updated density field remains strictly positive.

Instead of updating the density \(\rho\) directly, the update is applied to its logarithm \(\log \rho\). This transforms the multiplicative update behavior into an additive one in log-space.

The method computes the Lagrangian gradient \(dL\) as the sum of the compliance sensitivity \(dC\) and the derivative of the volume penalty term \(\lambda_v\). The update rule is:

\[\log \rho^{(t+1)} = \log \rho^{(t)} - \eta \cdot (dC + \lambda_v)\]

The density is then recovered by exponentiation:

\[\rho^{(t+1)} = \exp\left( \log \rho^{(t)} - \eta \cdot (dC + \lambda_v) \right) = \rho^{(t)} \cdot \exp\left( -\eta \cdot (dC + \lambda_v) \right)\]

Here:

  • \(\rho^{(t)}\) is the current density at iteration \(t\),

  • \(dC\) is the derivative of compliance with respect to density,

  • \(\lambda_v\) is the derivative of the volume constraint penalty,

  • \(\eta\) is a scalar step size (analogous to a learning rate).

Optional clipping is applied in log-space to limit excessive updates and preserve stability:

\[\log \rho^{(t+1)} = \text{clip}\left( \log \rho^{(t+1)},\ \log \rho_{\min},\ \log \rho_{\max} \right)\]

Finally, move limits can also be enforced using:

\[\rho^{(t+1)} = \text{clip}\left( \rho^{(t+1)},\ \rho^{(t)} - \Delta \rho_{\max},\ \rho^{(t)} + \Delta \rho_{\max} \right)\]

Advantages:

  • Naturally ensures \(\rho > 0\) without additional constraints.

  • Suitable for in-place and vectorized computation in large-scale problems.

  • Converts multiplicative effects into additive updates, improving numerical robustness.

Disadvantages:

  • Sensitive to step size \(\eta\) and penalty weight \(\lambda_v\).

  • Volume constraints are only enforced implicitly via penalty.

  • Requires careful initialization and parameter tuning to ensure convergence.