Ksenia Syrtseva
email: ksenia-k.syrtseva@edf.fr

Welington De Oliveira

Sophie Demassey

Wim Van Ackooij

Minimizing the difference of convex and weakly convex functions via bundle method

Keywords: Non-smooth Optimization, DC Programming, Bundle Methods, weakly convex Mathematics Subject Classification: 90C26, 65K05, 49J52, 49J53

We consider optimization problems with objective and constraint being the difference of convex and weakly convex functions. This framework covers a vast family of nonsmooth and nonconvex optimization problems, particularly those involving certain classes of composite and nonconvex value functions. We investigate several stationary conditions and extend the proximal bundle algorithm of [van Ackooij et al., Comput. Optim. Appl., 78 (2021), pp. 451-490 ] to compute critical points for problems in this class. Our modifications on that algorithm boil down to a different master program and an original rule to update the proximal parameter to ensure convergence in this more general setting. Thanks to this new rule, no pre-estimation of the underlying weakly-convex moduli is needed, opening the way to deal with optimization problems for which no practical and mathematically sound algorithms exist. Numerical experiments on some nonconvex stochastic problems illustrate the practical performance of the method.

1 Introduction.

This work presents a bundle method for nonsmooth and nonconvex optimization problems of the form min x∈X f (x) s.t. c(x) ≤ 0, (1.1a) where X is a nonempty bounded polyhedron contained in an open convex set O ⊂ R n , and functions f : O → R and c : O → R are decomposable as the difference of convex and (locally) weakly convex functions. More specifically, we assume that the following convex-weakly convex (CwC) decompositions are available:

f (x) = f 1 (x) -f 2 (x) and c(x) = c 1 (x) -c 2 (x), (1.1b)
with f 1 , c 1 : O → R convex and f 2 , c 2 : O → R weakly convex functions on some neighbourhood of each x ∈ O. We adopt the more general definition of weakly convex functions (see Definition 2.2 below) given in [START_REF] Vial | Strong and weak convexity of sets and functions[END_REF]Def. 4.2] so that we can exploit the equivalence between the families of locally weakly convex and Lower-C 2 functions [33, Thm. 1.3, Cor. 1.3] to highlight the breadth of our approach. In particular, we have in mind the following settings for f 2 (as well as for c 2 ): i) f 2 (x) = ϕ(x) is a (possibly nonsmooth) convex function;

ii) f 2 (x) = -h(x) with h having Lipschitz continuous gradient;

iii) f 2 (x) = ϕ(x) -h(x), with ϕ and h as given above; iv) f 2 (x) is the optimal value of max t∈T F (t, x), with T a (possibly nonconvex) compact set and F of class C 2 ; v) f 2 (x) = ϕ(G(x)), with ϕ : R m → R convex and Lipschitz and G : R n → R m a smooth mapping with Lipschitz Jacobian.

1

Analogous settings for c 2 , and their combinations with the ones for f 2 , are covered by our analysis (see Section 2 below for details).

Weakly convex functions enjoy favorable properties in so much as that they can be recast as Differenceof-Convex (DC) functions [START_REF] Hartman | On functions representable as a difference of convex functions[END_REF]. Hence, problem (1.1) can, in theory, be recast as a DC-constrained DC program, a setting that proves practical if explicit DC decompositions are available; see for instance [START_REF] Le Thi | Exact penalty and error bounds in DC programming[END_REF][START_REF] Khalaf | A DC optimization-based clustering technique for edge detection[END_REF][START_REF] Strekalovsky | A local search method for optimisation problem with d.c. inequality constraints[END_REF][START_REF] Montonen | Bundle-based descent method for nonsmooth multiobjective DC optimization with inequality constraints[END_REF][START_REF] Gaudioso | Minimizing nonsmooth DC functions via successive DC piecewise-affine approximations[END_REF][START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF][START_REF] Van Ackooij | Addendum to the paper 'nonsmooth DC-constrained optimization: constraint qualification and minimizing methodologies[END_REF][START_REF] Sun | Algorithms for difference-of-convex programs based on difference-of-moreauenvelopes smoothing[END_REF] and references therein. However, if no DC decomposition is known for f or c, the DC machinery is unsuitable, and the methods proposed in these references are not applicable. This is already the case for the more straightforward items ii) and iii) above if the underlying Lipschitz constant is unknown and no upper bound is readily available 1 . The situation becomes even more complicated for items iv) and v): in general, there are no formulae, rules, or practical insights to obtain a DC decomposition for f 2 in these cases (see Example 1.1 below for a particular case of iv). A strategy to handle problem (1.1) via DC programming algorithms is to replace functions f i and c i (i = 1, 2) with f i (x) + µ 2 ∥x∥ 2 and c i (x) + µ 2 ∥x∥ 2 for a large parameter µ > 0 estimating upper bounds on the unknown weakly-convex moduli µ f and µc of f 2 and c 2 (see Proposition 2.4), hoping that f 2 (x) + µ 2 ∥x∥ 2 and c 2 (x) + µ 2 ∥x∥ 2 are convex on X. As, in general, there is no reliable way to assert the convexity of these latter functions, DC programming algorithms applied in this context must be understood as heuristics. Remarkably, the work [START_REF] Van Ackooij | Nonsmooth and nonconvex optimization via approximate difference-of-convex decompositions[END_REF] exploits such a strategy by combining a dynamic rule to update µ with a nonconvexity test so that convergence is achieved, but only in a probability sense. Differently, for a class of nonconvex two-stage stochastic problems, the authors of [START_REF] Li | A decomposition algorithm for two-stage stochastic programs with nonconvex recourse[END_REF] exploit an implicitly convex-concave structure of the objective function and propose an algorithm based on the so-called partial Moreau envelope that disregards DC decompositions at the price of nonnegligible computational costs.

In contrast to the above references, this work investigates a bundle method approach for tackling (1.1), which neither requires explicit DC decompositions of the involved functions (in particular, bounds on the weakly-convex moduli µ f and µc need not be known), nor relies on (often costly) Moreau envelopes. For the method to work, it suffices to dispose of a difference of convex and weakly convex (CwC) decomposition of the involved functions, as in (1.1b). Compared to DC, the latter structure appears more naturally in applications (see [15, § 7.5]) and has yet to be exploited to design optimality conditions and numerical algorithms. This work aims to fill this gap.

Our approach broadens and enhances the method proposed in [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] for dealing with DC-constrained DCproblems in the following two ways. First, the availability of DC decompositions is no longer needed, which makes our approach applicable to a larger scope of problems in practice. Second, it is ensured to compute critical points for the original problem without any additional assumption on weakly-convex components (f 2 and c 2 need not be continuously differentiable as assumed in [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF]Thm. 2]). In addition, it has a lower cost per iteration (the master subproblem has fewer constraints than the one of [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF]). Similarly to [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF], our approach builds upon a problem reformulation via improvement function, a well-known and successful strategy in the nonsmooth optimization literature [START_REF] Sagastizábal | An infeasible bundle method for nonsmooth convex constrained optimization without a penalty function or a filter[END_REF][START_REF] Apkarian | Mixed H 2 /H∞ Control via Nonsmooth Optimization[END_REF][START_REF] Montonen | Bundle-based descent method for nonsmooth multiobjective DC optimization with inequality constraints[END_REF]. However, due to the above modifications, the convergence analysis of our extension of the method proposed in [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] must be done anew. Furthermore, a new criticality definition for the reformulated problem links directly with (necessary) optimality conditions for the original problem (1.1), which makes it a major ingredient for these enhancements. Such a criticality concept is introduced and analyzed in Section 3 below, where we also extend the alternative characterization of Bouligand stationarity given in [START_REF] Pang | Computing B-stationary points of nonsmooth DC programs[END_REF] to our CwC setting. Before that, we motivate this work with the following example that presents a class of problems (of great practical appeal) where the CwC decomposition arises upon applying a well-known interior-penalty strategy. An example of a real-life (chance-constrained optimal power flow) problem fitting our CwC structure without any approximation can be found in the Ph.D. thesis [15, § 7.5]. (1.2)

Assume that:

-f 1 , c 1 , c 2 : R n → R are convex (possibly nonsmooth) functions;

-X ⊂ R n , Y ⊂ R n 2 are two (non-empty) convex and bounded polytopes;

-q, ψ i : R n × R n 2 × Ξ → R, i = 1, . . . , m, possess the following characteristics: q(•, •, ξ) and ψ i (•, •, ξ) are twice-continuously differentiable for every ξ ∈ Ξ fixed and, moreover, q(x, •, ξ) and ψ i (x, •, ξ) are convex for every x and ξ fixed;

-the constraints in the subproblem Q(x; ξ) satisfy the Slater condition: for every x ∈ X and ξ ∈ Ξ, there exists y • (x; ξ) ∈ Y such ψ i (x, y • (x; ξ), ξ) < 0, i = 1, . . . , m. As presented in [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF], the DC constraint c 1 (x) -c 2 (x) ≤ 0 above is particularly useful in this stochastic programming setting to model chance constraints.

Under the above assumptions, evaluating the recourse function Q(x; ξ) amounts to solving a well-defined convex optimization problem on variable y. Although this essential property is present, the recourse function itself fails to be convex on variable x (but Q(•; ξ) is continuous as a result of [START_REF] Bonnans | Perturbation Analysis of Optimization problems[END_REF]Prop. 4.4]). Furthermore, without further assumptions, computing a (generalized) subgradient of Q(•; ξ) at x as well as asserting additional properties about this function are challenging tasks. This could for instance be done if the constraints satisfy a further Aubin or Lipschitz like property upon exploiting [START_REF] Mordukhovich | Variational Analysis and Applications[END_REF]Chapter 4]. Still though, most likely, at best we would be dealing with subdifferentials inclusions -and concrete algorithms to handle such general "marginal functions" would be unavailable.

A possible manner to curtail these difficulties is to approximate the recourse function with a more tractable one. As explained in [START_REF] Borges | A regularized smoothing method for fully parameterized convex problems with applications to convex and nonconvex two-stage stochastic programming[END_REF], with the help of the log-barrier penalty function and a penalization parameter ε > 0, we may approximate Q(x; ξ) with

Q ε (x; ξ) := min y∈Y q(x, y; ξ) - 1 ε m i=1 log(-ψ i (x, y; ξ)). (1.3)
Given the above assumptions, it is well known that Q ε (x; ξ) ↓ Q(x; ξ) as ε ↓ 0 (e.g., [3, § 2.2] and [2, p. 266]), and thus the model

min x∈X f 1 (x) -f 2 (x) s.t. c 1 (x) -c 2 (x) ≤ 0 with f 2 (x) := S s=1 πs[-Q ε (x; ξ s )]
is an accurate approximation of (1.2) when ε > 0 is small enough. Furthermore, as

-Q ε (x; ξ) = max y∈Y 1 ε m i=1 log(-ψ i (x, y; ξ)) -q(x
, y; ξ) is a weakly convex function (c.f. item iv) above), this model fits the structure (1.1). We highlight that Q ε (x; ξ) is generally a nonsmooth (nonconvex) function; hence, the above problem is challenging. To our knowledge, no practical and mathematically sound optimization algorithm could tackle this class of problems before this work. Indeed, [START_REF] Mine | A minimization method for the sum of a convex function and a continuously differentiable function[END_REF] requires f 2 to be smooth, [START_REF] Li | A decomposition algorithm for two-stage stochastic programs with nonconvex recourse[END_REF] assumes q(•, •, ξ) and ψ i (•, •, ξ) to be concave-convex functions, and [START_REF] Borges | A regularized smoothing method for fully parameterized convex problems with applications to convex and nonconvex two-stage stochastic programming[END_REF] requires another degree of approximation by adding a Tikhonov regularization to (1.3) to force Q ε (x; ξ) be smooth. In all these references, function c 2 is absent. Being nonsmooth, we mention in passing that a (generalized) subgradient of f 2 at x can be computed and seen to be S s=1 πsg(y(x; ξ s )), where g(•)

:= ∇x[ 1 ε m i=1 log(-ψ i (x, •; ξ s )) -q(x, •; ξ s )]
is the gradient w.r.t. x of the objective function of (1.3) multiplied by -1, and y(x; ξ) is an arbitrary optimal solution of (1.3) (see Proposition 2.1 and [START_REF] Rockafellar | Generalized subgradients in mathematical programming[END_REF]Thm. 7.3]). □

Our approach is still applicable in more general case, where the probability vector π in the above example is a function (of class C 2 ) of the first-stage variable x, i.e., πs(x), i = 1, . . . , S. Hence, this work's class of optimization problems includes the challenging family of stochastic programming recourse models with decision-dependent uncertainty considered, for instance, in [START_REF] Hellemo | Decision-dependent probabilities in stochastic programs with recourse[END_REF] and [START_REF] Li | A decomposition algorithm for two-stage stochastic programs with nonconvex recourse[END_REF].

The remaining of this manuscript is organized as follows. Section 2 recalls essential definitions, key elements, and well-known concepts from variational analysis. Necessary optimality conditions for problem (1.1) are presented in Section 3 as well as the problem reformulation via an improvement function. Once the link between the reformulated and the original problem is established in the same section, Section 4 focuses on an improvement-function-based bundle method for problem (1.1). Section 5 presents the method's convergence analysis to critical points, and finally, Section 6 illustrates the practical performance of our approach on some nonconvex stochastic optimization problems.

Notation

The following notation is employed throughout the text. For a real number a, we denote by [a] + the value max{a, 0}. For any points x, y ∈ R n , ⟨x, y⟩ stands for the Euclidean inner product, and ∥ • ∥ for the associated norm, i.e., ∥x∥ = ⟨x, x⟩. For a convex set X, N X (x) stands for its normal cone at the point x, i.e., the set {y : ⟨y, z -x⟩ ⩽ 0 for all z ∈ X} if x ∈ X and the empty set otherwise. The Bouligand tangent cone to a (possibly nonconvex) set W ⊂ R n at a point w ∈ W is the set T W (w) of all tangent directions in the following sense: d ∈ T W (w) if there exist a sequence of vectors {w k } ⊂ W converging to w and a sequence of positive scalars

t k → 0 such that d = lim k→∞ (w k -w)/t k . The indicator function of X ⊂ R n is defined as i X (x) = 0 if x ∈ X and i X (x) = +∞

Definition and prerequisites

This section starts by recalling the concept of (generalized) directional derivatives and subdifferentials. Basic subdifferential calculus is summarized in Proposition 2.1 below, followed by the definitions of weakly convex and lower-C 2 functions. The section closes with Proposition 2.4 asserting that the definition of (locally) weakly convex function can be globally extended to the whole convex and compact set X. Such a property is of crucial importance in this work.

A function f : O → R is said to be locally Lipschitz continuous if for each x ′ ∈ O there is a neighbourhood

V x ′ ⊂ O of x ′ such that, for some L x ′ ≥ 0, |f(x) -f(y)| ≤ L x ′ ∥x -y∥ ∀ x, y ∈ V x ′ .
The function f is said to be Lipschitz continuous on O if L x ′ = L can be taken independent of x ′ ∈ O, and V x ′ in the above inequality is replaced with O.

Directional derivatives and subdifferentials. Let ϕ : O → R be a convex function. Then ϕ is locally Lipschitz continuous and, for each x ∈ O, the directional derivative ⟨s, d⟩,

ϕ ′ (x; d) := lim τ ↓0 ϕ(x + τ d) -ϕ(x)
where ∂ϕ(x) is the subdifferential of ϕ at point x:

∂ϕ(x) := {s ∈ R n : ϕ(y) ≥ ϕ(x) + ⟨s, y -x⟩ ∀ y ∈ R n } . (2.1)
The elements of ∂ϕ(x) are referred to as the subgradients of ϕ at x. The approximate subdifferential is defined, for ϵ ≥ 0, by

∂ϵϕ(x) := {s ∈ R n : ϕ(y) ≥ ϕ(x) + ⟨s, y -x⟩ -ϵ ∀ y ∈ R n } .
Let f : O → R be a locally Lipschitz continuous function. Then the generalized directional derivative defined by

f • (x; d) := lim sup x ′ → x, τ ↓ 0 f(x ′ + τ d) -f(x ′ ) τ is finite for all x ∈ O in every direction d ∈ R n [4, Prop. 2.1.1(a)
]. Such a mathematical concept permits us to define the Clarke subdifferential of f at x ∈ O,

∂ C f(x) := {g ∈ R n : ⟨g, d⟩ ≤ f • (x; d) for all d ∈ R n }, (2.2)
which is a nonempty, convex, and compact subset of R n [4, Prop. 2.1.2(a)] satisfying

f • (x; d) = max g∈∂ C f(x)
⟨g, d⟩.

The elements of ∂ C f(x) are referred to as generalized (or Clarke) subgradients, as they are the usual subgradients, i.e.,

∂ C f = ∂f, when f is convex [4, Prop. 2.2.7]. Furthermore, when f is continuously differentiable, ∂ C f(x) reduces to the singleton {∇f(x)}. An alternative representation, in finite dimensions, of ∂ C f(x) is (see [4, Thm. 2.5.1]) ∂ C f(x) := conv lim ι→∞ ∇f(xι), xι → x, f differentiable at xι .
A fundamental result, often evoked in this work, is the following one [4, Prop. 2.1.2]: the mapping

∂ C f is locally bounded in the interior of Dom(f) := {x ∈ R n : f(x) < ∞}. As a result, the image ∂ C f(X) of every bounded set X ⊂ O (⊂ Dom(f)) is bounded in R n .
Useful calculus rules of subdifferentials are listed in Proposition 2.1 below and rely on the concept of regularity.

A locally Lipschitz continuous function f : O → R is subdifferentially regular (or simply regular) at x ∈ O if for every d ∈ R n the ordinary directional derivative at x exists and coincides with the generalized one:

f ′ (x; d) = f • (x; d) ∀ d ∈ R n .
It holds that smooth functions, as well as convex ones, are regular at every point in the interior of their domains. Moreover, a finite linear combination (by nonnegative scalars) of regular functions at x is regular [4, Prop. The last item can be strengthened when more structure is present, such as in the case of weakly convex functions (see Eq. (2.3) below).

Weakly convex functions: definition and main properties Definition 2.2 (Def. 4.2 [START_REF] Vial | Strong and weak convexity of sets and functions[END_REF]). A function f : O → R is said to be (locally) weakly convex on O if, on some neighbourhood

V x ′ ⊂ O of each x ′ ∈ O, there exists µ x ′ ≥ 0 such that, for all µ ≥ µ x ′ ϕ(x) := f(x) + µ 2 ∥x∥ 2 is finite and convex on V x ′ .
Furthermore, f is said to be weakly convex in the global sense on O if the above property holds for

V x ′ = O and µ x ′ = μ ≥ 0 regardless of x ′ ∈ O. □
Clearly, a convex function on O is weakly convex in the global sense on O: it suffices to take µ x ′ = 0 and

V x ′ = O for all x ′ ∈ O.
f(x) = max t∈T ft(x).
in which the functions ft are of differentiability class C 2 on V x ′ and the index set T is a compact space such that ft(x), ∇ft(x), and ∇ 2 ft(x) depend continuously not just on x ∈ V x ′ but jointly on (t, x) ∈ T × V x ′ . □

In particular, if f is given by f(x) = max{f 1 (x), . . . , fm(x)} and all functions f 1 , ..., fm are of class C 2 , then f is Lower-C 2 /weakly convex. Furthermore, the functions of item iv) are also weakly convex, since they are Lower-C 2 by definition.

An important property of LC 2 /weakly convex functions is regularity [28, Thm. 1]: for every x ∈ O, the equality f ′ (x; d) = f • (x; d) holds in every direction d ∈ R n . Furthermore, Theorem 7.3 in [START_REF] Rockafellar | Generalized subgradients in mathematical programming[END_REF] gives the following characterization of the Clarke subdifferential of f at x ∈ O: for I(x) = arg max t∈T ft(x),

∂ C f(x) = conv {∇xft(x) : t ∈ I(x)} for all x ∈ O. (2.3)
When constrained to a compact convex set X ̸ = ∅, we can say more about weakly convex functions. Indeed, the local property in Definition 2.2 globally extends to the whole X, and we have the following key result (whose proof can be found in the Appendix A).

Proposition 2.4. Let f : O → R n be a weakly convex function, and X ⊂ O a compact and convex set.

Then there exist a real number

µ f ≥ 0 and an open convex set O ′ satisfying X ⊂ O ′ ⊂ O such that, for all µ ≥ µ f : i) the function ϕ(x) := f(x) + µ 2 ∥x∥ 2 is convex on O ′ and ∂ϕ(x) = ∂ C f(x) + µx for all x ∈ O ′ ;
ii) for all s f ∈ ∂ C f (x) with x ∈ O ′ , the following inequality holds

f(y) ≥ f(x) + ⟨s f , y -x⟩ - µ 2 ∥y -x∥ 2 ∀y ∈ X. (2.4)
Concerning the setting of this work where X is compact, the appealing DC decomposition f(x) = ϕ(x) -µ f 2 ∥x∥ 2 is, unfortunately, unavailable: the threshold µ f in Proposition 2.4 is in general unknown. This fact precludes the application of DC techniques to optimization problems featuring general Lower-C 2 /weakly convex functions. Interested readers are referred to [START_REF] Van Ackooij | Nonsmooth and nonconvex optimization via approximate differenceof-convex decompositions[END_REF] for a strategy that uses approximated DC decompositions based on item i) of Proposition 2.4.

Necessary optimality conditions and problem reformulation

Let f, c : O → R be given by (1.1b). We highlight that some properties of their components can be transferred to these functions. (To ease the presentation, let us focus only on f (x) = f 1 (x) -f 2 (x), as the same conclusions hold for c.) For instance, f is locally Lipschitz continuous because f 1 and f 2 are so. Furthermore, as f 1 is convex and f 2 is (locally) weakly convex, they are both directional differentiable and these properties extend to f as well: for every x ∈ O, the directionally derivative of f is finite in every direction d ∈ R n as result of the following relation:

f ′ 1 (x; d) -f ′ 2 (x; d) = lim τ ↓0 f 1 (x + τ d) -f 1 (x) τ -lim τ ↓0 f 2 (x + τ d) -f 2 (x) τ = lim τ ↓0 f 1 (x + τ d) -f 1 (x) τ - f 2 (x + τ d) -f 2 (x) τ = lim τ ↓0 f (x + τ d) -f (x) τ = f ′ (x; d).
However, the important regularity condition of both f 1 and f 2 does not extend to f as a mere fact that the latter is not a linear combination with nonnegative coefficients of the two former functions (see Proposition 2.1.i)). Hence, we cannot expect to have equality in the following inclusion

∂ C f (x) ⊂ ∂f 1 (x) -∂ C f 2 (x)
unless one of the component functions is smooth at x. (Here we have used the fact that ∂ C f 1 = ∂f 1 due to convexity of f 1 .) Such an inclusion impacts stationary concepts as we will now discuss. Let us first consider the convexly-constrained problem

min x∈X f (x), with f (x) = f 1 (x) -f 2 (x). (3.1)
A point x ∈ X is said to be directional (d)-stationary for this problem if f ′ (x; d) ≥ 0 for all d ∈ T X (x). Here, T X (x) is the Bouligand tangent cone to X at a point x ∈ X. The following result generalizes [START_REF] Pang | Computing B-stationary points of nonsmooth DC programs[END_REF]Prop. 5], where a specific case of problem (3.1) with f 2 (x) = max{ψ 1 (x), . . . , ψm(x)} and convex ψ 1 , . . . , ψm is considered.

Proposition 3.1. A point x ∈ X is d-stationary of problem (3.1) if, and only if,

x ∈ arg min x∈X f 1 (x) -[f 2 (x) + ⟨s f 2 , x -x⟩] ∀ s f 2 ∈ ∂ C f 2 (x).
Proof. Observe that T X (x) = cl {d ∈ R n : d = t(x -x), x ∈ X, t ∈ R + } due to convexity of X. Therefore, the definition of d-stationarity can be equivalently written as f ′ (x; x -x) ≥ 0 for all x ∈ X. Recall that f 2 is (locally) weakly convex and thus regular:

f ′ 2 (x; x -x) = f • 2 (x; x -x), which implies that f ′ 2 (x; x -x) = max s f 2 ∈∂ C f 2 (x) ⟨s f 2 , x -x⟩. Hence, f ′ (x; x -x) = f ′ 1 (x; x -x) -f ′ 2 (x; x -x) ≥ 0 ∀ x ∈ X ⇔ f ′ 1 (x; x -x) ≥ f ′ 2 (x; x -x) ∀ x ∈ X ⇔ f ′ 1 (x; x -x) ≥ ⟨s f 2 , x -x⟩ ∀ s f 2 ∈ ∂ C f 2 (x), ∀ x ∈ X ⇔ f ′ 1 (x; x -x) -⟨s f 2 , x -x⟩ ≥ 0 ∀ s f 2 ∈ ∂ C f 2 (x), ∀ x ∈ X ⇔ x ∈ arg min x∈X f 1 (x) -⟨s f 2 , x -x⟩ ∀ s f 2 ∈ ∂ C f 2 (x). A point x ∈ X is said to be Clarke-stationary of problem (3.1) if 0 ∈ ∂ C f (x) + N X (x). (3.2)
Furthermore, by following the lead of DC programming (see for instance [6,

§3.1]), x ∈ X is said to be a critical point if 0 ∈ ∂f 1 (x) -∂ C f 2 (x) + N X (x). (3.3) It is not difficult to see that this inclusion means that x ∈ arg min x∈X f 1 (x) -[f 2 (x) + ⟨s f 2 , x -x⟩] for some s f 2 ∈ ∂ C f 2 (x).
Note that the concept of criticality is weaker than that of Clarke-stationarity, which in turn is weaker than dstationarity (because f ′ (•; d) ≤ f • (•; d) for all d ∈ R n ). However, criticality and Clarke-stationarity coincide when at least one component function is smooth (in which case f is regular). Furthermore, we can see from Proposition 3.1 and the above alternative characterization of criticality that the three concepts coincide when f 2 is continuously differentiable at x.

For the more general problem (1.1), x ∈ X is said to be a Bouligand (B)-stationary point of (1.1) if f ′ (x; d) ≥ 0 for all d ∈ T X c (x), with X c the feasible set of (1.1). B-stationary boils down to d-stationarity if the nonconvex constraint c(x) ≤ 0 is absent or if the considered point satisfies strictly the nonlinear constraint (in the latter case, T X c (x) = T X (x)). Necessary and sufficient conditions for B-stationarity are given in [START_REF] Pang | Computing B-stationary points of nonsmooth DC programs[END_REF]Prop. 4] for the case of DC-constrained DC problems. The next result deals with a more general case: we assume that only c 2 is convex, while f 2 remains a weakly-convex function.

Proposition 3.2. In addition to our assumptions on problem (1.1), let c 2 : O → R be a convex function and x ∈ X c := {x ∈ X : c(x) ≤ 0} such that c(x) = 0. Moreover, assume that the following constraint qualification (CQ) holds

cl d ∈ T X (x) : c ′ (x; d) < 0 = d ∈ T X (x) : c ′ (x; d) ≤ 0 . (3.4)
Then, x is a B-stationary point of problem (1.1) if and only if x solves the convex problems

min x∈X f 1 (x) -[f 2 (x) + ⟨s f 2 , x -x⟩] s.t. c 1 (x) -[c 2 (x) + ⟨sc 2 , x -x⟩] ≤ 0 ∀ s f 2 ∈ ∂ C f 2 (x), ∀ sc 2 ∈ ∂c 2 (x). (3.5) Proof. Denote Ȳ (x) = {x ∈ X : c 1 (x) ≤ c 2 (x) + c ′ 2 (x; x -x)}. As the CQ (3.4) holds, Proposition 2.1 of [41] ensures that T X c (x) = T Ȳ (x) (x) = cl d ∈ R n : d = t(x -x), x ∈ Ȳ (x), t ∈ R + . Thus, the B-stationary definition is equivalent to f ′ (x; x -x) = f ′ 1 (x; x -x) -f ′ 2 (x; x -x) ≥ 0 ∀ x ∈ Ȳ (x) ⇔ f ′ 1 (x; x -x) ≥ ⟨s f 2 , x -x⟩ ∀ s f 2 ∈ ∂ C f 2 (x), ∀ x ∈ Ȳ (x). (3.6)
The stated result follows upon establishing the equivalence between (3.6) and (3.5).

[(3.6) ⇒ (3.5)]. Suppose (3.6) holds and let sc 2 ∈ ∂c 2 (x) be arbitrary. As

Y (sc 2 ) := {x ∈ X : c 1 (x) ≤ c 2 (x) + ⟨sc 2 , x -x⟩} ⊂ Ȳ (x) due to convexity of c 2 , we conclude that f ′ 1 (x; x -x) ≥ ⟨s f 2 , x -x⟩ for all s f 2 ∈ ∂ C f 2 (x)
and all x ∈ Y (sc 2 ). Convexity of the latter set implies that x minimizes f 1 (x) -⟨s f 2 , x -x⟩ over Y (sc 2 ) for all s f 2 ∈ ∂ C f 2 (x). Thus, condition (3.5) holds because sc 2 ∈ ∂c 2 (x) was taken arbitrarily.

[(3.5) ⇒ (3.6)]. To show the reverse implication, we proceed with a proof by contrapositive. Suppose that there exist s

′ f 2 ∈ ∂ C f 2 (x) and x ′ ∈ Ȳ (x) such that f ′ 1 (x; x ′ -x) < ⟨s ′ f 2 , x ′ -x⟩ (and hence x ′ ̸ = x), i.e., (3.6) does not hold. Let s ′ c 2 ∈ ∂c 2 (x) be such that c ′ 2 (x; x ′ -x) = ⟨s ′ c 2 , x ′ -x⟩. Therefore, x ′ is feasible for the convex problem min x∈X f 1 (x) -[f 2 (x) + ⟨s ′ f 2 , x -x⟩] s.t. c 1 (x) -[c 2 (x) + ⟨s ′ c 2 , x -x⟩] ≤ 0.
Together with our assumption

f ′ 1 (x; x ′ -x) < ⟨s ′ f 2 , x ′ -x⟩, we have that d = x ′ -
x is a feasible descent direction for the above problem, and thus x can not be one of its solution. Hence, x does not satisfy (3.5). The proof is thus complete.

Note that convexity of c 2 plays an important role in the above proposition. Indeed, if c 2 is not convex, then the set {x ∈ X : c 1 (x) ≤ c 2 (x) + ⟨sc 2 , x -x⟩} is not necessarily a subset of X c and when solving the linearized subproblem (3.5) we may get a point that is infeasible for the original problem (1.1).

Example 3.3. Let f 1 = x, f 2 = 0, c 1 = 0, c 2 = x 3
3 and X = [-2, 2]. We are not in the framework of Proposition 3.2, since c 2 is not convex on [-2, 2], but weakly convex (with mudulus µ = 4). At x = 0, which globally solves (1.1), the convex problem (3.5) becomes min x∈ [-2,2] x (because we have dropped the trivial constraint 0 ≤ 0), and thus does not provide a feasible point for the original problem. □ However, if the mudulus µ is known for the weakly convex function c 2 , adding a quadratic term in the constraint of the convex problem (3.5) makes the corresponding set feasible for the original problem. Moreover, Proposition 3.2 is generalized in case of weakly convex c 2 .

Corollary 3.4. Let c 2 : O → R be a weakly convex function and µc ≥ 0 be a real number from Proposition 2.4 corresponding to c 2 . Moreover, assume that the CQ (3.4) holds. Then, x is a B-stationary point of problem (1.1) if and only if, for any given µ ≥ µc, x solves the convex problems

min x∈X f 1 (x) -[f 2 (x) + ⟨s f 2 , x -x⟩] s.t. c 1 (x) -[c 2 (x) + ⟨sc 2 , x -x⟩] + µ 2 ∥x -x∥ 2 ≤ 0 ∀ s f 2 ∈ ∂ C f 2 (x), ∀ sc 2 ∈ ∂ C c 2 (x). (3.7)
Proof. Consider the convex functions c1 (x) = c 1 (x) + µ 2 ∥x∥ 2 and c2 (x) = c 2 (x) + µ 2 ∥x∥ 2 with µ ≥ µc. The result follows from Proposition 3.2 by using instead the DC decomposition c(x) = c1 (x) -c2 (x) and by noting that, for an arbitrary sc 2 ∈ ∂c 2 (x), we obtain

c1 (x) -[c 2 (x) + ⟨sc 2 , x -x⟩] = c 1 (x) -[c 2 (x) + ⟨sc 2 , x -x⟩] + µ 2 ∥x -x∥ 2 with sc 2 = sc 2 -µx, sc 2 ∈ ∂ C c 2 (x).
Except for some particular cases, checking B-stationarity numerically is out of reach. Therefore, weaker stationarity concepts need to come into play: x ∈ X is said to be Clarke-stationary for (1.1) if there exists a Lagrange multiplier λ such that

0 ∈ ∂ C f (x) + λ ∂ C c(x) + N X (x) c(x) ≤ 0, λ c(x) = 0, λ ≥ 0, x ∈ X. (3.8)
Analogously, x is a critical point of (1.1) if there exists a Lagrange multiplier λ such that

0 ∈ ∂f 1 (x) -∂ C f 2 (x) + λ [∂c 1 (x) -∂ C c 2 (x)] + N X (x) c(x) ≤ 0, λ c(x) = 0, λ ≥ 0, x ∈ X. (3.9)
Observe that if f 1 or f 2 and c 1 or c 2 are smooth, then criticality boils down to Clarke stationarity. Next, we revisit the proximal bundle method of [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] and extend it to the more general setting of problem (1.1).

To this end, the method must be modified, and its convergence analysis must be done anew.

Problem reformulation via improvement function

Nonsmooth and noncovex constraints in optimization problems are in general dealt with numerically via exact penalization [START_REF] Le Thi | DC programming and DCA for general DC programs[END_REF][START_REF] Lipp | Variations and extension of the convex-concave procedure[END_REF][START_REF] Strekalovsky | A local search method for optimisation problem with d.c. inequality constraints[END_REF][START_REF] Jara-Moroni | A study of the difference-of-convex approach for solving linear programs with complementarity constraints[END_REF], linearization of certain components [START_REF] Pang | Computing B-stationary points of nonsmooth DC programs[END_REF][START_REF] Van Ackooij | Nonsmooth and nonconvex optimization via approximate difference-of-convex decompositions[END_REF], and improvement functions [START_REF] Sagastizábal | An infeasible bundle method for nonsmooth convex constrained optimization without a penalty function or a filter[END_REF][START_REF] Apkarian | Mixed H 2 /H∞ Control via Nonsmooth Optimization[END_REF][START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF]. The latter has a recognized good practical performance, does not require the additional assumptions normally assumed in exact penalization methods, and employs parameters that are simple to set. For these reasons, we handle problem (1.1) via the improvement function (3.10b) Observe that if x is a global solution of (1.1), then H(x; x) ≥ 0 for all x ∈ X and H(x; x) = 0.

H : O × O → R given by H(x; y) = max f (x) -τ f (y), c(x) -τc(y) , (3.10a
Improvement functions (also known as progress functions) have been considered within bundle methods in [START_REF] Sagastizábal | An infeasible bundle method for nonsmooth convex constrained optimization without a penalty function or a filter[END_REF][START_REF] Van Ackooij | Constrained bundle methods for upper inexact oracles with application to joint chance constrained energy problems[END_REF][START_REF] Van Ackooij | Level bundle methods for constrained convex optimization with various oracles[END_REF] for convex problems, in [START_REF] Apkarian | Mixed H 2 /H∞ Control via Nonsmooth Optimization[END_REF] for a class of (nonconvex) optimal control problems, and in [START_REF] Montonen | Bundle-based descent method for nonsmooth multiobjective DC optimization with inequality constraints[END_REF][START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] for DC-constrained DC programs. In what follows we exploit some relevant mathematical properties of (3.10) and its link to the original problem (1.1). To this end, we need to consider necessary conditions for a point x to be a local solution of the reformulated problem

min x∈X H(x; x) . (3.11) As the second argument of H is fixed, it follows from (3.2) that x ∈ X is a Clarke-stationary point of (3.11) if 0 ∈ ∂ C 1 H(x; x) + N X (x), (3.12
) where ∂ C 1 H stands for the generalized subdifferential of H with respect to the first argument. Proposition 2.1 ii) yields

∂ C 1 H(x; x) ⊂    ∂ C c(x) if f (x) -τ f (x) < c(x) -τc(x) conv ∂ C f (x), ∂ C c(x) if f (x) -τ f (x) = c(x) -τc(x) ∂ C f (x) if f (x) -τ f (x) > c(x) -τc(x
). Since we do not work with generalized subgradients of either f or c, but only with subgradients of the functions yielding their CwC decompositions (1.1b), we must consider a weaker stationary definition: we say that x ∈ X is a critical point of the composite problem (3.11)

with CwC decompositions (1.1b) if 0 ∈ N X (x) +    ∂c 1 (x) -∂ C c 2 (x) if f (x) -τ f (x) < c(x) -τc(x) conv ∂f 1 (x) -∂ C f 2 (x), ∂c 1 (x) -∂ C c 2 (x) if f (x) -τ f (x) = c(x) -τc(x) ∂f 1 (x) -∂ C f 2 (x) if f (x) -τ f (x) > c(x) -τc(x). (3.13)
Note that if both f and c are regular, then the above condition coincides with that of Clarke stationarity: recall (3.2), Proposition 2.1 i), and observe that the set defined by the expressions in the curly brackets above is nothing but ∂ C 1 H(x; x). The following result, inspired by both [1, Lemma 5.1] that deals with the (stronger) Clarke stationarity and [49, Thm. 2] that works with the (weaker) criticality definition from DC programming, links condition (3.13) with criticality of the original problem. Theorem 3.5. Let x ∈ X be a point satisfying condition (3.13). Then, the following hold: i) If c(x) > 0, then x is a critical point (in the sense of (3.3)) of the optimization problem

min x∈X c 1 (x) -c 2 (x). (3.14)
ii) If c(x) = 0 and x is not a critical point of (3.14), then x satisfies (3.9) for some λ > 0.

iii) If c(x) < 0, then x satisfies (3.9) with λ = 0.

Proof. With the τ function defined in (3.10b), note that

f (x) -τ f (x) -[c(x) -τc(x)] =    -[ρ + (1 -σ)]c(x) < 0 if c(x) > 0, 0 if c(x) = 0, -c(x) > 0 if c(x) < 0.
Hence,

c(x) > 0 ⇔ f (x) -τ f (x) < c(x) -τc(x), c(x) = 0 ⇔ f (x) -τ f (x) = c(x) -τc(x), c(x) < 0 ⇔ f (x) -τ f (x) > c(x) -τc(x),
and items i) and iii) follow directly from (3.13). To show item ii), recall that c(x) = 0 and condition (3.13) ensures the existence of λ

∈ [0, 1] such that 0 ∈ λ[∂f 1 (x) -∂ C f 2 (x)] + (1 -λ)[∂c 1 (x) -∂ C c 2 (x)] + N X (x).
By assumption, x is not a critical point of (3.14). Then 0 ̸ ∈ ∂c 1 (x)-∂ C c 2 (x)+N X (x), implying that λ above must be strictly positive. Dividing the displayed inclusion by λ > 0 we obtain the criticality condition (3.9) with λ = (1 -λ)/λ > 0.

At item ii) above, the assumption that x is not a critical point of (3.14) can be seen as a constraint qualification, which turns out to be more restrictive than (3.4). Indeed, the latter excludes d-stationary points of (3.14), but not necessarily critical ones. The following example gives a critical point x of (3.14) that satisfies (3.4) but not the criticality condition (3.9) for the nonlinearly-constrained problem (1.1).

Example 3.6. Take c 1 (x) = max{x, 2x}, c 2 (x) = max{2x, 4x}, X = [-2, 2] and x = 0. Then T X (x) = R, N X (x) = {0}, and x is a critical point of (3.14)

because 0 ∈ ∂c 1 (x) -∂c 2 (x) = [1, 2] -[2, 4] = [-3, 0]. Furthermore, note that c ′ 1 (x; d) = max s∈[1, 2] sd = 2d if d ≥ 0 d if d ≤ 0 and c ′ 2 (x; d) = max s∈[2, 4] sd = 4d if d ≥ 0 2d if d ≤ 0 , thus c ′ (x; d) = min{-d, -2d}. We conclude that {d ∈ T X (x) : c ′ (x; d) < 0} = R + , whereas {d ∈ T X (x) : c ′ (x; d) ≤ 0} = R + ∪ {0},
showing that x = 0 satisfies the CQ (3.4). However, if we take f 1 (x) = 0 and f 2 (x) = -1 2 x 2 + x, the following system does not have a solution:

0 ∈ ∂f 1 (0) -∂ C f 2 (0) + λ[∂c 1 (0) -∂c 2 (0)] λ ≥ 0 ≡ 0 ∈ -1 + λ[-3, 0] λ ≥ 0 ≡ 0 ∈ [-3 λ -1, -1] λ ≥ 0,
i.e., x does not satisfy (3.9). Figure 1 illustrates the objective and constraint function in this example: it is clear that x is indeed a global maximizer of f (x) under the constraints x ∈ X and c(x) ≤ 0. □ This example shows that, at item ii) of Theorem 3.5, we cannot replace the assumption that x is not a critical point of (3.14) with the CQ (3.4).

The DC setting.

In the DC setting, functions f 2 and c 2 are convex and the improvement function (3.10) is DC. Indeed, for x fixed, we can write Proof. Let x ∈ X be a point satisfying (3.13). Let us first observe that since f 2 , c 2 are convex and thus regular, we have ∂G(x) = ∂f 2 (x) + ∂c 2 (x). A similar observation can be made concerning the computation for F . Our analysis splits into three possible cases.

H(x; x) = F (x; x) -G(x), with F (x; x) = max f 1 (x) + c 2 (x) -τ f (x), f 2 (x) + c 1 (x) -τc(x) , G(x) = f 2 (x) + c 2 (x). ( 3
a) f (x) -τ f (x) < c(x) -τc(x). It follows from (3.13) that 0 ∈ N X (x) + ∂c 1 (x) -∂c 2 (x) ⊂ N X (x) + ∂f 2 (x) + ∂c 1 (x) -[∂f 2 (x) + ∂c 2 (x)].
We claim that this inclusion implies (3.16). To see that, observe that the above inequality implies

f 1 (x) + c 2 (x) -τ f (x) < f 2 (x) + c 1 (x) -τc(x), which in turn gives ∂ 1 F (x; x) = ∂f 2 (x) + ∂c 1 (x). Therefore, the right-hand-side of the above inclusion is (3.16). b) f (x) -τ f (x) = c(x) -τc(x). It follows from (3.13) that there exists λ ∈ [0, 1] such that 0 ∈ N X (x) + λ[∂f 1 (x) -∂f 2 (x)] + (1 -λ)[∂c 1 (x) -∂c 2 (x)] = N X (x) + λ[∂f 1 (x) + ∂c 2 (x)] -λ∂f 2 (x) + (1 -λ)∂c 1 (x) -∂c 2 (x) ⊂ N X (x) + λ[∂f 1 (x) + ∂c 2 (x)] + (1 -λ)[∂f 2 (x) + ∂c 1 (x)] -[∂f 2 (x)+∂c 2 (x)] ⊂ N X (x) + ∂ 1 F (x; x) -∂G(x). c) f (x) -τ f (x) > c(x) -τc(x). Again, (3.13) gives 0 ∈ N X (x) + ∂f 1 (x) -∂f 2 (x) ⊂ N X (x) + ∂f 1 (x) + ∂c 2 (x) -[∂f 2 (x) + ∂c 2 (x)].
The proof is complete because in this case ∂ 1 F (x; x) = ∂f 1 (x) + ∂c 2 (x) due to the fact that

f 1 (x) + c 2 (x) -τ f (x) > f 2 (x) + c 1 (x) -τc(x).
Remark 3.8. In the DC setting, the three concepts of criticality (3.12), (3.13), and (3.16) are equivalent when f 2 and c 2 are continuously differentiable at x. Indeed, in this case f and c are regular at x and (3.12) coincides with (3.13) (regardless of convexity of f 2 and c 2 ). Theorem 2 in [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] ensures that, under these assumptions, (3.16) is equivalent to (3.12). □

The following example shows that (3.16) does not necessarily imply (3.13) in the nondifferentiable DC case.

Example 3.9.

Let X = [-1, 1], f 1 = 2x, f 2 = |x|, c 1 = 4x and c 2 = 2|x|. At x = 0, f (x) = c(x) = 0 and thus, τ f (x) = τc(x) = 0 due to (3.10b). Furthermore, we have that ∂f 1 (0) = {2}, ∂f 2 (0) = [-1, 1], ∂c 1 (0) = {4}, and ∂c 2 (0) = [-2, 2]. As a result, ∂f 1 (0) -∂f 2 (0) = [1, 3], ∂c 1 (0) -∂c 2 (0) = [2, 6], ∂f 1 (0) + ∂c 2 (0) = [0, 4], and ∂f 2 (0) + ∂c 1 (0) = [3, 5]. As in N X (0) = {0}, we conclude that 0 ̸ ∈ [1, 6] = conv{∂f 1 (0) -∂f 2 (0), ∂c 1 (0) -∂c 2 (0)} + N X (0), whereas conv{∂f 1 (0) + ∂c 2 (0), ∂c 1 (0) + ∂f 2 (0)} + N X (0) = [0, 5] and ∂f 2 (0) + ∂c 2 (0) = [-3, 3],
showing that (3.16) is satisfied but not (3.13). □

The paper [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] proposes a bundle method for DC-constrained DC programs employing the DC decomposition H = F -G above. Once a critical point satisfying (3.16) is computed, the link with criticality of the original problem is adequate if f 2 and c 2 are continuously differentiable at x. In the next section we modify that method to compute a point satisfying the stronger criticality condition (3.13). As a result, the link with criticality of the original problem is nicely established by Theorem 3.5 without any additional assumption. In fact, f 2 and c 2 need not even be convex, but weakly convex on some neighbourhood of each x ∈ O. We, therefore, strengthen the analysis provided in [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] even though significantly fewer assumptions are required: [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] works in the DC configuration, whereas here, we deal with the more general CwC structure. These improvements, together with the optimality conditions presented above, feature the main contributions of this work.

Proximal bundle method with improvement function

This section extends the proximal bundle method of [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] for computing a critical point of problem (1.1). The main tool in our analysis is the improvement function H given in (3.10). In the DC setting, the algorithm of [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] works with the explicit DC decomposition (3.15) of H and computes a point x ∈ X satisfying the classic criticality condition in DC programming (3.16). In this section we do not decompose H and consider the milder assumption that f 2 and c 2 are weakly convex and target the stronger criticality condition (3.13).

The method's main ingredients: model, subproblem, and descent test

The algorithm requires four oracles (black-boxes) providing, for every given x ∈ X, i ∈ {1, 2}, the function values f i (x), c i (x), arbitrary subgradients s f 1 ∈ ∂f 1 (x), sc 1 ∈ ∂c 1 (x) (c.f., (2.1)) and arbitrary generalized subgradients

s f 2 ∈ ∂ C f 2 (x), sc 2 ∈ ∂ C c 2 (x) (c.f. (2.
2)). We do not impose any assumption on these (generalized) subgradients, as they are assumed to be computed by (external) oracles that do not accept any intervention from the algorithm. (This is particularly useful in industrial applications where companies do not want or cannot share information on the underlying functions with optimizers.) At iteration k ∈ N, given a trial point x k ∈ X, we construct a linearization of every component (here

s k f i , s k c i , i ∈ {1, 2}
, denote the respective -generalized -subgradients at x k ):

f k 1 (x) := f 1 (x k ) + ⟨s k f 1 , x -x k ⟩ (4.1a) f k 2 (x) := f 2 (x k ) + ⟨s k f 2 , x -x k ⟩ (4.1b) ck 1 (x) := c 1 (x k ) + ⟨s k c 1 , x -x k ⟩ (4.1c) ck 2 (x) := c 2 (x k ) + ⟨s k c 2 , x -x k ⟩. (4.1d)
Due to convexity of f 1 and c 1 , we have the following inequalities

f k 1 (x) ≤ f 1 (x) and ck 1 (x) ≤ c 1 (x) for all x ∈ R n . (4.2)
Since X is compact and components f 2 and c 2 are assumed to be only weakly convex on some neighbourhood of each x ∈ O, we have weaker inequalities for these functions. Let O ′ ⊂ O be an open convex set and µ f , µc real numbers ensured by Proposition 2.4. As x k ∈ X ⊂ O ′ , the following inequalities are due to Proposition 2.4, item ii)

f k 2 (x) ≤ f 2 (x) + μ 2 ∥x -x k ∥ 2 and ck 2 (x) ≤ c 2 (x) + μ 2 ∥x -x k ∥ 2 for all x ∈ X, (4.3)
where μ := max{µ f , µc} > 0. Observe that the threshold μ is in general unknown, and the inequalities in (4.3) are only supposed to hold for x in X, in contrast with the (subgradient) inequalities in (4.2). Let B k f and B k c be two index sets gathering the bundle of information (function values and subgradients) given by the oracles. In general, B k f , B k c ⊂ {0, . . . , k} but other possibilities exist making it possible to design a limited-memory method (see Remark 5.5 below). These index sets are useful to define the following individual cutting-plane models for the convex functions f 1 and c 1 :

f k 1 (x) := max j∈B k f f j 1 (x) ≤ f 1 (x) for all x ∈ R n čk 1 (x) := max j∈B k c cj 1 (x) ≤ c 1 (x) for all x ∈ R n .
Furthermore, let ℓ k ∈ {0, . . . , k} be the iteration index of the best candidate solution (stability center, in the parlance of bundle methods) among the trial points {x 0 , . . . , x k }: whenever a better candidate solution x k+1 is computed by the algorithm, at a so-called serious step, such a point becomes the new stability center and the counter ℓ is increased by one: for κ ∈ (0, 1 2 ), we declare a serious step and let

ℓ k+1 := k + 1 if x k+1 ̸ = x ℓ k and the inequality H(x k+1 ; x ℓ k ) ≤ H(x ℓ k ; x ℓ k ) - κ 2 ∥x k+1 -x ℓ k ∥ 2 (4.5)
holds, and declare a null step and let ℓ k+1 := ℓ k otherwise. Since the descent test is independent of the model, the following result from [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] also holds in our framework.

Lemma 4.1 (Lemma 1 in [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF]). Let x ℓ k ∈ X be the stability center at iteration k. Then H(x ℓ k ; x ℓ k ) ≥ 0 and if inequality (4.5) holds, we have that either

i) f (x k+1 ) ≤ f (x ℓ k ) -κ 2 ∥x k+1 -x ℓ k ∥ 2 and c(x k+1 ) ≤ 0 when c(x ℓ k ) ≤ 0; or ii) c(x k+1 ) ≤ c(x ℓ k ) -κ 2 ∥x k+1 -x ℓ k ∥ 2 when c(x ℓ k ) > 0.
The rationale of serious iterates is to ensure sufficient decrease on one component function of H(•; x ℓ k ) while maintaining feasibility for (1.1) once reached.

Having all these ingredients at our disposal, we can now define our convex model for the improvement function (3.10) at iteration k :

Ȟk (x; x ℓ k ) = max f k 1 (x) - f ℓ k 2 (x) -τ f (x ℓ k ), čk 1 (x) - cℓ k 2 (x) -τc(x ℓ k ) . (4.6)
(Even in the particular setting where f 2 and c 2 are convex functions, this model differs from the one employed in [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] and is crucial to obtain convergence results stronger than the ones in that paper.) Given a proxparameter µ k > 0 estimating the threshold μ in (4.3), the next iterate is the solution of the following strict convex subproblem

x k+1 = arg min x∈X Ȟk (x; x ℓ k ) + µ k 2 ∥x -x ℓ k ∥ 2 , (4.7)
which can be transformed into a QP (provided X is a polyhedron) by adding an extra variable r ∈ R

           min x, r r + µ k 2 ∥x -x ℓ k ∥ 2 s.t. f j 1 (x) - f ℓ k 2 (x) -r ≤ τ f (x ℓ k ) ∀ j ∈ B k f cj 1 (x) - cℓ k 2 (x) -r ≤ τc(x ℓ k ) ∀ j ∈ B k c x ∈ X, r ∈ R. (4.8)
The optimality condition for (4.7) gives

x k+1 = x ℓ k - 1 µ k [p k+1 + s k+1 X ], with p k+1 ∈ ∂ 1 Ȟk (x k+1 ; x ℓ k ) s k+1 X ∈ N X (x k+1 ). (4.9)
As usual in bundle methods, we may remove from the model the inactive linearizations to keep (4.8) small. To this end, we denote by Bk f ⊂ B k f and Bk c ⊂ B k c the index set of active linearizations in the QP subproblem (4.8), i.e.,

Bk f := j ∈ B k f : α j f > 0 and Bk c := j ∈ B k c : α j c > 0 (4.10)
where α j f ≥ 0, j ∈ B k f , denote the Lagrange multipliers associated with the first set of constraints and α j c ≥ 0, j ∈ B k c , the ones associated with the second family of constraints. We mention in passing that the index sets B k f and B k c can be kept bounded at the price of including artificial (aggregate) linearizations. We postpone this discussion to Remark 5.5, right after the analysis of null steps (the only place in the convergence analysis where bundle management plays a role.)

We can now present the following proximal bundle method algorithm for CwC-constrained CwC programs (1.1), which modifies [49, Alg. 1] in two ways. First, the convex model (4.6) for the improvement function is distinct. On the one hand, it is a key element to obtain the stronger criticality condition (3.13), and on the other hand, it leads to a simpler/smaller strongly convex QP (4.8) (more details can be found in Subsection 4.2 below). Second, Algorithm 1 employs an ad-hoc rule to update the proximal parameter µ k so that no pre-estimation of the underlying weakly-convex moduli is needed. The proposed rule employs the following value

ν k := 2 max f ℓ k 2 (x k+1 ) -f 2 (x k+1 ) ∥x k+1 -x ℓ k ∥ 2 , cℓ k 2 (x k+1 ) -c 2 (x k+1 ) ∥x k+1 -x ℓ k ∥ 2 , 0 . (4.11)
Algorithm 1 Proximal Bundle Method for CwC-constrained CwC programs -CwC-PBM

Step 0 (Initialization) Let ), and Tol ≥ 0 be given. Call the oracles to compute f i (x 0 ), c i (x 0 ), and (generalized) subgradients s 0 fi ,

x 0 ∈ X, κ ∈ (0, 1 2 ), κ ≤ µ 0 , ρ ≥ 0, δ ∈ [0, 1
s 0 ci , i = 1, 2. Define k := ℓ k = 0 and B 0 f = B 0 c := {0}.
Step 1 (Trial point) Compute x k+1 by solving the QP (4.8).

Step 2 (Stopping test) If ∥x k+1 -x ℓ k ∥ ≤ Tol, then stop and return x ℓ k .

Step 3 (Oracles call) Compute f i (x k+1 ), c i (x k+1 ), and subgradients s k+1 fi , s k+1 ci , i = 1, 2.

Step 4 (Descent test)

ν k ≥ µ k -2κ, set µ k+1 = ν k + 1; otherwise µ k+1 = µ k .
Step 5 (Loop) Set k := k + 1 and go back to Step 1.

A drawback of the rule for updating the prox-parameter is that µ k only increases after a null step when the inequality ν k ≥ µ k -2κ is verified. As a result, µ k may never increase: this is, for instance, the case when f 2 and c 2 are convex (thus ν k = 0 for all k). The motivation for this rule is to eventually keep the prox-parameter fixed if the algorithm performs an infinite sequence of null steps after a last serious step (see Lemma 5.1). This is a condition necessary to prove Proposition 5.4 below. We care to mention that increasing µ k after a null step is a simple strategy that pays off in practice: it helps the algorithm to either stop or produce a new serious step, and thus accelerate the numerical performance.

The DC setting: a comparison with the earlier bundle method for DC probgrams

In the DC setting, both functions f 2 and c 2 are convex and the improvement function (3.10) is DC. The DC decomposition given in (3.15), with x replaced with x ℓ k , was exploited in the bundle method of [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] through the following model for the improvement function H(•; x ℓ k ) (see Eq. ( 18) therein):

max f k 1 (x) + čk 2 (x) -τ f (x ℓ k ), f k 2 (x) + čk 1 (x) -τc(x ℓ k ) -[ f ℓ k 2 (x) + cℓ k 2 (x)]. (4.12)
Differently from our model (4.6), the above gathers also cutting-planes for f 2 and c 2 and, although gathering more information, only the weaker criticality condition (3.16) is ensured by the method of [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF]. Hence, the proposed model (4.6) is more advantageous than (4.12) from both practical and theoretical point of view:

• the quadratic program (QP) issued by our model has only half of the linearizations, and is thus simpler to solve; • convexity of f 2 and c 2 are required in (4.12), but not in (4.6);

• both models (4.6) and (4.12) are iteratively updated to ensure that every cluster point x ∈ X of the sequence of stability centers satisfies a criticality condition. To show that such a point is also critical for (the DC counterpart of) (1.1), [49, Thm. 2] requires both f 2 and c 2 to be continuously differentiable at x. As we will show in Theorem 5.8 below, neither convexity nor differentiability of f 2 and c 2 are required to establish that x issued by Algorithm 1 is also critical for (1.1) in the sense of (3.9). Thus, Algorithm 1 strengthens the results of [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] even though significantly fewer assumptions are required.

Although the apparently small changes concerning [49, Alg. 1], the convergence analysis in that paper cannot be reused here. The reason is that the analysis in [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] strongly depends on the DC decomposition of the employed model for the improvement function. That reasoning is no longer valid for our new model, even if f 2 and c 2 were convex. Furthermore, our more general setting requires extra steps to cope with the weakly convex functions.

Convergence Analysis

The goal of this section is to show that every cluster point x of the sequence {x ℓ k } k ⊂ X generated by Algorithm 1 satisfies the necessary optimality condition (3.13). To this end, we first observe that the sequence of prox-parameters issued by Algorithm 1 is bounded.

Lemma 5.1. The value µmax := sup k∈N µ k is finite. Furthermore, if the algorithm produces an infinite sequence of null steps after a last serious step, then the prox-parameter becomes eventually constant.

Proof. Let μ := max{µ f 2 , µc 2 , µ 0 } > 0 be given, where µ f 2 and µc 2 are as in Proposition 2.4 for the weakly convex functions f 2 and c 2 , and µ 0 is the parameter given to the algorithm at initialization. Then, by taking y := x k+1 and x := x ℓ k in (2.4) it follows that

2 f ℓ k 2 (x k+1 ) -f 2 (x k+1 ) ∥x k+1 -x ℓ k ∥ 2 ≤ μ, and 2 cℓ k 2 (x k+1 ) -c 2 (x k+1 ) ∥x k+1 -x ℓ k ∥ 2 ≤ μ for all k with x k+1 ̸ = x ℓ k .
As a result, ν k ≤ μ for all k. Note that the prox-parameter is only increased after a null step such that ν k ≥ µ k -2κ. In this case, the rule employed in Step 4 of the algorithm sets µ k+1 = ν k + 1, which gives µ k+1 = ν k + 1 ≤ μ + 1. Since the algorithm does not increase the prox-parameter after a serious step or null step such that ν k < µ k -2κ, we conclude that µmax := sup k∈N µ k ≤ μ+1 is finite. Finally, note that the proxparameter is sharply increased after a null step such that ν k ≥ µ k -2κ:

µ k+1 = ν k +1 ≥ µ k -2κ+1 > µ k +δ because κ ∈ (0, 1 2 ), with δ = 1 2 -κ > 0.
As a result, if the algorithm produces an infinite sequence of null steps after a last serious step, then the inequality ν k < µ k -2κ will be satisfied for all k large enough and the prox-parameter will become constant (otherwise µ k would increase indefinitely, which contradicts that µmax is finite).

We now define the following function H : O × O → R, which is of key importance in our analysis:

H(x; y) := max f 1 (x) -[f 2 (y) + ⟨s f 2 , x -y⟩] -τ f (y), c 1 (x) -[c 2 (y) + ⟨sc 2 , x -y⟩] -τc(y) , (5.1)
with s f 2 ∈ ∂ C f 2 (y) and sc 2 ∈ ∂ C c 2 (y). As these subgradients are not specified, the above definition is ambiguous. However, when y is a point previously computed by the algorithm, say y = x j for j ≤ k, then

s j f 2 ∈ ∂ C f 2 (x j ) and s j c 2 ∈ ∂ C c 2 (x j
) are the subgradients provided by the oracles and ambiguity disappears:

H(x; x j ) := max f 1 (x) -f j 2 (x) -τ f (x j ), c 1 (x) -cj 2 (x) -τc(x j ) .
It follows from convexity of f 1 and c 1 that, for every y ∈ O fixed, the function H(•; y) is convex and satisfies H(•;

x ℓ k ) ≥ Ȟk (•; x ℓ k ) for all k. Furthermore, as ℓ k ∈ B k f ∩ B k c for all k, we have that f k i (x ℓ k ) = f i (x ℓ k ), čk i (x ℓ k ) = c i (x ℓ k ), i = 1, 2, and thus H(x ℓ k ; x ℓ k ) = Ȟk (x ℓ k ; x ℓ k ) = H(x ℓ k ; x ℓ k ). (5.2)
The following lemma is of particular interest in the remainder of this work.

Lemma 5.2. Suppose that x minimizes H(•; x) over X. Then, x satisfies the necessary optimality condition (3.13).

Proof. Convexity of H(•; x) in the first argument and assumption on x ∈ X imply that 0 ∈ ∂ 1 H(x; x)+N X (x).

The result follows by noting that, for some pair of generalized subgradients sf 2 ∈ ∂ C f 2 (x) and sc 2 ∈ ∂ C c 2 (x), the following set

∂ 1 H(x; x) =    ∂c 1 (x) -sc 2 if f (x) -τ f (x) < c(x) -τc(x) conv ∂f 1 (x) -sf 2 , ∂c 1 (x) -sc 2 if f (x) -τ f (x) = c(x) -τc(x) ∂f 1 (x) -sf 2 if f (x) -τ f (x) > c(x) -τc(x)
is contained in the one defined by the curly brackets in (3.13).

We begin the convergence analysis for the case Tol = 0 with the remark that the sequence of stability centers {x ℓ k } k has at least one cluster point, since it is contained in the compact set X. We split the analysis into three cases: the algorithm performs only finitely many steps; the algorithm performs infinitely many steps and the sequence {x ℓ k } k is either finite or infinite. Proposition 5.3 (Finitely many iterations). Assume that Algorithm 1 stops at iteration k with Tol = 0. Then, the last stability center x := x ℓ k = x k+1 satisfies condition (3.13).

Proof. It follows from the model's definition (4.6) and (5.1) that H(x; x ℓ k ) ≥ Ȟk (x; x ℓ k ) for all x ∈ O. Hence, as x ℓ k ∈ X we have that

H(x ℓ k ; x ℓ k ) ≥ min x∈X H(x; x ℓ k ) + µ k 2 ∥x -x ℓ k ∥ 2 ≥ min x∈X Ȟk (x; x ℓ k ) + µ k 2 ∥x -x ℓ k ∥ 2 = Ȟk (x k+1 ; x ℓ k ) + µ k 2 ∥x k+1 -x ℓ k ∥ 2 = Ȟk (x ℓ k ; x ℓ k ) = H(x ℓ k ; x ℓ k ),
where the first equality is due to (4.7), the second one follows by the fact that x k+1 = x ℓ k since the algorithm stops at iteration k with Tol = 0, and the last one is due to (5.2). Hence,

x ℓ k minimizes H(•; x ℓ k ) + µ k 2 ∥ • -x ℓ k ∥ 2
over X and the quadratic term vanishes in the corresponding optimality condition: 0 ∈ ∂ 1 H(x; x ℓ k ) + N X (x) and the stated result follows from Lemma 5.2.

If the algorithm performs finitely many serious steps and infinite number of null steps, the following result shows that the last stability center satisfies (3.13).

Proposition 5.4 (Finitely many serious steps). Suppose that Algorithm 1 with Tol = 0 does not stop but produces only finitely many serious steps. Then the last stability center x satisfies the condition (3.13).

Proof. Let ℓ ∈ N denote the last serious iteration, then x = x ℓ and note that, for all subsequent (null) iterations k > ℓ, ℓ k = ℓ and the linearizations f ℓ 2 and cℓ 2 are fixed in the model Ȟk (•; x), which is in this case a cutting-plane model for the convex function H(•; x). Here we take τ ℓ f = τ f (x ℓ ), τ ℓ c = τ ℓ c (x ℓ ), and function H(•; x) defined with the fixed linearizations f ℓ 2 and cℓ 2 , i.e., H(•;

x ℓ ) := max f 1 (•) -f ℓ 2 (•) -τ ℓ f , c 1 (•) -cℓ 2 (•) -τ ℓ c .
We highlight that the updating rule for µ k in Algorithm 1 ensures that the sequence {µ k } k>ℓ is nondecreasing and becomes constant at a certain value µ ′ ∈ (0, µmax] after finitely many steps k ′ > ℓ, as a consequence of Lemma 5.1. More precisely, the updating rule at Step 4(b) of Algorithm 1 ensures that

µ k = µ ′ and ν k + 2κ < µ ′ for all k > k ′ . (5.3)
Hence, from iteration k ′ on, Algorithm 1 becomes a cutting-plane procedure to compute the unique solution x of

min x∈X H(x; x) + µ ′ 2 ∥x -x∥ 2 .
(5.4)

As the algorithm keeps all the active linearizations in the bundles (Step 4(b)), standard arguments from the convex bundle methods' theory (see [START_REF] Correa | Convergence of some algorithms for convex minimization[END_REF]Prop. 4.3]) ensure that lim k→∞ x k = x and lim k→∞ [ Ȟk (x k+1 ; x) -H(x k+1 ; x)] = 0.

(The last inequality implies that the convex model asymptotically coincides with the function at the limit point.) We claim that x = x. To show that, let us assume the opposite, i.e., x ̸ = x, and arrive to a contradiction. In this case, for some δ > 0, we may find and index k 1 such that ∥x k+1 -x∥ 2 > δ for all k ≥ k 1 . We may furthermore find an index k 2 such that Ȟk (x k+1 ; x) -H(x k+1 ; x) ≥ -κ 2 δ for all k ≥ k 2 as the left-hand side vanishes. Therefore, for

k ′′ ≥ max{k 1 , k 2 , k ′ }, we have Ȟk (x k+1 ; x) -H(x k+1 ; x) ≥ - κ 2 ∥x k+1 -x∥ 2 ̸ = 0 for all k > k ′′ .
The following chain of inequalities holds at every iteration k > k ′′ :

Ȟk (x; x) ≥ Ȟk (x k+1 ; x) + µ ′ 2 ∥x k+1 -x∥ 2
(by (4.7) and ( 5.

3)) = [ Ȟk (x k+1 ; x) -H(x k+1 ; x)] + H(x k+1 ; x) + µ ′ 2 ∥x k+1 -x∥ 2 ≥ -κ 2 ∥x k+1 -x∥ 2 + H(x k+1 ; x) + µ ′ 2 ∥x k+1 -x∥ 2 ≥ -κ 2 ∥x k+1 -x∥ 2 + max f 1 (x k+1 ) -f ℓ 2 (x k+1 ) -τ ℓ f + µ ′ 2 ∥x k+1 -x∥ 2 c 1 (x k+1 ) -cℓ 2 (x k+1 ) -τ ℓ c + µ ′ 2 ∥x k+1 -x∥ 2
(by (5.1))

> -κ 2 ∥x k+1 -x∥ 2 + max f 1 (x k+1 ) -f ℓ 2 (x k+1 ) -τ ℓ f + ν k +2κ 2 ∥x k+1 -x∥ 2 c 1 (x k+1 ) -cℓ 2 (x k+1 ) -τ ℓ c + ν k +2κ 2 ∥x k+1 -x∥ 2 (by (5.3)) ≥ -κ 2 ∥x k+1 -x∥ 2 + max f 1 (x k+1 ) -f 2 (x k+1 ) -τ ℓ f + 2κ 2 ∥x k+1 -x∥ 2 c 1 (x k+1 ) -c 2 (x k+1 ) -τ ℓ c + 2κ 2 ∥x k+1 -x∥ 2 (by (4.11)) = H(x k+1 ; x) + κ 2 ∥x k+1 -x∥ 2 .
(by (3.10))

As x = x ℓ and Ȟk (x ℓ ; x ℓ ) = H(x ℓ ; x ℓ ) due to (5.2), we have shown that the descent test (4.5) is satisfied at

x k+1 ̸ = x ℓ : H(x k+1 ; x ℓ ) ≤ H(x ℓ ; x ℓ ) - κ 2 ∥x k+1 -x ℓ ∥ 2 ,
contradicting thus the assumption that only null steps are performed for k > ℓ. Hence, x = x and the last stability center solves (5.4). This allows us to conclude (thanks to convexity of H(•; x)) that x = x ℓ solves min x∈X H(x; x ℓ ). Lemma 5.2 then concludes the proof.

Remark 5.5 (Bundle compression). It is worth mentioning that the index sets B k f and B k c gathering the information bundle can be kept bounded; each one having at most Mmax indices, for a chosen integer Mmax ≥ 3. Indeed, it suffices to keep in the bundles the linearizations issued by the stability center x ℓ k , the new trial point x k+1 and the so-called aggregate linearization as in [START_REF] Correa | Convergence of some algorithms for convex minimization[END_REF]Eq. 4.5]. When transcribed to our setting, the aggregate linearizations for f 1 and c 1 read as = {k + 1, ℓ k , a k c }. Indeed, Proposition 5.4 is still valid if the algorithm employs the above economical rule for updating the bundles: the key Proposition 4.3 from [START_REF] Correa | Convergence of some algorithms for convex minimization[END_REF] still applies and thus the displayed equations right after (5.4) hold. As it can be noted in the sequel, no bundle management restriction (besides the requirement that

f a k f 1 (x) := f k 1 (x k+1 ) + ⟨p k f , x -x k+1 ⟩ ≤ f 1 (x) ∀x ∈ R n ca k c 1 (x) := čk 1 (x k+1 ) + ⟨p k c , x -x k+1 ⟩ ≤ c 1 (x) ∀x ∈ R n , with p k f := j∈B k f α j f s j f 1 , p k c := j∈B k c α j c s j c
k + 1 ∈ B k+1 f ∩ B k+1 c
) is required after a serious steps. □

We consider now the case of infinitely many serious steps. To this end, we need the following auxiliary result.

Lemma 5.6. There exist constants L, M > 0 such that, for all k ∈ N, the three following conditions hold for

p k+1 ∈ ∂ 1 Ȟk (x k+1 ; x ℓ k ), s k+1 X ∈ N X (x k+1 ), and e k+1 = L||x k+1 -x ℓ k ||: ||p k+1 + s k+1 X || ≤ µmax ||x k+1 -x ℓ k || ≤ M, (5.5a
)

p k+1 + s k+1 X ∈ ∂ e k+1 H(x ℓ k ; x ℓ k ) + i X (x ℓ k ) , (5.5b
)

p k+1 ∈ ∂ e k+1 H(x ℓ k ; x ℓ k ).
(5.5c)

Proof. As µ k ∈ (0, µmax] (c.f. Lemma 5.1), expression (4.9) yields the first inequality in (5.5a). Recall that the iterates x k+1 and x ℓ k are contained in the bounded set X for all k. The second inequality in (5.5a) then follows. Convexity of the function Ȟk + i X and (4.9) gives that, for all x ∈ R n , Ȟk (x;

x ℓ k ) + i X (x) ≥ Ȟk (x k+1 ; x ℓ k ) + ⟨p k+1 + s k+1 X , x -x k+1 ⟩ ≥ Ȟk (x k+1 ; x ℓ k ) + ⟨p k+1 + s k+1 X , x -x ℓ k ⟩ + ⟨p k+1 + s k+1 X , x ℓ k -x k+1 ⟩ ≥ Ȟk (x k+1 ; x ℓ k ) + ⟨p k+1 + s k+1 X , x -x ℓ k ⟩ -M ||x ℓ k -x k+1 ||, (5.6)
where the last inequality is due to (5.5a) and Cauchy-Schwarz inequality. Definition (4.6) of Ȟk (•; x ℓ k ) as well as the fact that ℓ k ∈ B k f ∩ B k f give the following chain of inequalities:

Ȟk (x k+1 ; x ℓ k ) ≥ max f ℓ k 1 (x k+1 ) - f ℓ k 2 (x k+1 ) -τ f (x ℓ k ), cℓ k 1 (x k+1 ) - cℓ k 2 (x k+1 ) -τc(x ℓ k ) ≥ max f 1 (x ℓ k ) -f 2 (x ℓ k ) -τ f (x ℓ k ), c 1 (x ℓ k ) -c 2 (x ℓ k ) -τc(x ℓ k ) + min ⟨s ℓ k f 1 -s ℓ k f 2 , x k+1 -x ℓ k ⟩, ⟨s ℓ k c 1 -s ℓ k c 2 , x k+1 -x ℓ k ⟩ .
Since X ⊂ O is compact, we have that ∂f 1 (X), ∂c 1 (X), ∂ C f 2 (X), and ∂ C c 2 (X) are bounded sets (see Section 2). Hence, there exist K f > 0 and Kc > 0 such that, ∥s

ℓ k f 1 -s ℓ k f 2 ∥ ≤ K f and ∥s ℓ k c 1 -s ℓ k c 2 ∥ ≤ Kc for all k.
Applying the Cauchy-Schwarz inequality to the inequalities above and recalling that max

f 1 (x ℓ k ) -f 2 (x ℓ k ) -τ f (x ℓ k ), c 1 (x ℓ k ) -c 2 (x ℓ k ) -τc(x ℓ k ) = H(x ℓ k ; x ℓ k ) = H(x ℓ k ; x ℓ k ) by definition, we get Ȟk (x k+1 ; x ℓ k ) ≥ H(x ℓ k ; x ℓ k ) -L 0 ||x k+1 -x ℓ k ||, with L 0 = max{K f , Kc}. (5.7)
Recall that H(x; x ℓ k ) ≥ Ȟk (x; x ℓ k ) for all x ∈ R n and combine (5.6) with (5.7) to obtain

H(x; x ℓ k ) + i X (x) ≥ Ȟk (x; x ℓ k ) + i X (x) ≥ Ȟk (x k+1 ; x ℓ k ) + ⟨p k+1 + s k+1 X , x -x ℓ k ⟩ -M ||x ℓ k -x k+1 || ≥ H(x ℓ k ; x ℓ k ) -(L 0 + M )||x k+1 -x ℓ k || + ⟨p k+1 + s k+1 X , x -x ℓ k ⟩.
We have thus shown (5.5b) with L = M + L 0 . To prove the last inclusion (5.5c), observe that this chain of inequalities remains true if the term i X (x) is excluded together with corresponding subdifferential s k+1 X : for all

x ∈ R n , H(x; x ℓ k ) ≥ Ȟk (x; x ℓ k ) ≥ Ȟk (x k+1 ; x ℓ k ) + ⟨p k+1 , x -x ℓ k ⟩ -M ||x ℓ k -x k+1 || ≥ H(x ℓ k ; x ℓ k ) + ⟨p k+1 , x -x ℓ k ⟩ -(L 0 + M )||x k+1 -x ℓ k ||.
Observe that Ȟk (•; x ℓ k ) given in (4.6) is the pointwise maximum of finitely many affine functions. Hence, its subdifferential is the convex hull of the "active" linearization slopes, i.e., Proposition 2.1 ii) asserts that

∂ 1 Ȟk (x k+1 ; x ℓ k ) := conv s j f 1 -s ℓ k f 2 j∈ Bk f , s j c 1 -s ℓ k c 2 j∈ Bk c , (5.8)
with Bk f and Bk c given in (4.10). Since X ⊂ O is compact, we have that ∂f 1 (X), ∂c 1 (X), ∂ C f 2 (X), and ∂ C c 2 (X) are bounded sets (see Section 2). Thus, (4.9) and (5.8) certificate that the sequence of model's subgradients p k is bounded. This property is used in the proof of the following proposition. Proposition 5.7 (Infinitely many serious steps). Assume that the algorithm performs infinitely many serious steps. Then, any cluster point x ∈ X of the sequence {x ℓ k } k satisfies the necessary optimality condition (3.13).

Proof. We first show that lim

k→∞ ∥x ℓ k+1 -x ℓ k ∥ = 0 . (5.9)
To this end, we must analyze the two cases of Lemma 4.1. In case i), Algorithm 1 produces a feasible point for (1.1) after finitely many serious steps and all subsequent points are feasible. Let x ℓ k 1 be the first feasible serious iterate. Then, Lemma 4.1 i) yields

f (x ℓ k+1 ) ≤ f (x ℓ k ) -κ 2 ∥x ℓ k+1 -x ℓ k ∥ 2 and c(x ℓ k+1 ) ≤ 0 for all k ≥ ℓ k 1 .
The telescopic sum of the first inequality above yields

∞ k=k 1 ∥x ℓ k+1 -x ℓ k ∥ 2 ≤ 2 κ ∞ k=k 1 f x ℓ k -f x ℓ k+1 ≤ 2 κ f x ℓ k 1 -lim k→∞ f x ℓ k+1 .
Since f is finite-valued and continuous over the bounded set X, the right-hand side of the above inequality is finite. Hence, (5.9) holds. Assume now that the sequence {x ℓ k } is infeasible for (1.1). Lemma 4.1 ii) yields

0 < c(x ℓ k+1 ) ≤ c(x ℓ k ) -κ 2 ∥x ℓ k+1 -x ℓ k ∥ 2 for all ℓ.
Once again, by using the telescopic sum we get (5.9).

As X ⊂ O is compact, with O an open set contained in the domains of component functions, and the generalized subdifferential is locally compact, we conclude that for s

ℓ k f 2 ∈ ∂ C f 2 (x ℓ k ) and s ℓ k c 2 ∈ ∂ C c 2 (x ℓ k ) x ℓ k , s ℓ k f 2 and s ℓ k c 2
are bounded sequences.

By taking subsequences, we can define an index set K ⊂ {0, 1, 2, . . .} such that lim K∋k→∞

x ℓ k = x ∈ X, lim K∋k→∞ s ℓ k f 2 = sf 2 ∈ ∂ C f 2 (x) and lim K∋k→∞ s ℓ k c 2 = sc 2 ∈ ∂ C c 2 (x),
where the two last limits are due to the fact that the generalized subdifferential is outer-semicontinuous [4, Prop. 2.1.5(b)]. Let us now define φ k (•) = H(•; x ℓ k ) and φ(•) = H(•; x), with the latter defined in (5.1) with y = x and the pair of generalized subgradients (s f 2 , sc 2 ) above. For every x ∈ R n fixed, the above limits imply lim K∋k→∞ H(x;

x ℓ k ) = H(x; x), i.e., {φ k } k∈K converges pointwise to φ. With k + 1 = ℓ k+1 , recall that p ℓ k+1 ∈ ∂ e ℓ k+1 φ ℓ (x ℓ k ) due to
Lemma 5.6, Eq. (5.5c). Furthermore, as the sequence {p ℓ k+1 } is bounded (see the paragraph right before this proposition), we may take another subsequence indexed by K ′ ⊂ K so that lim K ′ ∋ℓ→∞ p ℓ k+1 = p ∈ R n and lim K ′ ∋ℓ→∞ e ℓ k+1 = 0 in view of (5.9) and definition of e k+1 given in Lemma 5.6. With these conditions at hand, Lemma A.1 (in Appendix A) ensures that p ∈ ∂φ(x), i.e., p ∈ ∂ 1 H(x; x). Next, observe that {s ℓ k+1 X } is a bounded sequence as the inequality

||p ℓ k+1 + s ℓ k+1 X || ≤ µmax||x ℓ k+1 -x ℓ k || (5.10)
holds due to Lemma 5.6 (with µmax finite due to Lemma 5.1). By definition of the convex normal cone, it follows that there exists a suitable subsequence of {s

ℓ k+1 X } k∈K ′′ , with K ′′ ⊂ K ′ converging to a cluster point s ∈ N X (x) = ∂i X (x). Hence, since X is polyhedral and ri(Dom( H(•; x))) = O ̸ = ∅, p + s ∈ ∂ 1 H(x; x) + ∂i X (x) = ∂ 1 H(x; x) + i X (x) . (5.11)
Finally, inequality (5.10) combined with (5.9) yield p + s = 0. Hence, 0 ∈ ∂ 1 H(x; x) + i X (x) , showing that x minimizes H(•; x) over X. Lemma 5.2 thus concludes the proof.

The following theorem sums up the algorithm's convergence analysis. Furthermore, concerning the original problem (1.1), the following holds:

i) If c(x) > 0, then x is a critical point of min x∈X c 1 (x) -c 2 (x).
ii) If c(x) = 0 and x is not a critical point of min x∈X c 1 (x) -c 2 (x), then x satisfies the criticality condition (3.9) with λ > 0.

iii) If c(x) < 0, then x satisfies the criticality condition (3.9) with λ = 0. If Tol > 0, then the algorithm stops after finitely many steps k ∈ N with an approximate critical point x ℓ k of (3.13).

Proof. For the case Tol = 0, condition (3.13) follows directly from Proposition 5.3 if {x ℓ k } is finite, from Proposition 5.4 if the the algorithm produces only finitely many serious steps followed by an infinite sequence of null steps, and from Proposition 5.7 if infinitely many serious steps are produced.

Furthermore, the connection with the necessary optimality condition for the original problem (1.

Illustrative numerical examples

A deep analysis of the numerical performance of Algorithm 1 is out of the scope of this paper. Instead, this section aims to illustrate our approach to solving some challenging test problems. Here we have two goals:

show that it provides good-quality critical points (examples of Subsection 6.1, Subsection 6.2 and Subsection 6.4) and is able to solve problems that, to the best of our knowledge, could not be resolved with other solvers (example of Subsection 6.3). We consider three nonconvex stochastic optimization problems and one coming from signal processing. Notice that the stochastic problems do not have explicit DC decompositions, and thus the algorithm programming algorithms are not directly applicable.

We have coded Algorithm 1 in MATLAB using Gurobi for solving the master QP subproblem (4.8). Our implementation allows for simpler problems without non-linear constraints (as in the case of the problem in Subsection 6.1). We invite the interested reader to check Appendix B for a brief presentation of how the approach can be simplified in this case. As for the algorithm's parameters, we have chosen Tol = 10 -4 , ρ = σ = 1 2 , and κ = 0.3. Parameter µ 0 > 0 has been selected differently (in the range [10 -2 , 10 2 ]) depending on the problem.

Investment like problems

We will see how the structure of Example 1.1 can appear in practice. Here we will follow the general discussion in [START_REF] Van Ackooij | On investment in power systems[END_REF]. We are interested in the situation wherein we dispose of a set of different technologies i = 1, ..., m capable of generating electricity. Each technology comes with a specific and detailed set of constraints P i , cost function c i attributing to p i ∈ R T the cost of generation. Altogether, the various technologies are meant to ensure the satisfaction of a given customer load d ∈ R T . We are interested in finding the optimal mix. Thus for i = 1, ..., m, we are given θ i ∈ N, the number of "units" of a given type we would like to invest in. The vector θ comes with an investment cost F (θ). In a deterministic setting this would amount to solving min θ∈Θ,p

j i ∈P i F (θ) + m i=1 θ i j=1 c i (p j i ) s.t. m i=1 θ i j=1 p j i ≥ d.
Now should for each i, the mappings c i as well as the feasible sets P i be convex, then it must be so that the averaged solution:

p * i = 1 θ i θ i
j=1 (p j i ) * , in which each power plant of technology i produces exactly this amount is also optimal. This follows from using convexity of P i showing feasibility of p * i and through using Jensen's inequality for c i . This is also exactly what would happen if we would solve the subproblems of the Lagrangian dual (w.r.t. the load constraint) for a given fixed investment vector. The convexifying effect of the Lagrangian is well known, e.g., [START_REF] Lemaréchal | A geometric study of duality gaps, with applications[END_REF][START_REF] Van Ackooij | Decomposition algorithm for large-scale two-stage unit-commitment[END_REF] and thus for this dual setting convexity of c i or P i would not be essential. Either way, as a result we may thus assume that each power plant of the same technology produces the same amount. This would thus lead to the simpler problem (less variables):

min θ∈Θ,p i ∈P i F (θ) + m i=1 θ i c i (p i ) s.t. m i=1 θ i p i ≥ d.
We will investigate a two-stage stochastic version of the last problem, wherein d is uncertain. We thus define:

Q(θ, ξ) = min p i ∈P i m i=1 θ i c i (p i ) s.t. m i=1 θ i p i ≥ d. (6.1)
and consider the optimization problem

min θ∈Θ F (θ) + E[Q(θ, ξ)], (6.2)
where for the sake of simplicity we will assume θ to be allowed to take continuous values (Θ is a polytope).

We will also assume that the feasible set P i is convex, although one could investigate problem (6.2) without this assumption -for instance by arguing through Lagrangian duality.

Let us now look at a concrete instance. We will consider a time horizon of t = 1, ..., T time steps where each time step is considered to be ∆t hours long. The problem disposes of m types of technology, having the following characteristics. Each technology type has a maximum power output level p mx i , proportional cost c i and gradient condition g i . Additionally, each unit is assumed to dispose of a carbon emission rate e i , and the system subject to a carbon cost f . Concretely this means that the proportional cost gets updated through the formula c i ← c i + f e i .

The system is moreover endowed with a given customer load that we will assume to be multivariate Gaussian with a given mean and positive definite Covariance matrix. We refer to [36, § 5] for the description of P i (polyhedral). Furthermore for the various technologies we will assume that F (θ) = F ⊤ θ. The purpose of our experiment is to showcase how concretely the new algorithm can process specifically structured problems such as these.

Following the description of Example 1.1, we need to compute Q ε (θ, ξ) at each given θ. The latter involves the solution of a convex optimization problem, wherein ψt is given by the t-th component of d -m i=1 θ i p i . As a result of the logarithm, the objective function defining Q ε is convex non-linear in y. We will therefore use a cutting plane approach to internally compute Q ε , as well as it's gradient. The inner optimization is initialized from the optimal solution y 0 of the inner optimization problem of Q(x, ξ). The latter can be computed by solving a linear program. This will give us the oracle for f 2 (in the notation of Subsection 4.1).

Table 1 provides the concrete data. is more a demonstration of the capabilities of the algorithm rather than an attempt to provide practical insights into investment regarding the electrical system. The nominal set of assets consists of one asset of each type. This is the initial vector chosen to start the algorithm. The first run considers a situation wherein f = 0, i.e. with zero carbon cost. The algorithm in a total of 13 iterations performs a single serious step and essentially keeps the nominal investment vector fixed -the total installed capacity is 4000 MW. In a second run we have set f = 100 to see the potential impact of such a penalized setting for emitting technologies. In this case the algorithm updates the set of installed capacities quite significantly, by shifting essentially all generation from carbon emitting technologies to technologies 1 -3 not emitting CO2 at all. The algorithm performs a total of 37 iterations and 28 serious steps. The resulting total amount of installed capacity is 3500, slightly less that the nominal vector. This can be explained since the original setting was slightly overcapacitated -and as a result of the introduction of the "fictive" imbalance unit -infeasibility is no longer an issue. The maximum load over the considered scenarios being roughly 3500 as well. The overall objective function value is of the same magnitude as the earlier run. This numerical experiment thus clearly shows that the algorithm has potential to provide meaningful solutions for this type of problem.

Decision dependent probability constraints in two stage problems

In this section we consider the following stochastic problem having different random vectors:

     min x∈X f 1 (x) + S s=1 πsQ(x; ξ s ) s.t. P[A 1 x + b 1 ≥ ω 1 ] ≥ p 1 with Q(x; ξ) := min y∈Y q(x, y; ξ) s.t. Px[A 2 (ξ)y + b 2 (ξ) ≥ ω 2 (x)] ≥ p 2 .
In this problem, ξ ∈ Ξ := {ξ 1 , . . . , ξ S }, ω 1 ∼ N (µ 1 , Σ 1 ), and ω 2 (x) ∼ N (µ 2 (x), Σ 2 (x)). The latter random vector depends on the first-stage decisions. We assume that the covariance matrices Σ 1 and Σ 2 (x) are positive definite for all x ∈ X. As a result, the probability functions are twice-differentiable [START_REF] Henrion | A gradient formula for linear chance constraints under Gaussian distribution[END_REF][START_REF] Van Ackooij | Second-order differentiability of probability functions[END_REF]. Furthermore, as the multivariate Gaussian distribution is log-concave, we get that c 1 (x) = log(p 1 ) -log(P

[A 1 x + b 1 ≥ ω 1 ]
) is a convex function and so is the objective of the penalized subproblem

Q ε (x; ξ) = min y∈Y q(x, y; ξ) -1 ε log Px[A 2 (ξ)y + b 2 (ξ) ≥ ω 2 (x)] -p 2 .
We are thus in the setting of Example 1.1 with f 2 (x) = S s=1 πs[-Q ε (x; ξ s )]. We can observe that the just given optimization problem is a version of two-stage stochastic program having unhedgeable, or post-decision random realizations. Now in order to compute the gradient of both of the involved probability functions, we can rely on two different formulae for the gradients. The mapping c 1 is continuously differentiable and its gradient can be evaluated by employing the results shown in [START_REF] Van Ackooij | On probabilistic constraints induced by rectangular sets and multivariate normal distributions[END_REF]. The second stage probability function is also differentiable under a mild regularity condition, its gradient can be evaluated using the formulae from [START_REF] Van Ackooij | Sub-) Gradient formulae for probability functions of random inequality systems under Gaussian distribution[END_REF]Thm. 5.1]. Indeed with L 2 (x) the matrix resulting from the Cholesky decomposition Σ 2 (x) = L 2 (x)L 2 (x) ⊤ , we may write

c i 2 (x, y) = P[-A 2 (ξ i )y -b 2 (ξ i ) + µ 2 (x) + L 2 (x)ω 2 ≤ 0],
where ω 2 ∈ R m , ω 2 ∼ N (0, I). Hence we can observe that:

∇c i 2 (x, y) = v∈S m-1 :J * (v)̸ =∅,|J * * (v)|=1 - χ(ρ(v)) (L 2 (x)) j(v) v ρ(v)∇(L 2 (x)) j(v),. v) + ∇µ 2 (x), -(A 2 ) ⊤ j(v),: dµ ζ (v) (6.3) with J * (v) = {j = 1, ..., r : (L 2 (x)) j v > 0} ρ(v) = min j∈J * (v) A 2 y + b 2 -µ 2 (x) (L 2 (x)) j v J * * (v) = j ∈ J * (v) : ρ(v) = A 2 y + b 2 -µ 2 (x) (L 2 (x)) j v .
and j(v) being the unique element of J * * (v). In this case since L 2 has linearly independent rows -which is the case since Σ 2 is positive definite -the aforementioned regularity condition (R2CQ) holds true. In fact (R2CQ) holds true locally and as a consequence it is indeed so that both c 1 and c 2 are twice continuously differentiable. This was already clear for c 1 upon using well known classic arguments.

Let us now consider the following concrete example of a problem of this kind. We are interested in a situation considering a manufacturer capable of producing two different products. The first-stage decision variables of the problem consist of setting prices for the products and an advertisement levels. The price will be assumed to be in relation to the average second-stage demand for the given product. We will use the following rule µ 2 (x) = (μ 1 /x 1 , μ2 /x 2 ), with x 1 , x 2 being the price levels for product 1 and 2 respectively. Advertisement is assumed to have a beneficial effect on the variance of the demand, but simultaneous advertisement for both products will be counterproductive. In other words:

Σ 2 (x) := (0.1μ 1 /x 1 x 3 ) 2 -0.4(0.01μ 1 /x 1 x 3 μ2 /x 2 x 4 ) -0.4(0.01μ 1 /x 1 x 3 μ2 /x 2 x 4 ) (0.1μ 2 /x 2 x 4 ) 2 .
The second stage decision y involves the production of the goods. The production process of the goods is subject to some possible unreliability as the amount of actually produced goods are concerned. The matrix A 2 is thus a diagonal matrix, where the first entry is a uniform random variable over the interval [0.9, 1]on average only 95% of the commissioned products actually get manufactured. The second diagonal entry is uniform over the interval [0.8, 1] -the process of production here is more unreliable. However producing with the more unreliable process is slightly cheaper. Any products that are manufactured but not sold, will incur a penalty. The second stage cost function is thus given by

q(x, y, ξ) = (2 -x 1 )y 1 + (1 -x 2 )y 2 + 12E[max {y 1 -(ω 2 ) 1 , 0}] + 12E[max {y 2 -(ω 2 ) 2 , 0}].
The last two terms correspond to the penalization of produced, but not sold goods. It turns out that the latter expectations can be computed "analytically" as they are related to the computation of an expectation of a truncated Gaussian random variable. Therefore, we can observe that the following identity holds true:

E[max {y 1 -(ω 2 ) 1 , 0}] = Φ((y 1 -(µ 2 (x)) 1 )x 1 /(0.1μ 1 x 3 ))(y 1 -(µ 2 (x)) 1 ) - 1 √ 2π e -1 2 ((y 1 -(µ 2 (x)) 1 )x 1 /(0.1 μ1 x 3 )) 2 0.1μ 1 /x 1 x 3 .
The second formula is of course immediately deduced as it is analogous. Both products require a different setup of the factory, so y 1 + y 2 ≤ 10. Furthermore, the first-stage cost function is related to the cost of advertisement f 1 (x) = q 1 x 2 3 + q 2 x 2 4 . Furthermore all first-stage variables are bounded. The implementation of this example requires first the implementation of the formulae for the gradient of the probability function. Here we can exploit the earlier given formula immediately. It can be observed (see the more extensive discussion in [START_REF] Van Ackooij | A discussion of probability functions and constraints from a variational perspective[END_REF]) that the probability value itself can be computed with exactly the same cost. Subsequently the algorithm scheme is very similar to the one of the investment problem. In particular, combining the computations for probability function value and subgradient (6.3) with the reasoning of the previous example, we will obtain the oracle for f 2 component (in the notation of Subsection 4.1), while the oracle for f 1 is straightforward from the formula of the advertisement cost.

Therefore, we have also run this case through the algorithm and found an approximate critical solution after a total of 12 iterations (1600 seconds on personal laptop). The found solution is x = (3.37, 3.21, 0.096, 0.784), showing that there is interest in balancing the prices, i.e., not taking maximal prices, while also investing in advertisement. We have done the same test with IPOPT solver. The computation was aborted after 50 000 seconds with the resulting infeasible point x = (9.997, 9.983, 0.0004, 0.167) slowly approaching the bound (10, 10, 0, 0). This example thus shows that the new algorithm allows us to consider settings beyond classic convexity, even when dealing with probability functions -in this case with decision dependent random vectors.

Highly nonconvex chance-constrained problem

In this section we investigate the following optimization problem (having weakly-convex constraint):

min x∈R n c ⊤ x s.t. P 1 2 ξ ⊤ Q j (x)ξ + q j (x) ⊤ ξ + d j (x) ≤ 0, j = 1, ..., k ≥ p (6.5) x ≤ x ≤ x.
We first note that as a result of [START_REF] Van Ackooij | Gradient formulae for probability functions depending on a heterogenous family of constraints[END_REF] and the upfollowing concrete data, that the probability function is continuously differentiable. Furthermore the underlying feasible set is compact and so we are in case ii) of the introduction: c 1 (x) = p and -c 2 indicating the probability function. The underlying data is not convex in the parameter replaced by the random vector. As a result, the underlying feasible set is not expected to be convex. Concretely we will consider the following data, for k = 2: as well as

Q 1 (x) = 3(x 1 -1) -x 2 -x 2 3(x 1 -1) , Q 2 (x) = -2x 2 x 1 -1 x 1 -1 -2x 2
q 1 (x) = 3 1 (x 1 -1), q 2 (x) = 1 4 x 2 d 1 (x) = -2, d 2 (x) = -2.
We have also picked p = 0.7 together with c = (-1, -1), x = (-2, -2), x = (2, 2). The random vector is taken to be multivariate Gaussian with mean vector 0 and covariance matrix

Σ = 2 -1 -1 2 .
This optimization problem is quite challenging. First it can be observed that an alternative sample average approximation along the lines of [START_REF] Luedtke | A sample approximation approach for optimization with probabilistic constraints[END_REF] would be a MILP. It is thus tempting to first try to solve the resulting optimization problem with such a formulation. We have done so for the following sample sizes {100, 1000, 5000, 10000, 50000}, using CPLEX 12.10. The resulting computation times are 0.5, 4, 22, 100, > 8000 seconds respectively. The last computation was aborted still showcasing a 25.9% gap after more than 2 hours of computation. Unfortunately, none of the obtained solutions turn out to be feasible, quite to the contrary: the typically obtained final probability value is roughly 0.04 being far from the required 0.7. We have also performed a run of a sampled problem with 10000 random realizations, but with a significantly higher probability level of 0.9. In this case, after roughly one hour of computation, the resulting solution being at a MIPgap of 3.3 %, is still highly infeasible having only a probability value of 0.02. Furthermore, we have tested IPOPT solver for the problem resolution: tests have been performed for the six initial points listed in Table 6.3. After at most 26 seconds of computation, IPOPT halted with highly infeasible points with a probability constraint value equal to 0. The difficulty of generating feasible points might come from the form of probability distribution as the probability level sharply raises from zero (blue, Figure 2) to nearly 1 (yellow, Figure 2) when approaching the feasible area from most directions, which causes the loss of gradient information in a large zone of zero probability. But of course this gradient information is not exploited by the MILP formulation at all.

In contrast, as Table 6.3 shows, our algorithm manages to improve the probability level if the starting point is infeasible, and to improve the objective function value while satisfying probability constraint for a feasible initial case. Moreover, for one of the tested starting points we have managed to generate the near (globally) optimal solution (1.2400; -0.1126). Since the problem is indeed very difficult, a precise internal sampling scheme for the probability function is required. This amounts to the number of samples used to compute a formula of the type (6.3), which subsequently leads to design of the oracle for c 2 component. We can prematurely end further sampling, very much like the implementation of Genz' code [START_REF] Genz | Numerical computation of multivariate normal probabilities[END_REF], by checking if sampling variance -in fact the confidence interval bounds -are sufficiently small. Unlike Genz' code a crude antithetic Monte-Carlo scheme has been used for sampling, thus leaving much room for significant improvement in terms of speed and precision, by using for instance QMC as, for instance, in [START_REF] Heitsch | On probabilistic capacity maximization in a stationary gas network[END_REF].

Compressed sensing problem

In this section we focus on the problem of compressed sensing considered in [START_REF] Yu | Convergence rate analysis of a sequential convex programming method with line search for a class of constrained difference-of-convex optimization problems[END_REF]:

min x∈R n ∥x∥ 1 -∥x∥ (6.6) s.t. ∥Ax -b∥ LL 2 ,γ ≤ δ,
where A ∈ R q×n is a full row rank matrix and b ∈ R q . For given γ > 0, Lorentzian norm ∥ • ∥ LL 2 ,γ of a vector y ∈ R q is defined as

∥y∥ 2 LL 2 ,γ = q i=1 log 1 + y 2 i γ 2 .
As discussed in [START_REF] Yu | Convergence rate analysis of a sequential convex programming method with line search for a class of constrained difference-of-convex optimization problems[END_REF], the problem (6.6) is DC with twice continuously differentiable constraint, whose modulus of Lipschitz gradient is known. This allows us to construct the oracle for the constraint component. To compute subgradient of the component f 1 (x) = ∥x∥ 1 , we have chosen the sign function.

As in [START_REF] Yu | Convergence rate analysis of a sequential convex programming method with line search for a class of constrained difference-of-convex optimization problems[END_REF], we have taken q = 3 600, n = 12 800 and generated A ∈ R q×n with normally distributed random entries normalizing it so that each column has a unit norm. To set the original point, we have chosen a subset of size s 0 = [ q 9 ] among basis vectors and generated a s 0 -sparse vector x orig with i.i.d. normally distributed random entries. We have taken b = Ax orig + 0.01η, each η i having a standard Cauchy distribution, and δ = 1.1∥0.01η∥ LL 2 ,γ with γ = 0.02.

The CwC-bundle algorithm manages to obtain a feasible solution after 60 iterations and 52 serious steps, as well as to recover a global solution within tolerance 9.5 × 10 -4 after 2400 iterations and 1613 serious steps with execution time of 8800 seconds (on personal laptop). This example shows the ability of our algorithm to solve problems with an explicit DC structure, even if its performance is not as good as for the method in [START_REF] Yu | Convergence rate analysis of a sequential convex programming method with line search for a class of constrained difference-of-convex optimization problems[END_REF] developed for more specific framework.

Conclusion

In this manuscript, we have considered nonsmooth and nonconvex optimization problems where the objective function and nonlinear constraint are represented as the difference of convex and weakly convex functions (CwC). Our work studies various stationary conditions and a bundle method approach enabling to compute critical (generalized KKT) points. The latter broadens and enhances the algorithm developed in [START_REF] Van Ackooij | A bundle method for nonsmooth DC programming with application to chance-constrained problems[END_REF] for the case of Difference-of-Convex (DC)-constrained DC-problems, and likewise relies on problem reformulation via an improvement function. To the best of our knowledge, proposed method is the first one that directly exploits the CwC-structure of the involved functions and does not require additional assumptions or transformations as, for instance, explicit Difference-of-Convex decompositions or Moreau envelopes. We have illustrated the method performance with a few stochastic problems, including two-stage and chanceconstrained problems and a compressed sensing problem with nonlinear constraint. Preliminary results are meaningful and show that the algorithm enables tackling settings beyond the classic Difference-of-Convex setting.

Appendices A Some mathematical results

Proof of Proposition 2.4. Since f : O → R n is weakly convex, it follows by definition that, relative to some neighbourhood V x ′ of each point x ′ ∈ O, there exist µ x ′ > 0 such that for all µ ≥ µ x ′ the function ϕ(x) = f(x) + µ 2 ∥x∥ 2 is finite and convex on V x ′ . In such a representation, there is no loss of generality in assuming that V x ′ ⊂ O (if necessary we can define a new/smaller neighbourhood as V x ′ ∩ O for which the above conclusion obviously stands). By considering all the points in X, let V := {V To show item ii), let us now define φ(x) = f(x) + µ 2 ∥x∥ 2 + i O ′ (x), an extended real-valued convex function: φ : R n → R ∪ {∞}. Note that for each x ∈ O ′ , there exists a neighbourhood Vx ⊂ O ′ such that φ(x ′ ) = ϕ(x ′ ) for all x ′ ∈ Vx. This fact permits us to conclude that ∂ φ(x) = ∂ϕ(x) for all x ∈ O ′ . It thus follows from item i) that, for every x ∈ O ′ and every s ∈ ∂ φ(x), there exists s f ∈ ∂ C f (x) such that Proof. First, let us prove that lim inf ℓ φ ℓ (x ℓ ) ≥ φ(x). Since dom(φ ℓ ) = dom(φ) = R n , it follows from [32, Cor. 2C] that the pointwise convergence of {φ ℓ } ℓ∈N is equivalent to epi-convergence, which in turn is equivalent (see [START_REF] Salinetti | On the relations between two types of convergence for convex functions[END_REF]Eq. (4.2)]) to epi-convergence of {φ * ℓ } ℓ∈N , the sequence of conjugate functions to φ ℓ . Hence, it follows that lim ℓ φ * ℓ (x) = φ * (x) for every given x ∈ R n . Now consider the following development: Recall that g ℓ ∈ ∂ ϵ ℓ φ ℓ (x ℓ ). Then, φ ℓ (x) ≥ φ ℓ (x ℓ ) + ⟨g ℓ , x -x ℓ ⟩ -ϵ ℓ for all x ∈ R n . By taking the limit when ℓ goes to infinity we get

φ ℓ (x ℓ ) = (φ * ℓ ) * (x ℓ ) = sup

B Simplified algorithm for the case without nonlinear constraints

This section describes how Algorithm 1 can be simplified to deal with the simpler convexly-constrained problem min x∈X f (x), with f (x) = f 1 (x) -f 2 (x). (B.1)

In this case, the problem's model (4.6) reduces to Ȟk (x; x ℓ k ) = f k 1 (x) -f ℓ k 2 (x), and the descent test (4.5) becomes f (x k+1 ) ≤ f (x ℓ k ) -κ 2 ∥x -x ℓ k ∥ 2 . Hence, Algorithm 1 boils down to the following plainer scheme. Convergence analysis for Algorithm 2 follows from that of Algorithm 1 upon several simplifications. Instead of doing this exercise, we simply state the following result. If Tol > 0, then the algorithm stops after finitely many steps k ∈ N with an approximate critical point x ℓ k of (3.3).

To have an intuition of why the above theorem is valid, the reader may think of adding a dummy convex nonlinear convex function c(x) = c 1 (x) -0 to (B.1) and rely on the results from Sections 4 and 5. Indeed, by selecting a constant M > 0 large enough and function c such that c(x) ≤ -M < 0 for all x ∈ X, we can see that Algorithm 1 applied to (B.1) with the additional and superfluous constraint c(x) ≤ 0 boils down to Algorithm 2. Furthermore, in this artificial setting, the above convergence result follows directly from Theorem 5.8, item iii).

Example 1 . 1 (c 1

 111 Nonconvex two-stage programming). Let Ξ := {ξ 1 , . . . ξ S } be a set of scenarios and πs > 0 the probability of occurrence of event ξ s , s = 1, . . . , S. Consider the following nonconvex two-stage program (x) -c 2 (x) ≤ 0 with Q(x; ξ) := min y∈Y q(x, y; ξ) s.t. ψ i (x, y; ξ) ≤ 0, i = 1, . . . , m.

  otherwise. The convex hull of a set X is convX and the relative interior is denoted by ri X. The domain of a function φ : R n → (-∞, +∞] is represented by Dom(φ) = {x ∈ R n : φ(x) < +∞}. Notation O stands for an open convex set of the Euclidean space R n and, given the definitions of f and c, we have that O ⊂ Dom(f ) and O ⊂ Dom(c). The component functions of f and c are f 1 , f 2 , and c 1 , c 2 respectively: f 1 and c 1 are convex, whereas f 2 and c 2 are weakly convex on some neighbourhood of every x ∈ O. Finally, f * stands for the Legendre-Fenchel transform of a function f : R n → (-∞, +∞].

  ) with τ f (y) = f (y) + ρ[c(y)] + and τc(y) = σ[c(y)] + , for ρ ≥ 0 and σ ∈ [0, 1).

Figure 1 :

 1 Figure 1: Function f (x) = 1 2 x 2 -x in red and c(x) = max{x, 2x} -max{2x, 4x} in blue.

  .15) Since F and G are convex functions, the criticality condition (3.3) for (3.11) (under this DC decomposition) reads as 0 ∈ ∂ 1 F (x; x) -∂G(x) + N X (x), (3.16) where ∂ 1 F stands for the subdifferential of F with respect to the first argument. It turns out that our new condition (3.13) is stronger than (3.16), used in [49]. Lemma 3.7. In addition to our basic assumptions on problem (1.1), suppose that f 2 and c 2 are convex. Then the necessary optimality condition (3.13) implies (3.16).

  (a) If (4.5) holds, then declare a serious step: define ℓ k+1 := k + 1, choose B k+1 f , B k+1 c ⊂ {0, . . . , k + 1} with {k + 1} ∈ B k+1 f ∩ B k+1 c and arbitrarily select µ k+1 ∈ (0, µ k ]. (b) Else, declare a null step: define ℓ k+1 := ℓ k and choose B k+1 f , B k+1 c ⊂ {0, . . . , k + 1} with Bk f ∪ {k + 1, ℓ k } ⊂ B k+1 f and Bk c ∪ {k + 1, ℓ k } ⊂ B k+1 c ( Bk f and Bk c as in (4.10)). Compute ν k by (4.11). If

1

 1 and multipliers α f , αc as in (4.10). We claim that the following economical rule for managing B k f and B k c (in Step 4 of Algorithm 1) is enough to ensure convergence: Serious step: set B k+1 f = {k + 1} and B k+1 c = {k + 1}; Null step: set B k+1 f = {k + 1, ℓ k , a k f } and B k+1 c

Theorem 5 . 8 (

 58 Convergence analysis). Let X ̸ = ∅ be a bounded polyhedron contained in the open set O, f 1 , c 1 : O → R convex, and f 2 , c 2 : O → R weakly convex functions on some neighbourhood of each x ∈ O. If in Algorithm 1 the stopping test tolerance Tol is set to zero, then any cluster point x of the sequence of stability centers {x ℓ k } satisfies the necessary optimality condition (3.13).

  1) is established by Theorem 3.5. Proposition 5.4 ensures that lim k→∞ ||x k+1 -x ℓ || = 0 if x ℓ is the last stability center. Otherwise, lim k→∞ ||x ℓ k+1 -x ℓ k || = 0, as shown in the proof of Proposition 5.7. Thus, Algorithm 1 stops after finitely many steps provided Tol > 0.

Figure 2 :

 2 Figure 2: Probability distribution associated with the chance constraint in (6.5)

Figure 3 :

 3 Figure 3: Computed solution (marked by circle) of (6.6) and x orig (marked by asterisk)

  x ′ : x ′ ∈ X} be the collection of all such neighbourhoods. Then, by construction, V is an open cover of the compact set X and, by definition of compactness, it has a finite open subcover, i.e., there exists finitely many points {x′ 1 , . . . , x ′ m } ⊂ X such that O ′ := ∪ m i=1 V x ′ i ⊃ X,and by construction O ′ is an open subset of O. The first part of item i) thus follows by taking µ f := max i=1,...,m µ x ′ i < ∞. By writing f(x) = ϕ(x) -µ 2 ∥x∥ 2 and recalling Proposition 2.1 i) we get ∂ C f(x) = ∂ϕ(x) -µx for all x ∈ O ′ . This concludes item i).

s

  = s f + µx and the subgradient inequality reads as φ(y) ≥ φ(x) + ⟨s f + µx, y -x⟩ ∀y ∈ R n , i.e., f(y)+ µ 2 ∥y∥ 2 + i O ′ (y) ≥ f(x) + µ 2 ∥x∥ 2 + i O ′ (x) + ⟨s f + µx, y -x⟩ for all y ∈ R n . The latter simplifies to f(y) + i O ′ (y) ≥ f(x) + ⟨s f , y -x⟩ -µ 2 ∥y -x∥ 2 ∀y ∈ R n .By restricting y to the set X and recalling thats f = s -µx ∈ ∂ C f(x)is an arbitrary subgradient (because no restriction was imposed to s ∈ ∂ φ(x)), the above inequality becomes (2.4). □ Lemma A.1. Let φ : R n → R be a convex function, and {φ ℓ } ℓ∈N a sequence of convex functions φ ℓ : R n → R converging pointwise to φ, i.e., lim ℓ→∞ φ ℓ (x) = φ ℓ (x) for every given point x. Furthermore, let {x ℓ } ⊂ R n be such that lim ℓ→∞ x ℓ = x and {ϵ ℓ } ⊂ R + satisfy lim ℓ→∞ ϵ ℓ = 0. If g ℓ ∈ ∂ ϵ ℓ φ ℓ (x ℓ ) for all ℓ and lim ℓ→∞ g ℓ = ḡ, then ḡ ∈ ∂φ(x).

  y∈R n [⟨y, x ℓ ⟩ -φ * ℓ (y)] ≥ ⟨y, x ℓ ⟩ -φ * ℓ (y) ∀y ∈ R n . Accordingly, lim inf ℓ φ ℓ (x ℓ ) ≥ lim inf ℓ [⟨y, x ℓ ⟩ -φ * ℓ (y)] = ⟨y, x⟩ -φ * (y) for all y ∈ R n , showing that lim inf ℓ φ ℓ (x ℓ ) ≥ sup y [⟨y, x⟩ -φ * (y)] = (φ * ) * (x) = φ(x).

  φ(x) = lim ℓ φ ℓ (x) = lim inf ℓ φ ℓ (x) ≥ lim inf ℓ [φ ℓ (x ℓ ) + ⟨g ℓ , x -x ℓ ⟩ -ϵ ℓ ] ≥ lim inf ℓ φ ℓ (x ℓ ) + lim inf ℓ ⟨g ℓ , x -x ℓ ⟩ -lim sup ℓ ϵ ℓ ≥ φ(x) + ⟨ḡ, x -x⟩,showing that ḡ ∈ ∂φ(x).

Theorem B. 1 .

 1 Consider problem (B.1) with X ̸ = ∅ a bounded polyhedron contained in the open set O, f 1 : O → R convex, and f 2 : O → R weakly convex on some neighbourhood of each x ∈ O. If in Algorithm 2 the stopping test tolerance Tol is set to zero, then any cluster point x of the sequence of stability centers {x ℓ k } satisfies the necessary optimality condition (3.3).

  Proposition 2.1. Let ft : O → R, t = 1, 2, . . . , m, be locally Lipschitz functions and x ∈ O an arbitrary point. Then i) ∂ C [ m t=1 atft](x) ⊂ m t=1 at∂ C ft(x) for all a ∈ R m , and equality holds if • all but one of ft are smooth [4, Prop. 2.3.3 and Cor. 2];

	• or if every ft is regular at x and a ∈ R m + [4, Cor. 3];
	ii) ∂ C f(x) ⊂ conv{∂ C ft(x) : t ∈ I(x)}, for f(x) = max t=1,...,m ft(x) and I(x) := arg max t=1,...,m ft(x),
	and equality holds and f is regular if every ft is regular at x. [4, Prop. 2.3.12].

2.3.6].

  When f is a smooth function with a Lipschitz continuous gradient, then f is weakly convex in the global sense with µ = L the Lipschitz constant of ∇f [6, Prop. 1]. Moreover, it follows from [7, Lemma 4.2] that the family of composite functions given in item v) of the Introduction is also weakly convex in the global sense. Definition 2.2 implies that weakly convex functions are locally DC : the decomposition f(x) = ϕ(x)µ 2 ∥x∥ 2 holds on some neighbourhood of every x ′ ∈ O. As a result, [29, Thm. 10.33] ensures that the class of weakly convex functions coincide with that of Lower-C 2 functions; see also [33, Thm. 1.3, Cor. 1.3]. Definition 2.3 (Def. 10.29 [29]). (LC 2 functions). A function f : O → R is said to be Lower-C 2 or LC 2 on O if, on some neighbourhood V x ′ ⊂ O of each x ′ ∈ O, there is a representation

Table 1 :

 1 Data for the stylized investment problem

		1	2	3	4	5	6	7	8	9	10
	p mx (MW)	900	900	900	300	300	200	200	200	100	10000
	g (MW/h)	100	100	100	30	30	20	20	70	70	5000
	c (e/MWh) 30	35	37	45	47	60	100	110	150	10000
	F inv (e)	493151 493151 493151 41096 41096 32877 32877 32877 21918 0
	e (t/MW)	0	0	0	1	1	0.5	0.5	0.5	1.1	0

Table 2 :

 2 Results obtained with Algorithm 1 for problem (6.5) depending on the initial point

	Initial point	Time (s)	Iterations	Initial	Final	Objective
				probability	probability	value
				level	level
	(1.2, -0.1) ⊤	47.6	17	0.80	0.70	-1.1097
	(1.5, -0.4) ⊤	78.4	28	0.43	0.70	-1.1346
	(1.6, -0.3) ⊤	36.9	14	0.42	0.70	-1.0815
	(1.7, -0.4) ⊤	59.8	20	0.36	0.70	-1.1443
	(1.1, -0.4) ⊤	21.3	10	0.43	0.43	-0.7001
	(1, -0.8) ⊤	10.8	9	0.18	0.18	-0.2001

It is worth noting that in many practical problems, mainly those from data science, an upper bound on such a constant can be computed.

We can observe that the last unit described in the previous table is an imbalance unit -a computational trick to ensure that one can always meet the load, in this case even despite a potentially completely unbalanced set of invested assets. In terms of constraints on investment, we do not allow investment in this last unit, the capacity will remain at 1. The cost of investment was set up taking inspiration from https://www.e-education.psu.edu/eme801/node/530, upon rescaling to match T and while accounting for life span of the various technologies. The data of the case is stylized and the general purpose of the study

Acknowledgments.

W. Oliveira acknowledges financial support from the Gaspard-Monge Program for Optimization and Operations Research (PGMO) project "Scalable Optimization for Learning and Energy Management."

Algorithm 2 Proximal Bundle Method for Convexly-Constrained CwC programs

Step 0 (Initialization) Let x 0 ∈ X, κ ∈ (0, 1 2 ), κ ≤ µ 0 , and Tol ≥ 0 be given. Call the oracles to compute f i (x 0 ) and (generalized) subgradients s 0 fi , i = 1, 2. Define k := ℓ k = 0 and B 0 f := {0}.

Step 1 (Trial point) Compute x k+1 the (x-part) solution of the QP

Step 2 (Stopping test) If ∥x k+1 -x ℓ k ∥ ≤ Tol, then stop and return x ℓ k .

Step 3 (Oracles call) Compute f i (x k+1 ), and subgradients s k+1 fi , i = 1, 2.

Step 4 (Descent test)

Step 5 (Loop) Set k := k + 1 and go back to Step 1.