Search

The implicit trace estimation problem asks for an approximation of the trace of a square matrix, accessed via matrix-vector products (matvecs). This paper designs new randomized algorithms, XTrace and XNysTrace, for the trace estimation problem by exploiting both variance reduction and the exchangeability principle. For a fixed budget of matvecs, numerical experiments show that the new methods can achieve errors that are orders of magnitude smaller than existing algorithms, such as the Girard--Hutchinson estimator or the Hutch++ estimator. A theoretical analysis confirms the benefits by offering a precise description of the performance of these algorithms as a function of the spectrum of the input matrix. The paper also develops an exchangeable estimator, XDiag, for approximating the diagonal of a square matrix using matvecs. [ABSTRACT FROM AUTHOR]

This is the second paper in a series of two. The goal of the series is to give a polynomial-time algorithm for the 4-coloring problem and the 4-precoloring extension problem restricted to the class of graphs with no induced six-vertex path, thus proving a conjecture of Huang. Combined with previously known results, this completes the classification of the complexity of the 4-coloring problem for graphs with a connected forbidden induced subgraph. In this paper we give a polynomial time-algorithm that starts with a 4-precoloring of a graph with no induced six-vertex path and outputs a polynomial-sized collection of so-called excellent precolorings. Excellent precolorings are easier to handle than general ones, and, in addition, in order to determine whether the initial precoloring can be extended to the whole graph, it is enough to answer the same question for each of the excellent precolorings in the collection. The first paper in the series deals with excellent precolorings, thus providing a complete solution to the problem. [ABSTRACT FROM AUTHOR]

This is the first paper in a series whose goal is to give a polynomial-time algorithm for the 4-coloring problem and the 4-precoloring extension problem restricted to the class of graphs with no induced six-vertex path, thus proving a conjecture of Huang. Combined with previously known results, this completes the classification of the complexity of the 4-coloring problem for graphs with a connected forbidden induced subgraph. In this paper we give a polynomial-time algorithm that determines if a special kind of precoloring of a P6-free graph has a precoloring extension, and constructs such an extension if one exists. Combined with the main result of the second paper of the series, this gives a complete solution to the problem. [ABSTRACT FROM AUTHOR]

In the first part of this paper we develop a theory for image restoration with a learned regularizer that is analogous to that of Meyer's geometric characterization of solutions of the classical variational method of Rudin-Osher-Fatemi (ROF). The learned regularizer we use is a Kantorovich potential for an optimal transport problem of mapping a distribution of noisy images onto clean ones, as first proposed by Lunz, \" Oktem, and Sch\"onlieb. We show that the effect of their restoration method on the distribution of the images is an explicit Euler discretization of a gradient flow on probability space, while our variational problem, dubbed Wasserstein ROF (WROF), is the corresponding implicit Euler discretization. We obtain our geometric characterization of the solution in this learned regularizer setting by first proving a much more general convex analysis theorem for variational problems having solutions characterized by projections. We then use optimal transport arguments to obtain the corresponding theorem for WROF from this general result, as well as a natural decomposition of a transport map into large scale "features" and small scale "details," where scale refers to the magnitude of the transport distance. In the second part of the paper we leverage our theory for restoration with learned regularizers to analyze two algorithms which iterate WROF. We refer to these as iterative regularization and multiscale transport. For the former we obtain a proof of convergence to the clean data. For the latter we produce successive approximations to the target distribution that match it up to finer and finer scales. These two algorithms are in complete analogy to well-known effective methods based on ROF for iterative denoising, respectively hierarchical image decomposition. We also obtain an analogue of the Tadmor-Nezzar-Vese energy identity, which decomposes the Wasserstein 2 distance between two measures into a sum of nonnegative terms that correspond to transport costs at different scales. [ABSTRACT FROM AUTHOR]

This paper presents an interior point method for pure state and mixed-constraint optimal control problems for dynamics, mixed constraints, and cost function all affine in the control variable. This method relies on resolving a sequence of two-point boundary value problems of differential and algebraic equations. This paper establishes a convergence result for primal and dual variables of the optimal control problem. A primal and a primal-dual solving algorithm are presented, and a challenging numerical example is treated for illustration. [ABSTRACT FROM AUTHOR]

In this paper we develop and analyze a variational data assimilation method with efficient decoupled iterative numerical algorithms for the Stokes–Darcy equations with the Beavers–Joseph interface condition. By using Tikhonov regularization and formulating the variational data assimilation into an optimization problem, we establish the existence, uniqueness, and stability of the optimal solution. Based on the weak formulation of the Stokes–Darcy equations, the Lagrange multiplier rule is utilized to derive the first order optimality system for both the continuous and discrete variational data assimilation problems, where the discrete data assimilation is based on a finite element discretization in space and the backward Euler scheme in time. By rescaling the optimality system and then analyzing its corresponding bilinear forms, we prove the optimal finite element convergence rate with special attention paid to recovering uncertainties missed in the optimality system. To solve the discrete optimality system efficiently, three decoupled iterative algorithms are proposed to address the computational cost for both well-conditioned and ill-conditioned variational data assimilation problems, respectively. Finally, numerical results are provided to validate the proposed methods. [ABSTRACT FROM AUTHOR]

In this paper we consider finding a second-order stationary point (SOSP) of nonconvex equality constrained optimization when a nearly feasible point is known. In particular, we first propose a new Newton-conjugate gradient (Newton-CG) method for finding an approximate SOSP of unconstrained optimization and show that it enjoys a substantially better complexity than the Newton-CG method in [C. W. Royer, M. O'Neill, and S. J. Wright, Math. Program., 180 (2020), pp. 451-488]. We then propose a Newton-CG based augmented Lagrangian (AL) method for finding an approximate SOSP of nonconvex equality constrained optimization, in which the proposed Newton-CG method is used as a subproblem solver. We show that under a generalized linear independence constraint qualification (GLICQ), our AL method enjoys a total inner iteration complexity of O(ε7/peration complexity of O(ε7/2 min\{ n, ε3/4\}) for finding an (ε, ε)-SOSP of nonconvex equality constrained optimization with high probability, which are significantly better than the ones achieved by the proximal AL method in [Y. Xie and S. J. Wright, J. Sci. Comput., 86 (2021), pp. 1-30]. In addition, we show that it has a total inner iteration complexity of O(ε11/2) and an operation complexity of O(ε11/2 min\{ n, ε5/4\}) when the GLICQ does not hold. To the best of our knowledge, all the complexity results obtained in this paper are new for finding an approximate SOSP of nonconvex equality constrained optimization with high probability. Preliminary numerical results also demonstrate the superiority of our proposed methods over the other competing algorithms. [ABSTRACT FROM AUTHOR]

The Euler elastica (EE) model with surface curvature can generate artifact-free results compared with the traditional total variation regularization model in image processing. However, strong nonlinearity and singularity due to the curvature term in the EE model pose a great challenge for one to design fast and stable algorithms for the EE model. In this paper, we propose a new, fast, hybrid alternating minimization (HALM) algorithm for the EE model based on a bilinear decomposition of the gradient of the underlying image and prove the global convergence of the minimizing sequence generated by the algorithm under mild conditions. The HALM algorithm comprises three subminimization problems and each is either solved in the closed form or approximated by fast solvers, making the new algorithm highly accurate and efficient. We also discuss the extension of the HALM strategy to deal with general curvature-based variational models, especially with a Lipschitz smooth functional of the curvature. A host of numerical experiments are conducted to show that the new algorithm produces good results with much-improved efficiency compared to other state-of-the-art algorithms for the EE model. As one of the benchmarks, we show that the average running time of the HALM algorithm is at most one-quarter of that of the fast operator-splitting-based Deng--Glowinski--Tai algorithm. [ABSTRACT FROM AUTHOR]

The constrained pressure residual (CPR) algorithm is a family of well-known and industry-standard preconditioners for large-scale reservoir simulation. The CPR algorithm is a twostage preconditioner to deal with different blocks stage-by-stage, and is often able to effectively improve the robustness behavior and the convergence speed of linear iterations. Nonetheless, as a linear preconditioner, it is hard for the traditional CPR method to be capable of action at the nonlinear level for effectively solving large sparse nonlinear systems of equations with high nonlinearity. In this paper, we present and study the extension of this linear method to the nonlinearly CPR (NCPR) case for solving fractured reservoir problems, to deal with the difficulty of the slow convergence or stagnation from the nonlinear level. In the proposed nonlinear preconditioning, a subspace nonlinear block system is first built and solved to remove the unbalanced nonlinearities of the pressure and the saturation fields, and the fast convergence can then be restored when a variant of semismooth Newton methods is called after the subspace nonlinear block system is solved. Experiments on two or three dimensional porous media applications are presented to demonstrate the applicability and parallel scalability of the aforementioned NCPR method. We also show that the proposed algorithm is superior to the commonly used nonlinear algorithm in terms of the robustness and the positivity-preserving property. [ABSTRACT FROM AUTHOR]

Proximal mappings are essential in splitting algorithms for both convex and noncon-vex optimization. In this paper, we show that proximal mappings of every prox-bounded function are nonexpansive if and only if they are firmly nonexpansive if and only if they are averaged if and only if the function is convex. Lipschitz proximal mappings of prox-bounded functions are also characterized via hypoconvex or strongly convex functions. Our results generalize a recent result due to Rockafellar [ABSTRACT FROM AUTHOR]

Given closed convex sets Ci, i=1, ... l, and some nonzero linear maps Ai, i=1, ... l, of suitable dimensions, the multiset split feasibility problem aims at finding a point in ... based on computing projections onto Ci and multiplications by Ai and AiT. In this paper, we consider the associated best approximation problem, i.e., the problem of computing projections onto ...; we refer to this problem as the best approximation problem in multiset split feasibility settings (BA-MSF). We adapt the Dykstra's projection algorithm, which is classical for solving the BA-MSF in the special case when all Ai - I, to solve the general BA-MSF. Our Dykstra-type projection algorithm is derived by applying (proximal) coordinate gradient descent to the Lagrange dual problem, and it only requires computing projections onto Ci and multiplications by Ai and AiT in each iteration. Under a standard relative interior condition and a genericity assumption on the point we need to project, we show that the dual objective satisfies the Kurdyka-Łojasiewicz property with an explicitly computable exponent on a neighborhood of the (typically unbounded) dual solution set when each Ci is C1,α.-cone reducible for some α∈ (0,1]: this class of sets covers the class of C²-cone reducible sets, which include all polyhedrons, second-order cone, and the cone of positive semidefinite matrices as special cases. Using this, explicit convergence rate (linear or sublinear) of the sequence generated by the Dykstra-type projection algorithm is derived. Concrete examples are constructed to illustrate the necessity of some of our assumptions. [ABSTRACT FROM AUTHOR]

This paper proposes a path-based approach for the minimization of a continuously differentiable function over sparse symmetric sets, which is a hard problem that exhibits a restrictiveness-hierarchy of necessary optimality conditions. To achieve the more restrictive conditions in the hierarchy, state-of-the-art algorithms require a support optimization oracle that must exactly solve the problem in smaller dimensions. The path-based approach developed in this study produces a path-based optimality condition, which is placed well in the restrictiveness-hierarchy, and a method to achieve it that does not require a support optimization oracle and, moreover, is projection-free. In the development process, new results are derived for the regularized linear minimization problem over sparse symmetric sets, which give additional means to identify optimal solutions for convex and concave objective functions. We complement our results with numerical examples. [ABSTRACT FROM AUTHOR]

In this paper, we consider a class of structured nonsmooth difference-of-convex (DC) minimization in which the first convex component is the sum of a smooth and a nonsmooth function, while the second convex component is the supremum of finitely many convex smooth functions. The existing methods for this problem usually have weak convergence guarantees or need to solve lots of subproblems per iteration. Due to this, we propose hybrid algorithms for solving this problem in which we first compute approximate critical points and then check whether these points are approximate D-stationary points. Under suitable conditions, we prove that there exists a subsequence of iterates of which every accumulation point is a D-stationary point. Some preliminary numerical experiments are conducted to demonstrate the efficiency of the proposed algorithms. [ABSTRACT FROM AUTHOR]

In this paper, we have studied a decomposition method for solving a class of non-convex two-stage stochastic programs, where both the objective and constraints of the second-stage problem are nonlinearly parameterized by the first-stage variables. Due to the failure of the Clarke regularity of the resulting nonconvex recourse function, classical decomposition approaches such as Benders decomposition and (augmented) Lagrangian-based algorithms cannot be directly generalized to solve such models. By exploring an implicitly convex-concave structure of the recourse function, we introduce a novel decomposition framework based on the so-called partial Moreau envelope. The algorithm successively generates strongly convex quadratic approximations of the recourse function based on the solutions of the second-stage convex subproblems and adds them to the first-stage mas-ter problem. Convergence has been established for both a fixed number of scenarios and a sequential internal sampling strategy. Numerical experiments are conducted to demonstrate the effectiveness of the proposed algorithm. [ABSTRACT FROM AUTHOR]

A lattice of integers is the collection of all linear combinations of a set of vectors for which all entries of the vectors are integers and all coefficients in the linear combinations are also integers. Lattice reduction refers to the problem of finding a set of vectors in a given lattice such that the collection of all integer linear combinations of this subset is still the entire original lattice and so that the Euclidean norms of the subset are reduced. The present paper proposes simple, efficient iterations for lattice reduction which are guaranteed to reduce the Euclidean norms of the basis vectors (the vectors in the subset) monotonically during every iteration. Each iteration selects the basis vector for which projecting off (with integer coefficients) the components of the other basis vectors along the selected vector minimizes the Euclidean norms of the reduced basis vectors. Each iteration projects off the components along the selected basis vector and efficiently updates all information required for the next iteration to select its best basis vector and perform the associated projections. [ABSTRACT FROM AUTHOR]

Structure-preserving doubling algorithms (SDAs) are efficient algorithms for solving Riccati-type matrix equations. However, breakdowns may occur in SDAs. To remedy this drawback, in this paper, we first introduce Ω -symplectic forms (Ω -SFs), consisting of symplectic matrix pairs with a Hermitian parametric matrix Ω. Based on Ω -SFs, we develop modified SDAs (MSDAs) for solving the associated Riccati-type equations. MSDAs generate sequences of symplectic matrix pairs in Ω -SFs and prevent breakdowns by employing a reasonably selected Hermitian matrix Ω. In practical implementations, we show that the Hermitian matrix Ω in MSDAs can be chosen as a real diagonal matrix that can reduce the computational complexity. The numerical results demonstrate a significant improvement in the accuracy of the solutions by MSDAs. [ABSTRACT FROM AUTHOR]

In the Bin Packing problem one is given n items with weights w(1),. . . ,w(n) and m bins with capacities c1, . . ., cm. The goal is to partition the items into sets S1, . . . ,Sm such that w(Sj 0 such that an instance of Bin Packing with m bins can be solved in O (2(1 y m)n) randomized time. Before our work, such improved algorithms were not known even for m= 4. A key step in our approach is the following new result in Littlewood--Offord theory on the additive combinatorics of subset Σ: For every γ > 0 there exists an >0 such that if | { X n{ 1, . . . ,n}: w(X) =v} | 2(1)n for some v, then { w(X) :X { 1, . . . ,n\} | &#947n. [ABSTRACT FROM AUTHOR]

Semidefinite programs (SDPs) arising from the Kalman—Yakubovich—Popov (KYP) lemma are frequently encountered in systems robustness analysis, filter design, and other control/signal processing related applications. These programs possess a special structure that can be exploited to construct specialized algorithms that substantially outperform general-purpose SDP solvers. In this paper, a new interior path-following algorithm that utilizes this structure is proposed. The main idea behind the algorithm is a new barrier function for these specially structured SDPs. Convergence of the new algorithm is shown and a measure of the accuracy of suboptimal solutions produced by the algorithm is provided. The algorithm is tested in numerical experiments and the results indicate that the new algorithm is indeed favorable against general-purpose SDP solvers in many circumstances. [ABSTRACT FROM AUTHOR]

Stochastic Parareal (SParareal) is a probabilistic variant of the popular parallel-in-time algorithm known as Parareal. Similarly to Parareal, it combines fine- and coarse-grained solutions to an ODE using a predictor-corrector (PC) scheme. The key difference is that carefully chosen random perturbations are added to the PC to try to accelerate the location of a stochastic solution to the ODE. In this paper, we derive superlinear and linear mean-square error bounds for SParareal applied to nonlinear systems of ODEs using different types of perturbations. We illustrate these bounds numerically on a linear system of ODEs and a scalar nonlinear ODE, showing a good match between theory and numerics. [ABSTRACT FROM AUTHOR]

We propose in this paper a new minimization algorithm based on a slightly modified version of the scalar auxialiary variable (SAV) approach coupled with a relaxation step and an adaptive strategy. It enjoys several distinct advantages over popular gradient based methods: (i) it is unconditionally energy diminishing with a modified energy which is intrinsically related to the original energy, and thus no parameter tuning is needed for stability; (ii) it allows the use of large step-sizes which can effectively accelerate the convergence rate. We also present a convergence analysis for some SAV based algorithms, which include our new algorithm without the relaxation step as a special case. We apply our new algorithm to several illustrative and benchmark problems and compare its performance with several popular gradient based methods. The numerical results indicate that the new algorithm is very robust, and its adaptive version usually converges significantly faster than those popular gradient decent based methods. [ABSTRACT FROM AUTHOR]

Policy gradient algorithms have been widely applied to Markov decision processes and reinforcement learning problems in recent years. Regularization with various entropy functions is often used to encourage exploration and improve stability. This paper proposes an approximate Newton method for the policy gradient algorithm with entropy regularization. In the case of Shannon entropy, the resulting algorithm reproduces the natural policy gradient algorithm. For other entropy functions, this method results in brand-new policy gradient algorithms. We prove that all these algorithms enjoy Newton-type quadratic convergence and that the corresponding gradient flow converges globally to the optimal solution. We use synthetic and industrial-scale examples to demonstrate that the proposed approximate Newton method typically converges in single-digit iterations, often orders of magnitude faster than other state-of-the-art algorithms. [ABSTRACT FROM AUTHOR]

This paper presents an extension of the Adaptive Antoulas-Anderson (AAA) algorithm for rational modeling. Specifically, our new stable multi-input multi-output AAA (smiAAA) algorithm builds rational approximations of multi-input signals with a common set of stable poles. A new methodology is presented for iteratively enforcing the stability constraint on the poles. We demonstrate the strengths of this approach compared to the stability enforcement in the FastAAA algorithm. Results using the smiAAA algorithm are compared with the commonly used Vector Fitting algorithm and the more recently published Rational Krylov fitting (RKFIT) algorithm. The Vector Fitting and RKFIT algorithms both require the user to input the number of poles to use in the approximations. If the final approximation is not accurate enough, the user must restart the Vector Fitting or RKFIT algorithm with a larger number of poles and/or a new starting location for the poles. In contrast, the smiAAA algorithm is designed to allow the user to simply input the desired accuracy of the approximations, and the necessary number of poles is detected automatically. This permits users to produce approximations of a desired accuracy with no knowledge about the underlying order of the system being approximated, preventing the algorithm from ever needing to be rerun. An additional feature for preventing extraneous poles from being returned by the AAA algorithm is also discussed. The cause of these extraneous poles is efficiently detected and removed by our presented methodology. The examples presented demonstrate that the smiAAA algorithm can efficiently produce approximations of similar or better accuracy than the Vector Fitting and RKFIT algorithms while requiring less input from the user. [ABSTRACT FROM AUTHOR]

In this paper we consider large-scale composite optimization problems having the objective function formed as a sum of two terms (possibly nonconvex); one has a (block) coordinatewise Lipschitz continuous gradient and the other is differentiable but nonseparable. Under these general settings we derive and analyze two new coordinate descent methods. The first algorithm, referred to as the coordinate proximal gradient method, considers the composite form of the objective function, while the other algorithm disregards the composite form of the objective and uses the partial gradient of the full objective, yielding a coordinate gradient descent scheme with novel adaptive stepsize rules. We prove that these new stepsize rules make the coordinate gradient scheme a descent method, provided that additional assumptions hold for the second term in the objective function. We present a complete worst-case complexity analysis for these two new methods in both convex and nonconvex settings, provided that the (block) coordinates are chosen random or cyclic. Preliminary numerical results also confirm the efficiency of our two algorithms for practical problems. [ABSTRACT FROM AUTHOR]

An adaptive regularization algorithm for unconstrained nonconvex optimization is presented in which the objective function is never evaluated but only derivatives are used. This algorithm belongs to the class of adaptive regularization methods, for which optimal worst-case complexity results are known for the standard framework where the objective function is evaluated. It is shown in this paper that these excellent complexity bounds are also valid for the new algorithm despite the fact that significantly less information is used. In particular, it is shown that if derivatives of degree one to p are used, the algorithm will find an ε 1-approximate first-order minimizer in at most O(ε (p+1)/p 1 ) iterations and an (ε 1, ε 2)-approximate second-order minimizer in at most O(max[ε¹ (p+1)/p 1, ε (p+1)/(p 1) 2 ]) iterations. As a special case, the new algorithm using first and second derivatives, when applied to functions with Lipschitz continuous Hessian, will find an iterate xk at which the gradient's norm is less than ε 1 in at most O(ε1 3/2) iterations. [ABSTRACT FROM AUTHOR]

Nonlinear filtering is fundamental to many engineering problems, as it involves inferring the state of a system given complicated dynamics and noisy observations. Some approaches to nonlinear filtering use the analysis of the underlying PDE or stochastic PDE governing the evolution of the posterior probability distribution over time, one approach, in particular, being the Yau-Yau method. In this paper, we give a maximum a posteriori (MAP) framework for the Yau-Yau method. Furthermore, we propose convex filtering, intermediate between linear and nonlinear filtering, which gives criteria under which the posterior preserves log-concavity. The key tool from the PDE is a result from Korevaar, giving criteria under which a quasilinear parabolic PDE preserves convexity. A bridge between the MAP estimator and the structure of the posterior is explained. Finally, we propose a novel numerical method based on iteration and apply this method to a tracking problem. [ABSTRACT FROM AUTHOR]

Motivated by an application in cancer genomics, Hajirasouliha and Raphael [Proceedings of the 14th International Workshop on Algorithms in Bioinformatics, 2014, pp. 354-367] proposed the split-row problem (SR). In this problem, an m\times n binary matrix M is given. A split-row operation on M is defined as replacing a row r by k > 1 rows r(1), r(2), . . ., r(k) whose bitwise OR is equal to r. The cost of the operation is the number of additional rows induced, that is, k 1. The goal is to find a sequence of split-row operations that transforms M into a matrix corresponding to a perfect phylogeny and the total cost is minimized. Recently, Hujdurovi\'c et al. [ACM Trans. Algorithms, 14 (2018), 26] proved the APX-hardness of SR and presented efficient exact and approximation algorithms. The parameterized study of SR was left as a direction for future work. Let\varepsilon (M) denote the minimum total cost. This paper gives an O*(2min(n,2x(M))-time exact algorithm for SR. This result indicates that SR is fixed-parameter tractable when parameterized by E (M). In addition, in the worst case, our algorithm requires O* (2n) time, significantly improving the previous upper bound of O\ast (nn). Hujdurovi\'c et al.'s exact algorithm can be modified to solve a variant of SR, called the distinct split-row problem (DSR). Our algorithm can be adapted to this variant as well. In addition, our algorithms can be extended to solve SR and DSR with the following additional constraint: only the rows in a given subset are allowed to be split. [ABSTRACT FROM AUTHOR]

For many algorithmic problems on graphs of treewidth t, a standard dynamic programming approach gives algorithms with time and space complexity 2O(t)\cdot nO(1). It turns out that when one considers the more restrictive parameter treedepth, it is often the case that a variation of this technique can be used to reduce the space complexity to polynomial, while retaining time complexity of the form 2O(d)\cdot nO(1), where d is the treedepth. This transfer of methodology is, however, far from automatic. For instance, for problems with connectivity constraints, standard dynamic programming techniques give algorithms with time and space complexity 2O(t)...nO(1) on graphs of treewidth t, but it is not clear how to convert them into time-efficient polynomial space algorithms for graphs of low treedepth. Cygan et al. [ACM Trans. Algorithms, 18 (2022), 17] introduced the Cut\&Count technique and showed that a certain class of problems with connectivity constraints can be solved in time and space complexity 2O(t) · nO(1). Recently, Hegerfeld and Kratsch (STACS'20) showed that, for some of those problems, the Cut\&Count technique can be also applied in the setting of treedepth, and it gives algorithms with running time 2O(d) · nO(1) and polynomial space usage. However, several important problems eluded such a treatment, with the most prominent examples being Hamiltonian Cycle and Longest Path. In this paper, we clarify the situation by showing that Hamiltonian cycle, Hamiltonian Path, Long Cycle, Long Path, and Min Cycle Cover all admit 5d · nO(1)-time and polynomial space algorithms on graphs of treedepth d. The algorithms are randomized Monte Carlo with only false negatives. [ABSTRACT FROM AUTHOR]

In this paper, we present the new "asynchronous truncated multigrid-reduction-intime"" (AT-MGRIT) algorithm for introducing time parallelism to the solution of discretized timedependent problems. The new algorithm is based on the multigrid-reduction-in-time (MGRIT) approach, which, in certain settings, is equivalent to another common multilevel parallel-in-time method, Parareal. In contrast to Parareal and MGRIT that both consider a global temporal grid over the entire time interval on the coarsest level, the AT-MGRIT algorithm uses truncated local time grids on the coarsest level, each grid covering certain temporal subintervals. These local grids can be solved completely in an independent way from each other, which reduces the sequential part of the algorithm and, thus, increases parallelism in the method. Here, we study the effect of using truncated local coarse grids on the convergence of the algorithm, both theoretically and numerically, and show, using challenging nonlinear problems, that the new algorithm consistently outperforms classical Parareal/MGRIT in terms of time to solution. [ABSTRACT FROM AUTHOR]

The tensor-train (TT) decomposition expresses a tensor in a data-sparse format used in molecular simulations, high-order correlation functions, and optimization. In this paper, we propose four parallelizable algorithms that compute the TT format from various tensor inputs: (1) Parallel-TTSVD for traditional format, (2) PSTT and its variants for streaming data, (3) Tucker2TT for Tucker format, and (4) TT-fADI for solutions of Sylvester tensor equations. We provide theoretical guarantees of accuracy, parallelization methods, scaling analysis, and numerical results. For example, for a d-dimension tensor in Rn× ... × n, a two-sided sketching algorithm PSTT2 is shown to have a memory complexity of O (n d/2), improving upon O (nd 1) from previous algorithms. [ABSTRACT FROM AUTHOR]

In this paper, an efficient ensemble domain decomposition algorithm is proposed for fast solving the fully mixed random Stokes--Darcy model with the physically realistic Beavers--Joseph interface conditions. We utilize the Monte Carlo method for the coupled model with random inputs to derive some deterministic Stokes--Darcy numerical models and use the idea of the ensemble to realize the fast computation of multiple problems. One remarkable feature of the algorithm is that multiple linear systems share a common coefficient matrix in each deterministic numerical model, which significantly reduces the computational cost and achieves comparable accuracy with the traditional methods. Moreover, by domain decomposition, we can decouple the Stokes--Darcy system into two smaller subphysics problems naturally. Both mesh-dependent and mesh-independent convergence rates of the algorithm are rigorously derived by choosing suitable Robin parameters. Optimized Robin parameters are derived and analyzed to accelerate the convergence of the proposed algorithm. Especially, for small hydraulic conductivity in practice, the almost optimal geometric convergence can be obtained by finite element discretization. Finally, two groups of numerical experiments are conducted to validate the exclusive features of the proposed algorithm. [ABSTRACT FROM AUTHOR]

We study exact algorithms for Metric TSP in Rd. In the early 1990s, algorithms with nO(\surd n) running time were presented for the planar case, and some years later an algorithm with nO(n1 1/d) running time was presented for any d \geqslant 2. Despite significant interest in subexponential exact algorithms over the past decade, there has been no progress on Metric TSP, except for a lower bound stating that the problem admits no 2o(n1 1/d) algorithm unless ETH fails. In this paper we settle the complexity of Metric TSP, up to constant factors in the exponent and under ETH, by giving an algorithm with running time 2O(n1 1/d). [ABSTRACT FROM AUTHOR]

In this paper we bring together some of the key ideas and methods of two disparate fields of mathematical research, frame theory, and optimal transport, using the methods of the second to answer questions posed in the first. In particular, we construct gradient flows in the Wasserstein space P2(\BbbR d) for a new potential, the tightness potential, which is a modification of the probabilistic frame potential. It is shown that the potential is suited for the application of a gradient descent scheme from optimal transport that can be used as the basis of an algorithm to evolve an existing frame toward a tight probabilistic frame. [ABSTRACT FROM AUTHOR]

The paper deals with stochastic difference-of-convex-functions (DC) programs, that is, optimization problems whose cost function is a sum of a lower semicontinuous DC function and the expectation of a stochastic DC function with respect to a probability distribution. This class of nonsmooth and nonconvex stochastic optimization problems plays a central role in many practical applications. Although there are many contributions in the context of convex and/or smooth stochastic optimization, algorithms dealing with nonconvex and nonsmooth programs remain rare. In deterministic optimization literature, the DC algorithm (DCA) is recognized to be one of the few algorithms able to effectively solve nonconvex and nonsmooth optimization problems. The main purpose of this paper is to present some new stochastic DCAs for solving stochastic DC programs. The convergence analysis of the proposed algorithms is carefully studied, and numerical experiments are conducted to justify the algorithms' behaviors. [ABSTRACT FROM AUTHOR]

Packing and covering semidefinite programs (SDPs) appear in natural relaxations of many combinatorial optimization problems as well as a number of other applications. Recently, several techniques were proposed, which utilize the particular structure of this class of problems, to obtain more efficient algorithms than those offered by general SDP solvers. For certain applications, such as those described in this paper, it may be desirable to obtain sparse dual solutions, i.e., those with support size (almost) independent of the number of primal constraints. In this paper, we give an algorithm that finds such solutions, which is an extension of a logarithmic-potential based algorithm of Grigoriadis et al. [SIAM J. Optim. 11 (2001), pp. 1081-1091] for packing/covering linear programs. [ABSTRACT FROM AUTHOR]

We study the random reshuffling (RR) method for smooth nonconvex optimization problems with a finite-sum structure. Though this method is widely utilized in practice, e.g., in the training of neural networks, its convergence behavior is only understood in several limited settings. In this paper, under the well-known Kurdyka-Łojasiewicz (KL) inequality, we establish strong limit-point convergence results for RR with appropriate diminishing step sizes; namely, the whole sequence of iterates generated by RR is convergent and converges to a single stationary point in an almost sure sense. In addition, we derive the corresponding rate of convergence, depending on the KL exponent and suitably selected diminishing step sizes. When the KL exponent lies in [0,1/2], the convergence is at a rate of O(t-1) with t counting the number of iterations. When the KL exponent belongs to (1/2,1), our derived convergence rate is of the form O(t-q) with q∈(0,1) depending on the KL exponent. The standard KL inequality-based convergence analysis framework only applies to algorithms with a certain descent property. We conduct a novel convergence analysis for the nondescent RR method with diminishing step sizes based on the KL inequality, which generalizes the standard KL framework. We summarize our main steps and core ideas in an informal analysis framework, which is of independent interest. As a direct application of this framework, we establish similar strong limit-point convergence results for the reshuffled proximal point method. [ABSTRACT FROM AUTHOR]

Finding a zero of a sum of maximally monotone operators is a fundamental problem in modern optimization and nonsmooth analysis. Assuming that the resolvents of the operators are available, this problem can be tackled with the Douglas-Rachford algorithm. However, when dealing with three or more operators, one must work in a product space with as many factors as there are operators. In groundbreaking recent work by Ryu and by Malitsky and Tam, it was shown that the number of factors can be reduced by one. A similar reduction was achieved recently by Campoy through a clever reformulation originally proposed by Kruger. All three splitting methods guarantee weak convergence to some solution of the underlying sum problem; strong convergence holds in the presence of uniform monotonicity. In this paper, we provide a case study when the operators involved are normal cone operators of subspaces and the solution set is thus the intersection of the subspaces. Even though these operators lack strict convexity, we show that striking conclusions are available in this case: strong (instead of weak) convergence and the solution obtained is (not arbitrary but) the projection onto the intersection. To illustrate our results, we also perform numerical experiments. [ABSTRACT FROM AUTHOR]

We propose a gradient-based block coordinate descent (BCD-G) framework to solve the joint approximate diagonalization of matrices defined on the product of the complex Stiefel manifold and the special linear group. Instead of the cyclic fashion, we choose a block optimization based on the Riemannian gradient. To update the first block variable in the complex Stiefel manifold, we use the well-known line search descent method. To update the second block variable in the special linear group, based on four different kinds of elementary transformations, we construct three classes, GLU, GQU and GU, and then get three BCD-G algorithms, BCD-GLU, BCD-GQU and BCDGU. We establish the global convergence and weak convergence of these three algorithms using the Lojasiewicz gradient inequality under the assumption that the iterates are bounded. We also propose a gradient-based Jacobi-type framework to solve the joint approximate diagonalization of matrices defined on the special linear group. As in the BCD-G case, using the GLU and GQU classes of elementary transformations, we focus on the Jacobi-GLU and Jacobi-GQU algorithms and establish their global convergence and weak convergence as well. All the algorithms and convergence results described in this paper also apply to the real case. [ABSTRACT FROM AUTHOR]

For coupling problems such as the Navier-Stokes-Darcy problem in porous and fluid layers, a natural question is that, while the orders of polynomials approximating the variables are different, could the lower order elements pollute the higher ones? In this paper, the numerical pollution of a decoupled algorithm for the unsteady convection model in superposed fluid and porous layers is studied. We first prove the existence and uniqueness of the weak solution to the convection model with Beavers-Joseph-Saffman interface conditions, by using vanishing viscosity methods together with backward-Euler time discretization. Second, an efficient decoupled algorithm of the model based on the Lagrange multiplier method is proposed, and its stability is proved under some restrictions on coefficients and initial values. Then the optimal error estimates are provided revealing that the optimal convergence orders can be achieved without using equal or similar order finite elements for the fluid velocities in fluid and porous layers. This means that there is no numerical pollution between lower and higher order elements when the matched finite element pairs are used. Finally, numerical tests are implemented to verify our theoretical results. [ABSTRACT FROM AUTHOR]

Efficient sampling from a high-dimensional Gaussian distribution is an old but high-stakes issue. Vanilla Cholesky samplers imply a computational cost and memory requirements that can rapidly become prohibitive in high dimensions. To tackle these issues, multiple methods have been proposed from different communities ranging from iterative numerical linear algebra to Markov chain Monte Carlo (MCMC) approaches. Surprisingly, no complete review and comparison of these methods has been conducted. This paper aims to review all these approaches by pointing out their differences, close relations, benefits, and limitations. In addition to reviewing the state of the art, this paper proposes a unifying Gaussian simulation framework by deriving a stochastic counterpart of the celebrated proximal point algorithm in optimization. This framework offers a novel and unifying revisiting of most of the existing MCMC approaches while also extending them. Guidelines to choosing the appropriate Gaussian simulation method for a given sampling problem in high dimensions are proposed and illustrated with numerical examples. [ABSTRACT FROM AUTHOR]

This paper focuses on efficient algorithms for finding the Dantzig selector which was first proposed by Candès and Tao as an effective variable selection technique in the linear regression. This paper first reformulates the Dantzig selector problem as an equivalent convex composite optimization problem and proposes a semismooth Newton augmented Lagrangian (Ssnal) algorithm to solve the equivalent form. This paper also applies a proximal point dual semismooth Newton (PpdSsn) algorithm to solve another equivalent form of the Dantzig selector problem. Comprehensive results on the global convergence and local asymptotic superlinear convergence of the Ssnal and PpdSsn algorithms are characterized under very mild conditions. The computational costs of a semismooth Newton algorithm for solving the subproblems involved in the Ssnal and PpdSsn algorithms can be cheap by fully exploiting the second order sparsity and employing efficient techniques. Numerical experiments on the Dantzig selector problem with synthetic and real data sets demonstrate that the Ssnal and PpdSsn algorithms substantially outperform the state-of-the-art first order algorithms even for the required low accuracy, and the proposed algorithms are able to solve the large-scale problems robustly and efficiently to a relatively high accuracy. [ABSTRACT FROM AUTHOR]

The terminal backup problems ([E. Anshelevich and A. Karagiozova, SIAM J. Comput., 40 (2011), pp. 678--708]) form a class of network design problems: Given an undirected graph with a requirement on terminals, the goal is to find a minimum-cost subgraph satisfying the connectivity requirement. The node-connectivity terminal backup problem requires a terminal to connect other terminals with a number of node-disjoint paths. This problem is not known whether is NP-hard or tractable. Fukunaga [SIAM J. Discrete Math., 30 (2016), pp. 777--800] gave a 4/3-approximation algorithm based on an LP-rounding scheme using a general LP-solver. In this paper, we develop a combinatorial algorithm for the relaxed LP to find a half-integral optimal solution in O(mlog(n-A · MF(/cn,m T/c2n)) time, where n is the number of nodes, m is the number of edges, k is the number of terminals, A is the maximum edge-cost, U is the maximum edge-capacity, and MF(nz,mz) is the time complexity of a max-flow algorithm in a network with nz nodes and m' edges. The algorithm implies that the 4/3-approximation algorithm for the node-connectivity terminal backup problem is also efficiently implemented. For the design of algorithm, we explore a connection between the node-connectivity terminal backup problem and a new type of a multiflow, which is called a separately capacitated multiflow. We show a min-max theorem which extends the Lovasz-Cherkassky theorem to the node-capacity setting. Our results build on discrete convexity in the node-connectivity terminal backup problem. [ABSTRACT FROM AUTHOR]

This paper presents smoothing schemes for obtaining approximate stationary points of unconstrained or linearly constrained composite nonconvex-concave min-max (and hence nonsmooth) problems by applying well-known algorithms to composite smooth approximations of the original problems. More specifically, in the unconstrained (resp., constrained) case, approximate stationary points of the original problem are obtained by applying, to its composite smooth approximation, an accelerated inexact proximal point (resp., quadratic penalty) method presented in a previous paper by the authors. Iteration complexity bounds for both smoothing schemes are also established. Finally, numerical results are given to demonstrate the efficiency of the unconstrained smoothing scheme. [ABSTRACT FROM AUTHOR]

Symmetries of a set are linear transformations that map the set to itself. A fundamental domain of a symmetric set is a subset that contains at least one representative from each of the symmetric equivalence classes (orbits) in the set. This paper contributes a novel polynomial algorithm for constructing minimal polytopic fundamental domains of polytopic sets. Our algorithm is applicable for generic linear symmetries of the set and has linear complexity in the number of facets and dimension of the symmetric polytope, e.g., the feasible region of an optimization problem. In addition, this paper contributes a novel polynomial algorithm for mapping an element of the polytope to its representative in the fundamental domain. The algorithms are demonstrated in four examples--two illustrative and two practical. In the first practical example, we show that a minimal fundamental domain of the hypercube under the symmetric group is the set of points with sorted elements. In the second practical example, we show how the construction algorithm can be applied to the max-cut problem. [ABSTRACT FROM AUTHOR]

The characterization of all the constraint satisfaction problems solvable by local consistency checking (also known as CSPs of bounded width) was proposed by Feder and Vardi [SIAM J. Comput., 28 (1998), pp. 57104]. It was confirmed by two independent proofs by Bulatov [Bounded Relational Width, manuscript, 2009] and Barto and Kozik [L. Barto and M. Kozik, 50th Annual IEEE Symposium on Foundations of Computer Science, 2009, pp. 595 603], [L. Barto and M. Kozik, J. ACM, 61 (2014), 3]. Later Barto [J. Logic Comput., 26 (2014), pp. 923 943] proved a collapse of the hierarchy of local consistency notions by showing that (2; 3) minimality solves all the CSPs of bounded width. In this paper we present a new consistency notion, jpq consistency, which also solves all the CSPs of bounded width. Our notion is strictly weaker than (2; 3) consistency, (2; 3) minimality, path consistency, and singleton arc consistency (SAC). This last fact allows us to answer the question of Chen, Dalmau, and Gruien [J. Logic Comput., 23 (2013), pp. 87 108] by confirming that SAC solves all the CSPs of bounded width. Moreover, as known algorithms work faster for SAC, the result implies that CSPs of bounded width can be, in practice, solved more efficiently. The definition of jpq consistency is closely related to a consistency condition obtained from the rounding of an SDP relaxation of a CSP instance. In fact, the main result of this paper is used by Dalmau et al. [Proceedings of the 28th Annual ACM-SIAM Symposium on Discrete Algorithms, SIAM, Philadelphia, ACM, New York, 2017, pp. 340{357] to show that CSPs with near unanimity polymorphisms admit robust approximation algorithms with polynomial loss. Finally, an algebraic characterization of some term conditions satisfied in algebras associated with templates of bounded width, first proved by Brady, is a direct consequence of our result. [ABSTRACT FROM AUTHOR]

In this paper, we develop new primal-dual algorithms to solve a class of nonsmooth and nonlinear convex-concave minimax problems, which covers many existing and brand-new models as special cases. Our approach relies on a combination of a generalized augmented Lagrangian function, Nesterov's accelerated scheme, and adaptive parameter updating strategies. Our algorithmic framework is single-loop and unifies two important settings: general convex-concave and convex-linear cases. Under mild assumptions, our algorithms achieve O (1/k) convergence rates through three different criteria: primal-dual gap, primal objective residual, and dual objective residual, where k is the iteration counter. Our rates are both ergodic (i.e., on a weighted averaging sequence) and nonergodic (i.e., on the last-iterate sequence). These convergence rates can be boosted up to O (1/k²) if only one objective term is strongly convex (or, equivalently, its conjugate is L-smooth). To the best of our knowledge, this is the first algorithm achieving optimal rates on the primal last-iterate sequence for convex-linear minimax problems. As a byproduct, we specify our algorithms to solve a general convex cone constrained program with both ergodic and nonergodic rate guarantees. We test our algorithms and compare them with two recent methods on two numerical examples. [ABSTRACT FROM AUTHOR]

In this paper, we propose two new solution schemes to solve the stochastic strongly monotone variational inequality (VI) problems: the stochastic extra-point solution scheme and the stochastic extra-momentum solution scheme. The first one is a general scheme based on updating the iterative sequence and an auxiliary extra-point sequence. In the case of a deterministic VI model, this approach includes several state-of-the-art first-order methods as its special cases. The second scheme combines two momentum-based directions: the so-called heavy-ball direction and the optimism direction, where only one projection per iteration is required in its updating process. We show that if the variance of the stochastic oracle is appropriately controlled, then both schemes can be made to achieve optimal iteration complexity of O(κ ln (1/∈)) to reach an ∈ -solution for a strongly monotone VI problem with condition number κ. As a specific application to stochastic VI, we demonstrate how to incorporate a zeroth-order approach for solving stochastic minimax saddle-point problems in our schemes, where only noisy and biased samples of the objective can be obtained, with a total sample complexity of O (κ/∈). [ABSTRACT FROM AUTHOR]

Conjugate gradient methods are important first-order optimization algorithms both in Euclidean spaces and on Riemannian manifolds. However, while various types of conjugate gradient methods have been studied in Euclidean spaces, there are relatively fewer studies for those on Riemannian manifolds (i.e., Riemannian conjugate gradient methods). This paper proposes a novel general framework that unifies existing Riemannian conjugate gradient methods such as the ones that utilize a vector transport or inverse retraction. The proposed framework also develops other methods that have not been covered in previous studies. Furthermore, conditions for the convergence of a class of algorithms in the proposed framework are clarified. Moreover, the global convergence properties of several specific types of algorithms are extensively analyzed. The analysis provides the theoretical results for some algorithms in a more general setting than the existing studies and new developments for other algorithms. Numerical experiments are performed to confirm the validity of the theoretical results. The experimental results are used to compare the performances of several specific algorithms in the proposed framework. [ABSTRACT FROM AUTHOR]

This paper proposes a variable time step, second-order time filter algorithm for the closed-loop geothermal system for the first time. The algorithm is easy to implement in programming, which greatly improves calculation efficiency, while improving the accuracy of convergence. Theoretically, we mainly analyze the stability and the convergence of the algorithm in the case of variable time steps. Finally, numerical experiments confirm the theoretical results and show that the algorithm can simulate geothermal system well. [ABSTRACT FROM AUTHOR]

Recently, in a series of papers [Y. Rahimi, C. Wang, H. Dong, and Y. Lou, SIAM J. Sci. Comput., 41 (2019), pp. A3649-A3672; C. Wang, M. Tao, J. Nagy, and Y. Lou, SIAM J. Imaging Sci., 14 (2021), pp. 749-777; C. Wang, M. Yan, and Y. Lou, IEEE Trans. Signal Process., 68 (2020), pp. 2660-2669; P. Yin, E. Esser, and J. Xin, Commun. Inf. Syst., 14 (2014), pp. 87-109], the ratio of ℓ1 and ℓ2 norms was proposed as a sparsity inducing function for noiseless compressed sensing. In this paper, we further study properties of such model in the noiseless setting, and propose an algorithm for minimizing ℓ1/ℓ2 subject to noise in the measurements. Specifically, we show that the extended objective function (the sum of the objective and the indicator function of the constraint set) of the model in [Y. Rahimi, C. Wang, H. Dong, and Y. Lou, SIAM J. Sci. Comput., 41 (2019), pp. A3649-A3672] satisfies the Kurdyka--Lojasiewicz (KL) property with exponent 1/2; this allows us to establish linear convergence of the algorithm proposed in [C. Wang, M. Yan, and Y. Lou, IEEE Trans. Signal Process., 68 (2020), pp. 2660-2669] (see equation 11) under mild assumptions. We next extend the ℓ1/ℓ2 model to handle compressed sensing problems with noise. We establish the solution existence for some of these models under the spherical section property [S. A. Vavasis, Derivation of Compressive Sensing Theorems from the Spherical Section Property, University of Waterloo, 2009; Y. Zhang, J. Oper. Res. Soc. China, 1 (2013), pp. 79-105] and extend the algorithm in [C. Wang, M. Yan, and Y. Lou, IEEE Trans. Signal Process., 68 (2020), pp. 2660-2669] (see equation 11) by incorporating moving-balls-approximation techniques [A. Auslender, R. Shefi, and M. Teboulle, SIAM J. Optim., 20 (2010), pp. 3232-3259] for solving these problems. We prove the subsequential convergence of our algorithm under mild conditions and establish global convergence of the whole sequence generated by our algorithm by imposing additional KL and differentiability assumptions on a specially constructed potential function. Finally, we perform numerical experiments on robust compressed sensing and basis pursuit denoising with residual error measured by \ell 2 norm or Lorentzian norm via solving the corresponding ℓ1/ℓ2 models by our algorithm. Our numerical simulations show that our algorithm is able to recover the original sparse vectors with reasonable accuracy. [ABSTRACT FROM AUTHOR]

This paper develops a novel limited-memory method to solve dynamic optimization problems. The memory requirements for such problems often present a major obstacle, particularly for problems with PDE constraints such as optimal flow control, full waveform inversion, and optical tomography. In these problems, PDE constraints uniquely determine the state of a physical system for a given control; the goal is to find the value of the control that minimizes an objective. While the control is often low dimensional, the state is typically more expensive to store. This paper suggests using randomized matrix approximation to compress the state as it is generated and shows how to use the compressed state to reliably solve the original dynamic optimization problem. Concretely, the compressed state is used to compute approximate gradients and to apply the Hessian to vectors. The approximation error in these quantities is controlled by the target rank of the sketch. This approximate first- and second-order information can readily be used in any optimization algorithm. As an example, we develop a sketched trust-region method that adaptively chooses the target rank using a posteriori error information and provably converges to a stationary point of the original problem. Numerical experiments with the sketched trust-region method show promising performance on challenging problems such as the optimal control of an advection-reaction-diffusion equation and the optimal control of fluid flow past a cylinder. [ABSTRACT FROM AUTHOR]