# Packages suitable for numerical optimization of functions with discontinuous gradient at the point of minimum

Are there packages for numerical optimization in julia or python, or in any other system for scientific computing, capable of taking into account the discontinuity of gradient at the minimum point? Functions that I need to minimize have form $$F(X) = \mbox{max}(f_1(X), \ldots, f_n(X))$$ where $$f_i(X)$$ are smooth functions of $$X\in R^d$$. Actually, in my particular problem, $$f_i(X)$$ are unordered eigenvalues of large real symmetric sparse matrices. So $$F(X)$$ is a maximum eigenvalue. My matrix family is such that level crossing definitely takes place under change of parameters. Hence my formulation of problem with discontinuous gradient.

I consider functions $$F(X)$$ described above as multidimensional generalizations of the following function of one variable $$F(x) = |x| = \mbox{max}(x, -x),\quad x\in R,$$ which has non-zero and even discontinuous derivative at the point of minimum, $$x_0 = 0$$. In this last example $$f_1(x) = x$$ and $$f_2(x) = -x$$. In general case, functions $$f_i$$ are non-linear.

• Are you trying to find the maximum or the minimum of $F(X)$? You only say "optimize", but it makes a difference whether you want the max or the min :-) Jul 18, 2023 at 15:55
• @WolfgangBangerth, I need to find the minimum value of F(X), while this function is defined by taking maximum.
– Gec
Jul 18, 2023 at 17:07
• Maybe packages like github.com/JuliaSmoothOptimizers/RegularizedOptimization.jl could be of interest here.
– Tanj
Jul 18, 2023 at 18:23
• What specifically is (the multidimensional generalization of) $f_i$? If it is a norm, then the problem can be handled by Convex.jl or CVXPY, among others. If it s a 2-norm , then specifically as a Second Order Cone Problem (SOCP). If the 1 or $\infty$ norm, then as a Linear Programming (LP) problem. Those optimization modeling systems can accept the max of norms (of affine arguments) n the objective of such a problem, and will automatically use an under the hood epigraph reformulation as needed. Jul 18, 2023 at 18:52
• Is $F(X)$ the max of the eigenvalues of the symmetric matrix? If so,, Convex.jl (eigmax) or CVXPY, (lambda_max) among others can handle, using the built-in "atomic" functions I put in parentheses. Jul 18, 2023 at 19:21

Based on the OP's comments, it appears that the problem is to minimize the maximum eigenvalue of a symmetric matrix, which I will assume to be an affine (linear) function of the optimization variable. If so, then Disciplined Convex Programming (DCP) tools, such as Convex.jl (Julia), CVXPY (Python), CVXR (R), and CVX (MATLAB) can be used to formulate the problem and call a conic solver supporting the semidefinite cone to solve it.

These optimization modeling systems automatically reformulate the maximum eigenvalue to a semidefinite constraint, and then call a senidefinite solver to solve the problem.

If the matrix is not an affine function of the optimization problem, it is likely a difficult non-convex (Nonlinear SDP) problem, for which good, reliable solvers are of limited availability, and the YALMUP optimization modeling system might be your best bet. If you don't even have an explicit algebraic formula for the matrix X as a function of the optimization variable, the problem is even more difficult.

If X is the affine symmetric matrix, the maximum eigenvalue of X is represented as: eigmax(X) in Convex.jl

lambda_max(X) in CVXPY

lambda_max(X) in CVXR

lambda_max(X) in CVX

Constraints, if any, would have to be compatible with these optimization modeling systems in order to use then. If so, the optimization problem is convex, so the solution found should be the (a) global minimum, barring numerical issues.

• Thank you very much! Your answer perfectly matches my problem. Now I see that my problem can be formulated as $\mbox{min}_X(\mbox{eigmax}(\sum_i x_i A_i))$, where $A_i$ are fixed matrices. And I didn't know that this problem is convex. Convex.jl with SCS optimizer really finds the answer in the case of not too large matrices. Alas, I need a solution for large matrices.
– Gec
Jul 19, 2023 at 7:18
• How large are the matrices? Mosek is likely a faster solver for your problem, presuming there is enough memory available. Jul 19, 2023 at 11:18
• Matrix size : "352716×352716 SparseMatrixCSC{Float64, Int64} with 15502706 stored entries". Dimension of $x$ is 12. SCS in Convex.jl have given "not enough memory" message. My laptop has 16Gb of RAM.
– Gec
Jul 19, 2023 at 13:30
• What is the dimension of the $x_i's$, i.e., the optimization variables? Jul 19, 2023 at 13:50
• Try using DualOptimizer, which should solve the dual problem instead of the primal. it might require much less memory and be much faster on this problem, given the dimensions of x and the matrices. Jul 19, 2023 at 14:01

This problem is easy to reformulate into a smooth one with slack variables. Finding the minimum of $$F(X)$$ is equivalent to finding the minimizer of the following problem: $$\min_{z,X} z \\ \text{subject to } z\ge f_i(X), i=1,\ldots,n.$$ That is because each of the constraints ensures that $$z$$ is bigger than $$f_i(X)$$, and consequently $$z$$ is bigger than the largest one -- that is, $$z\ge F(X)$$. Minimizing $$z$$ ensures that $$z=F(X)$$ and so the two optimization problems are equivalent.

This gets rid of the non-smoothness of the problem. But, if the functions $$f_i(X)$$ are not convex, then you still end up with a non-convex problem which is, like all non-convex problems difficult to solve because it may have multiple local solutions and you are in all likelihood interested in the global solution.

In other words, I can help you with the non-smoothness, but not with the non-convexity of the problem.

• Many thanks for your answer. I am not familiar with optimization theory approaches, so it is useful for me to see such reformulation. I've edited my question, so you can see that $f_i(X)$ are eigenvalues of large real symmetric sparse matrix dependent on $X$. Actually, I can only find $F(X)$ and it's gradient numerically. So I am afraid your approach can not help in solving my problem because of I don't have simple expressions for $f_i(X)$.
– Gec
Jul 18, 2023 at 22:42

It depends on the specifics of your problem, but one potential approach is to use gradient descent. From a pragmatic perspective, gradient descent does not actually require a continuous gradient and can often work fine even if the gradient is discontinuous. Choosing a good initialization (or repeating multiple times with multiple different initializations) may be helpful, depending on the structure of the optimization landscape.

Since your functions are non-convex and the gradient is not smooth, there are no guarantees. It might fail badly. However, in some applications gradient descent can work in practice even when the theory doesn't offer any guarantees. The only way to find out is to try it.

• In fact, I got a satisfactory numerical solution to my problem using the L-BFGS method from the Optim.jl package. Due to the discontinuity of the gradient at the minimum point, the gradient-based optimization process fluctuates at the end and does not stop. An obvious practical solution is to limit the number of steps in the optimization process. I started this topic because I thought there might be a more scientific solution to this problem.
– Gec
Jul 19, 2023 at 18:14
• @Gec, L-BFGS is different from gradient descent, and I expect L-BFGS might be more sensitive to discontinuity in the gradient.
– D.W.
Jul 20, 2023 at 7:21