Class SQPoptimizer

java.lang.Object
neqsim.process.util.optimizer.SQPoptimizer
All Implemented Interfaces:
Serializable
public class SQPoptimizer extends Object implements Serializable
Sequential Quadratic Programming (SQP) optimizer for constrained process optimization.

Implements a reduced SQP method for solving nonlinear programming (NLP) problems of the form: 

minimize   f(x)
subject to g_i(x) = 0     i = 1,...,m_e   (equality constraints)
           h_j(x) >= 0   j = 1,...,m_i   (inequality constraints)
           x_L <= x <= x_U              (variable bounds)

The algorithm iteratively solves Quadratic Programming (QP) sub-problems using the Lagrangian: 

L(x, lambda, mu) = f(x) - lambda ^ T * g(x) - mu ^ T * h(x)

At each iteration, a QP sub-problem is formed with a BFGS approximation of the Hessian of the Lagrangian, and the solution provides a search direction. An Armijo-backtracking line search on a merit function ensures global convergence.

Features

  • BFGS quasi-Newton Hessian approximation (damped update for positive definiteness)
  • Active-set QP solver for bound and linear inequality constraints
  • L1 exact penalty merit function for global convergence
  • Finite-difference gradient estimation (user can provide analytical gradients)
  • Variable bounds enforced via projection

Usage Example

SQPoptimizer sqp = new SQPoptimizer();
sqp.setObjectiveFunction(x -> computeNPV(x, process));
sqp.addEqualityConstraint(x -> massBalance(x));
sqp.addInequalityConstraint(x -> maxPressure - x[0]);
sqp.setVariableBounds(new double[] {50.0, 0.5}, new double[] {200.0, 1.0});
sqp.setInitialPoint(new double[] {100.0, 0.8});
SQPoptimizer.OptimizationResult result = sqp.solve();

References

  • Nocedal, J. & Wright, S.J., "Numerical Optimization", 2nd ed., Springer (2006), Ch. 18
  • Biegler, L.T., "Nonlinear Programming", SIAM (2010), Ch. 3-4
  • Boggs, P.T. & Tolle, J.W., "Sequential Quadratic Programming", Acta Numerica (1995)
Version:
1.0
Author:
NeqSim
See Also:

  • Field Details

    • serialVersionUID

      private static final long serialVersionUID
      Serialization version UID.
      See Also:

    • logger

      private static final org.apache.logging.log4j.Logger logger
      Logger object for class.

    • objectiveFunction

      private transient SQPoptimizer.ObjectiveFunc objectiveFunction
      Objective function: f(x) > double.

    • equalityConstraints

      private transient List<SQPoptimizer.ConstraintFunc> equalityConstraints
      Equality constraints: g_i(x) = 0.

    • inequalityConstraints

      private transient List<SQPoptimizer.ConstraintFunc> inequalityConstraints
      Inequality constraints: h_j(x) >= 0.

    • lowerBounds

      private double[] lowerBounds
      Lower bounds on variables.

    • upperBounds

      private double[] upperBounds
      Upper bounds on variables.

    • x0

      private double[] x0
      Initial point.

    • n

      private int n
      Number of decision variables.

    • maxIterations

      private int maxIterations
      Maximum SQP iterations.

    • tolerance

      private double tolerance
      Convergence tolerance on KKT conditions.

    • finiteDifferenceStep

      private double finiteDifferenceStep
      Step size for finite-difference gradient.

    • penaltyParameter

      private double penaltyParameter
      Penalty parameter for L1 merit function.

    • armijoC1

      private double armijoC1
      Armijo sufficient decrease parameter.

    • maxLineSearchIterations

      private int maxLineSearchIterations
      Maximum line search iterations.

    • hessian

      private double[][] hessian
      BFGS Hessian approximation (n x n matrix).

    • lambdaEq

      private double[] lambdaEq
      Lagrange multipliers for equality constraints.

    • lambdaIneq

      private double[] lambdaIneq
      Lagrange multipliers for inequality constraints.

  • Constructor Details

    • SQPoptimizer

      public SQPoptimizer()
      Default constructor.

    • SQPoptimizer

      public SQPoptimizer(int numVariables)
      Constructor with number of variables.
      Parameters:
      numVariables - number of decision variables

  • Method Details

    • setObjectiveFunction

      public void setObjectiveFunction(SQPoptimizer.ObjectiveFunc func)
      Set the objective function to minimize.
      Parameters:
      func - objective function

    • addEqualityConstraint

      public void addEqualityConstraint(SQPoptimizer.ConstraintFunc func)
      Add an equality constraint g(x) = 0.
      Parameters:
      func - constraint function

    • addInequalityConstraint

      public void addInequalityConstraint(SQPoptimizer.ConstraintFunc func)
      Add an inequality constraint h(x) >= 0.
      Parameters:
      func - constraint function

    • setVariableBounds

      public void setVariableBounds(double[] lower, double[] upper)
      Set variable bounds.
      Parameters:
      lower - lower bounds array (length n)
      upper - upper bounds array (length n)

    • setInitialPoint

      public void setInitialPoint(double[] x0)
      Set the initial point.
      Parameters:
      x0 - initial point array (length n)

    • setMaxIterations

      public void setMaxIterations(int maxIter)
      Set maximum number of SQP iterations.
      Parameters:
      maxIter - maximum iterations

    • setTolerance

      public void setTolerance(double tol)
      Set convergence tolerance.
      Parameters:
      tol - tolerance

    • setFiniteDifferenceStep

      public void setFiniteDifferenceStep(double step)
      Set the finite-difference step size for gradient estimation.
      Parameters:
      step - step size

    • getMaxIterations

      public int getMaxIterations()
      Get maximum iterations setting.
      Returns:
      maximum iterations

    • getTolerance

      public double getTolerance()
      Get tolerance setting.
      Returns:
      convergence tolerance

    • solve

      public SQPoptimizer.OptimizationResult solve()
      Solve the NLP problem using SQP.
      Returns:
      optimization result

    • projectToBounds

      private void projectToBounds(double[] x)
      Project point onto variable bounds.
      Parameters:
      x - point to project (modified in place)

    • computeGradient

      private double[] computeGradient(SQPoptimizer.ObjectiveFunc func, double[] x)
      Compute finite-difference gradient of a function.
      Parameters:
      func - function to differentiate
      x - evaluation point
      Returns:
      gradient vector

    • evaluateConstraints

      private double[] evaluateConstraints(List<SQPoptimizer.ConstraintFunc> constraints, double[] x)
      Evaluate a list of constraint functions.
      Parameters:
      constraints - list of constraint functions
      x - evaluation point
      Returns:
      array of constraint values

    • computeJacobian

      private double[][] computeJacobian(List<SQPoptimizer.ConstraintFunc> constraints, double[] x)
      Compute Jacobian of constraint functions via finite differences.
      Parameters:
      constraints - list of constraint functions
      x - evaluation point
      Returns:
      Jacobian matrix [m x n]

    • computeKKTError

      private double computeKKTError(double[] gradF, double[] gEq, double[] hIneq, double[][] jacEq, double[][] jacIneq, double[] x)
      Compute KKT optimality error.
      Parameters:
      gradF - gradient of objective
      gEq - equality constraint values
      hIneq - inequality constraint values
      jacEq - equality constraint Jacobian [mEq x n]
      jacIneq - inequality constraint Jacobian [mIneq x n]
      x - current point
      Returns:
      KKT error (inf-norm)

    • solveQPSubproblem

      private double[] solveQPSubproblem(double[] gradF, double[] gEq, double[] hIneq, double[][] jacEq, double[][] jacIneq, double[] x)
      Solve the QP sub-problem for the SQP search direction.

      Minimizes: 0.5 * d^T * H * d + grad_f^T * d subject to linearized constraints. Uses a simplified projected gradient approach with active-set handling.

      Parameters:
      gradF - gradient of objective
      gEq - equality constraint values
      hIneq - inequality constraint values
      jacEq - equality constraint Jacobian
      jacIneq - inequality constraint Jacobian
      x - current point
      Returns:
      search direction d

    • lineSearch

      private double lineSearch(double[] x, double[] dx, double f0, double[] gEq, double[] hIneq)
      L1 exact penalty merit function line search.
      Parameters:
      x - current point
      dx - search direction
      f0 - current objective value
      gEq - current equality constraint values
      hIneq - current inequality constraint values
      Returns:
      step length alpha

    • computeMerit

      private double computeMerit(double f, double[] gEq, double[] hIneq)
      Compute L1 exact penalty merit function.
      Parameters:
      f - objective value
      gEq - equality constraint values
      hIneq - inequality constraint values
      Returns:
      merit function value

    • updateBFGS

      private void updateBFGS(double[] x, double[] xPrev, double[] gradF, double[] gradPrev, double[] gEq, double[] hIneq)
      Damped BFGS update of the Hessian approximation.
      Parameters:
      x - new point
      xPrev - previous point
      gradF - new gradient
      gradPrev - previous gradient
      gEq - equality constraints (for Lagrangian gradient)
      hIneq - inequality constraints (for Lagrangian gradient)

    • solveLinearSystem

      private double[] solveLinearSystem(double[][] aMatrix, double[] b)
      Solve a linear system A*x = b using Gaussian elimination with partial pivoting.
      Parameters:
      aMatrix - coefficient matrix
      b - right-hand side
      Returns:
      solution vector x

    • invertMatrix

      private double[][] invertMatrix(double[][] matrix)
      Invert a matrix using Gauss-Jordan elimination.
      Parameters:
      matrix - square matrix to invert
      Returns:
      inverse matrix

    • negateVector

      private double[] negateVector(double[] v)
      Negate a vector.
      Parameters:
      v - input vector
      Returns:
      negated vector

    • dotProduct

      private double dotProduct(double[] a, double[] b)
      Dot product of two vectors.
      Parameters:
      a - first vector
      b - second vector
      Returns:
      dot product

    • matVecMult

      private double[] matVecMult(double[][] matrix, double[] vec)
      Matrix-vector multiplication.
      Parameters:
      matrix - the matrix
      vec - the vector
      Returns:
      result vector