1// Copyright 2010-2024 Google LLC
2// Licensed under the Apache License, Version 2.0 (the "License");
3// you may not use this file except in compliance with the License.
4// You may obtain a copy of the License at
6// http://www.apache.org/licenses/LICENSE-2.0
8// Unless required by applicable law or agreed to in writing, software
9// distributed under the License is distributed on an "AS IS" BASIS,
10// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
11// See the License for the specific language governing permissions and
12// limitations under the License.
16package operations_research.pdlp;
18option java_package = "com.google.ortools.pdlp";
19option java_multiple_files = true;
20option csharp_namespace = "Google.OrTools.PDLP";
22import "ortools/glop/parameters.proto";
25 OPTIMALITY_NORM_UNSPECIFIED = 0;
26 OPTIMALITY_NORM_L_INF = 1; // The infinity norm.
27 OPTIMALITY_NORM_L2 = 2; // The Euclidean norm.
29 // The infinity norm of component-wise relative errors offset by the ratio of
30 // the absolute and relative error tolerances, i.e., the l_∞ norm of
31 // [residual / (eps_ratio + |bound|)], where eps_ratio =
32 // eps_optimal_{X}_residual_absolute / eps_optimal_{X}_residual_relative
33 // where {X} is either primal or dual, and bound is the corresponding primal
34 // or dual bound (that is, the violated constraint bound for primal residuals,
35 // and the objective coefficient for dual residuals).
36 // Using eps_ratio in this norm means that if the norm is <=
37 // eps_optimal_{X}_residual_relative, then the residuals satisfy
38 // residual <= eps_optimal_{X}_residual_absolute
39 // + eps_optimal_{X}_residual_relative * |bound|
40 OPTIMALITY_NORM_L_INF_COMPONENTWISE = 3;
43// A description of solver termination criteria. The criteria are defined in
44// terms of the quantities recorded in IterationStats in solve_log.proto.
46// Relevant readings on infeasibility certificates:
47// (1) https://docs.mosek.com/modeling-cookbook/qcqo.html provides references
48// explaining why the primal rays imply dual infeasibility and dual rays imply
49// primal infeasibility.
50// (2) The termination criteria for Mosek's linear programming optimizer
51// https://docs.mosek.com/9.0/pythonfusion/solving-linear.html.
52// (3) The termination criteria for OSQP is in section 3.3 of
53// https://web.stanford.edu/~boyd/papers/pdf/osqp.pdf.
54// (4) The termination criteria for SCS is in section 3.5 of
55// https://arxiv.org/pdf/1312.3039.pdf.
56message TerminationCriteria {
57 // The norm that we are measuring the optimality criteria in.
58 optional OptimalityNorm optimality_norm = 1 [default = OPTIMALITY_NORM_L2];
60 // Define CombineBounds(l,u) as a function that takes two equal-length vectors
61 // and returns a vector whose elements are max{|l_i|, |u_i|} where non-finite
62 // values of l_i and u_i are replaced with zero.
63 // For example, CombineBounds([-2, -Inf, 5], [1, 10, 5]) == [2, 10, 5].
65 // Recalling the definitions from
66 // https://developers.google.com/optimization/lp/pdlp_math, c is the linear
67 // objective vector, l^c are the constraint lower bounds, and u^c are the
68 // constraint upper bounds. Define b^c = CombineBounds(l^c, u^c). Also,
69 // define violated_bound(l^c, u^c, i) to be the bound (l^c or u^c) which is
70 // violated by the current x[i].
72 // When using DetailedOptimalityCriteria the conditions to declare a solution
74 // | primal_objective - dual_objective | <=
75 // eps_optimal_objective_gap_absolute +
76 // eps_optimal_objective_gap_relative *
77 // ( | primal_objective | + | dual_objective | )
78 // If optimality_norm is OPTIMALITY_NORM_L_INF or OPTIMALITY_NORM_L2 (where
79 // norm(x, p) is the l_∞ or l_2 norm):
80 // norm(primal_residual, p) <=
81 // eps_optimal_primal_residual_absolute +
82 // eps_optimal_primal_residual_relative * norm(b^c, p)
83 // norm(dual_residual, p) <=
84 // eps_optimal_dual_residual_absolute +
85 // eps_optimal_dual_residual_relative * norm(c, p)
86 // Otherwise, if optimality_norm is OPTIMALITY_NORM_L_INF_COMPONENTWISE, then,
88 // primal_residual[i] <=
89 // eps_optimal_primal_residual_absolute +
90 // eps_optimal_primal_residual_relative * |violated_bound(l^c, u^c, i)|
91 // dual_residual[i] <=
92 // eps_optimal_dual_residual_absolute +
93 // eps_optimal_dual_residual_relative * |c[i]|
94 // It is possible to prove that a solution satisfying the above conditions
95 // for L_INF and L_2 norms also satisfies SCS's optimality conditions (see
97 // ϵ_pri = ϵ_dual = ϵ_gap = eps_optimal_*_absolute = eps_optimal_*_relative.
98 // (ϵ_pri, ϵ_dual, and ϵ_gap are SCS's parameters).
99 // When using SimpleOptimalityCriteria all the eps_optimal_*_absolute have the
100 // same value eps_optimal_absolute and all the eps_optimal_*_relative have the
101 // same value eps_optimal_relative.
103 message SimpleOptimalityCriteria {
104 // Absolute tolerance on the primal residual, dual residual, and objective
106 optional double eps_optimal_absolute = 1 [default = 1.0e-6];
107 // Relative tolerance on the primal residual, dual residual, and objective
109 optional double eps_optimal_relative = 2 [default = 1.0e-6];
112 message DetailedOptimalityCriteria {
113 // Absolute tolerance on the primal residual.
114 optional double eps_optimal_primal_residual_absolute = 1 [default = 1.0e-6];
115 // Relative tolerance on the primal residual.
116 optional double eps_optimal_primal_residual_relative = 2 [default = 1.0e-6];
117 // Absolute tolerance on the dual residual.
118 optional double eps_optimal_dual_residual_absolute = 3 [default = 1.0e-6];
119 // Relative tolerance on the dual residual.
120 optional double eps_optimal_dual_residual_relative = 4 [default = 1.0e-6];
121 // Absolute tolerance on the objective gap.
122 optional double eps_optimal_objective_gap_absolute = 5 [default = 1.0e-6];
123 // Relative tolerance on the objective gap.
124 optional double eps_optimal_objective_gap_relative = 6 [default = 1.0e-6];
127 // If none are set explicitly the deprecated eps_optimal_absolute and
128 // eps_optimal_relative are used to construct a SimpleOptimalityCriteria.
129 oneof optimality_criteria {
130 SimpleOptimalityCriteria simple_optimality_criteria = 9;
131 DetailedOptimalityCriteria detailed_optimality_criteria = 10;
134 // Absolute tolerance on primal residual, dual residual, and the objective
136 // Deprecated, use simple_optimality_criteria instead.
137 // TODO(b/241462829) delete this deprecated field.
138 optional double eps_optimal_absolute = 2
139 [deprecated = true, default = 1.0e-6];
141 // Relative tolerance on primal residual, dual residual, and the objective
143 // Deprecated, use simple_optimality_criteria instead.
144 // TODO(b/241462829) delete this deprecated field.
145 optional double eps_optimal_relative = 3
146 [deprecated = true, default = 1.0e-6];
148 // If the following two conditions hold we say that we have obtained an
149 // approximate dual ray, which is an approximate certificate of primal
151 // (1) dual_ray_objective > 0,
152 // (2) max_dual_ray_infeasibility / dual_ray_objective <=
153 // eps_primal_infeasible.
154 optional double eps_primal_infeasible = 4 [default = 1.0e-8];
156 // If the following three conditions hold we say we have obtained an
157 // approximate primal ray, which is an approximate certificate of dual
159 // (1) primal_ray_linear_objective < 0,
160 // (2) max_primal_ray_infeasibility / (-primal_ray_linear_objective) <=
161 // eps_dual_infeasible
162 // (3) primal_ray_quadratic_norm / (-primal_ray_linear_objective) <=
163 // eps_dual_infeasible.
164 optional double eps_dual_infeasible = 5 [default = 1.0e-8];
166 // If termination_reason = TERMINATION_REASON_TIME_LIMIT then the solver has
167 // taken at least time_sec_limit time.
168 optional double time_sec_limit = 6 [default = inf];
170 // If termination_reason = TERMINATION_REASON_ITERATION_LIMIT then the solver
171 // has taken at least iterations_limit iterations.
172 optional int32 iteration_limit = 7 [default = 2147483647];
174 // If termination_reason = TERMINATION_REASON_KKT_MATRIX_PASS_LIMIT then
175 // cumulative_kkt_matrix_passes is at least kkt_pass_limit.
176 optional double kkt_matrix_pass_limit = 8 [default = inf];
179message AdaptiveLinesearchParams {
180 // At the end of each iteration, regardless of whether the step was accepted
181 // or not, the adaptive rule updates the step_size as the minimum of two
182 // potential step sizes defined by the following two exponents.
184 // The step size reduction exponent defines a step size given by
185 // (1 - (total_steps_attempted + 1)^(-step_size_reduction_exponent)) *
186 // step_size_limit where step_size_limit is the maximum allowed step size at
187 // the current iteration. This should be between 0.1 and 1.
188 optional double step_size_reduction_exponent = 1 [default = 0.3];
190 // The step size growth exponent defines a step size given by (1 +
191 // (total_steps_attempted + 1)^(-step_size_growth_exponent)) * step_size_.
192 // This should be between 0.1 and 1.
193 optional double step_size_growth_exponent = 2 [default = 0.6];
196message MalitskyPockParams {
197 // At every inner iteration the algorithm can decide to accept the step size
198 // or to update it to step_size = step_size_downscaling_factor * step_size.
199 // This parameter should lie between 0 and 1. The default is the value used in
200 // Malitsky and Pock (2016).
201 optional double step_size_downscaling_factor = 1 [default = 0.7];
203 // Contraction factor used in the linesearch condition of Malitsky and Pock.
204 // A step size is accepted if primal_weight * primal_stepsize *
205 // norm(constraint_matrix' * (next_dual - current_dual)) is less
206 // than linesearch_contraction_factor * norm(next_dual - current_dual).
207 // The default is the value used in Malitsky and Pock (2016).
208 optional double linesearch_contraction_factor = 2 [default = 0.99];
210 // Malitsky and Pock linesearch rule permits an arbitrary choice of the first
211 // step size guess within an interval [m, M]. This parameter determines where
212 // in that interval to pick the step size. In particular, the next stepsize is
213 // given by m + step_size_interpolation*(M - m). The default is the value used
214 // in Malitsky and Pock (2016).
215 optional double step_size_interpolation = 3 [default = 1];
218// Parameters for PrimalDualHybridGradient() in primal_dual_hybrid_gradient.h.
219// While the defaults are generally good, it is usually worthwhile to perform a
220// parameter sweep to find good settings for a particular family of problems.
221// The following parameters should be considered for tuning:
222// - restart_strategy (jointly with major_iteration_frequency)
223// - primal_weight_update_smoothing (jointly with initial_primal_weight)
224// - presolve_options.use_glop
225// - l_inf_ruiz_iterations
226// - l2_norm_rescaling
227// In addition, tune num_threads to speed up the solve.
228message PrimalDualHybridGradientParams {
229 enum RestartStrategy {
230 RESTART_STRATEGY_UNSPECIFIED = 0;
231 // No restarts are performed. The average solution is cleared every major
232 // iteration, but the current solution is not changed.
234 // On every major iteration, the current solution is reset to the average
235 // since the last major iteration.
236 EVERY_MAJOR_ITERATION = 2;
237 // A heuristic that adaptively decides on every major iteration whether to
238 // restart (this is forced approximately on increasing powers-of-two
239 // iterations), and if so to the current or to the average, based on
240 // reduction in a potential function. The rule more or less follows the
241 // description of the adaptive restart scheme in
242 // https://arxiv.org/pdf/2106.04756.pdf.
243 ADAPTIVE_HEURISTIC = 3;
244 // A distance-based restarting scheme that restarts the algorithm whenever
245 // an appropriate potential function is reduced sufficiently. This check
246 // happens at every major iteration.
247 // TODO(user): Cite paper for the restart strategy and definition of the
248 // potential function, when available.
249 ADAPTIVE_DISTANCE_BASED = 4;
252 enum LinesearchRule {
253 LINESEARCH_RULE_UNSPECIFIED = 0;
254 // Applies the heuristic rule presented in Section 3.1 of
255 // https://arxiv.org/pdf/2106.04756.pdf (further generalized to QP). There
256 // is not a proof of convergence for it. It is usually the fastest in
257 // practice but sometimes behaves poorly.
258 ADAPTIVE_LINESEARCH_RULE = 1;
259 // Applies Malitsky & Pock linesearch rule. This guarantees an
260 // ergodic O(1/N) convergence rate https://arxiv.org/pdf/1608.08883.pdf.
261 // This is provably convergent but doesn't usually work as well in practice
262 // as ADAPTIVE_LINESEARCH_RULE.
263 MALITSKY_POCK_LINESEARCH_RULE = 2;
264 // Uses a constant step size corresponding to an estimate of the maximum
265 // singular value of the constraint matrix.
266 CONSTANT_STEP_SIZE_RULE = 3;
269 optional TerminationCriteria termination_criteria = 1;
271 // The number of threads to use. Must be positive.
272 // Try various values of num_threads, up to the number of physical cores.
273 // Performance may not be monotonically increasing with the number of threads
274 // because of memory bandwidth limitations.
275 optional int32 num_threads = 2 [default = 1];
277 // For more efficient parallel computation, the matrices and vectors are
278 // divided (virtually) into num_shards shards. Results are computed
279 // independently for each shard and then combined. As a consequence, the order
280 // of computation, and hence floating point roundoff, depends on the number of
281 // shards so reproducible results require using the same value for num_shards.
282 // However, for efficiency num_shards should a be at least num_threads, and
283 // preferably at least 4*num_threads to allow better load balancing. If
284 // num_shards is positive, the computation will use that many shards.
285 // Otherwise a default that depends on num_threads will be used.
286 optional int32 num_shards = 27 [default = 0];
288 // If true, the iteration_stats field of the SolveLog output will be populated
289 // at every iteration. Note that we only compute solution statistics at
290 // termination checks. Setting this parameter to true may substantially
291 // increase the size of the output.
292 optional bool record_iteration_stats = 3;
294 // The verbosity of logging.
295 // 0: No informational logging. (Errors are logged.)
296 // 1: Summary statistics only. No iteration-level details.
297 // 2: A table of iteration-level statistics is logged.
298 // (See ToShortString() in primal_dual_hybrid_gradient.cc).
299 // 3: A more detailed table of iteration-level statistics is logged.
300 // (See ToString() in primal_dual_hybrid_gradient.cc).
301 // 4: For iteration-level details, prints the statistics of both the average
302 // (prefixed with A) and the current iterate (prefixed with C). Also prints
303 // internal algorithmic state and details.
304 // Logging at levels 2-4 also includes messages from level 1.
305 optional int32 verbosity_level = 26 [default = 0];
307 // Time between iteration-level statistics logging (if `verbosity_level > 1`).
308 // Since iteration-level statistics are only generated when performing
309 // termination checks, logs will be generated from next termination check
310 // after `log_interval_seconds` have elapsed. Should be >= 0.0. 0.0 (the
311 // default) means log statistics at every termination check.
312 optional double log_interval_seconds = 31 [default = 0.0];
314 // The frequency at which extra work is performed to make major algorithmic
315 // decisions, e.g., performing restarts and updating the primal weight. Major
316 // iterations also trigger a termination check. For best performance using the
317 // NO_RESTARTS or EVERY_MAJOR_ITERATION rule, one should perform a log-scale
318 // grid search over this parameter, for example, over powers of two.
319 // ADAPTIVE_HEURISTIC is mostly insensitive to this value.
320 optional int32 major_iteration_frequency = 4 [default = 64];
322 // The frequency (based on a counter reset every major iteration) to check for
323 // termination (involves extra work) and log iteration stats. Termination
324 // checks do not affect algorithmic progress unless termination is triggered.
325 optional int32 termination_check_frequency = 5 [default = 64];
327 // NO_RESTARTS and EVERY_MAJOR_ITERATION occasionally outperform the default.
328 // If using a strategy other than ADAPTIVE_HEURISTIC, you must also tune
329 // major_iteration_frequency.
330 optional RestartStrategy restart_strategy = 6 [default = ADAPTIVE_HEURISTIC];
332 // This parameter controls exponential smoothing of log(primal_weight) when a
333 // primal weight update occurs (i.e., when the ratio of primal and dual step
334 // sizes is adjusted). At 0.0, the primal weight will be frozen at its initial
335 // value and there will be no dynamic updates in the algorithm. At 1.0, there
336 // is no smoothing in the updates. The default of 0.5 generally performs well,
337 // but has been observed on occasion to trigger unstable swings in the primal
338 // weight. We recommend also trying 0.0 (disabling primal weight updates), in
339 // which case you must also tune initial_primal_weight.
340 optional double primal_weight_update_smoothing = 7 [default = 0.5];
342 // The initial value of the primal weight (i.e., the ratio of primal and dual
343 // step sizes). The primal weight remains fixed throughout the solve if
344 // primal_weight_update_smoothing = 0.0. If unset, the default is the ratio of
345 // the norm of the objective vector to the L2 norm of the combined constraint
346 // bounds vector (as defined above). If this ratio is not finite and positive,
347 // then the default is 1.0 instead. For tuning, try powers of 10, for example,
348 // from 10^{-6} to 10^6.
349 optional double initial_primal_weight = 8;
351 message PresolveOptions {
352 // If true runs Glop's presolver on the given instance prior to solving.
353 // Note that convergence criteria are still interpreted with respect to the
354 // original problem. Certificates may not be available if presolve detects
355 // infeasibility. Glop's presolver cannot apply to problems with quadratic
356 // objectives or problems with more than 2^31 variables or constraints. It's
357 // often beneficial to enable the presolver, especially on medium-sized
358 // problems. At some larger scales, the presolver can become a serial
360 optional bool use_glop = 1;
362 // Parameters to control glop's presolver. Only used when use_glop is true.
363 // These are merged with and override PDLP's defaults.
364 optional operations_research.glop.GlopParameters glop_parameters = 2;
366 optional PresolveOptions presolve_options = 16;
368 // Number of L_infinity Ruiz rescaling iterations to apply to the constraint
369 // matrix. Zero disables this rescaling pass. Recommended values to try when
370 // tuning are 0, 5, and 10.
371 optional int32 l_inf_ruiz_iterations = 9 [default = 5];
373 // If true, applies L_2 norm rescaling after the Ruiz rescaling. Heuristically
374 // this has been found to help convergence.
375 optional bool l2_norm_rescaling = 10 [default = true];
377 // For ADAPTIVE_HEURISTIC and ADAPTIVE_DISTANCE_BASED only: A relative
378 // reduction in the potential function by this amount always triggers a
379 // restart. Must be between 0.0 and 1.0.
380 optional double sufficient_reduction_for_restart = 11 [default = 0.1];
382 // For ADAPTIVE_HEURISTIC only: A relative reduction in the potential function
383 // by this amount triggers a restart if, additionally, the quality of the
384 // iterates appears to be getting worse. The value must be in the interval
385 // [sufficient_reduction_for_restart, 1). Smaller values make restarts less
386 // frequent, and larger values make them more frequent.
387 optional double necessary_reduction_for_restart = 17 [default = 0.9];
389 // Linesearch rule applied at each major iteration.
390 optional LinesearchRule linesearch_rule = 12
391 [default = ADAPTIVE_LINESEARCH_RULE];
393 optional AdaptiveLinesearchParams adaptive_linesearch_parameters = 18;
395 optional MalitskyPockParams malitsky_pock_parameters = 19;
397 // Scaling factor applied to the initial step size (all step sizes if
398 // linesearch_rule == CONSTANT_STEP_SIZE_RULE).
399 optional double initial_step_size_scaling = 25 [default = 1.0];
401 // Seeds for generating (pseudo-)random projections of iterates during
402 // termination checks. For each seed, the projection of the primal and dual
403 // solutions onto random planes in primal and dual space will be computed and
404 // added the IterationStats if record_iteration_stats is true. The random
405 // planes generated will be determined by the seeds, the primal and dual
406 // dimensions, and num_threads.
407 repeated int32 random_projection_seeds = 28 [packed = true];
409 // Constraint bounds with absolute value at least this threshold are replaced
411 // NOTE: This primarily affects the relative convergence criteria. A smaller
412 // value makes the relative convergence criteria stronger. It also affects the
413 // problem statistics LOG()ed at the start of the run, and the default initial
414 // primal weight, since that is based on the norm of the bounds.
415 optional double infinite_constraint_bound_threshold = 22 [default = inf];
418 // https://developers.google.com/optimization/lp/pdlp_math#treating_some_variable_bounds_as_infinite
419 // for a description of this flag.
420 optional bool handle_some_primal_gradients_on_finite_bounds_as_residuals = 29
423 // When solving QPs with diagonal objective matrices, this option can be
424 // turned on to enable an experimental solver that avoids linearization of the
425 // quadratic term. The `diagonal_qp_solver_accuracy` parameter controls the
427 // TODO(user): Turn this option on by default for quadratic
428 // programs after numerical evaluation.
429 optional bool use_diagonal_qp_trust_region_solver = 23 [default = false];
431 // The solve tolerance of the experimental trust region solver for diagonal
432 // QPs, controlling the accuracy of binary search over a one-dimensional
433 // scaling parameter. Smaller values imply smaller relative error of the final
435 // TODO(user): Find an expression for the final relative error.
436 optional double diagonal_qp_trust_region_solver_tolerance = 24
439 // If true, periodically runs feasibility polishing, which attempts to move
440 // from latest average iterate to one that is closer to feasibility (i.e., has
441 // smaller primal and dual residuals) while probably increasing the objective
442 // gap. This is useful primarily when the feasibility tolerances are fairly
443 // tight and the objective gap tolerance is somewhat looser. Note that this
444 // does not change the termination criteria, but rather can help achieve the
445 // termination criteria more quickly when the objective gap is not as
446 // important as feasibility.
448 // `use_feasibility_polishing` cannot be used with glop presolve, and requires
449 // `handle_some_primal_gradients_on_finite_bounds_as_residuals == false`.
450 // `use_feasibility_polishing` can only be used with linear programs.
452 // Feasibility polishing runs two separate phases, primal feasibility and dual
453 // feasibility. The primal feasibility phase runs PDHG on the primal
454 // feasibility problem (obtained by changing the objective vector to all
455 // zeros), using the average primal iterate and zero dual (which is optimal
456 // for the primal feasibility problem) as the initial solution. The dual
457 // feasibility phase runs PDHG on the dual feasibility problem (obtained by
458 // changing all finite variable and constraint bounds to zero), using the
459 // average dual iterate and zero primal (which is optimal for the dual
460 // feasibility problem) as the initial solution. The primal solution from the
461 // primal feasibility phase and dual solution from the dual feasibility phase
462 // are then combined (forming a solution of type
463 // `POINT_TYPE_FEASIBILITY_POLISHING_SOLUTION`) and checked against the
464 // termination criteria.
466 optional bool use_feasibility_polishing = 30 [default = false];
468 reserved 13, 14, 15, 20, 21;