\documentclass[autumn_m,distiller,nototal,slideColor,colorBG]{prosper} \usepackage{graphics} \usepackage{graphicx} %\usepackage{amsmath,amssymb,amsthm} %\usepackage{subeqnarray} \usepackage{hyperref} \usepackage{epsfig} \usepackage{float} \usepackage{seceqn} \usepackage{amssymb} %%%%%%%%%% \newenvironment{pf}{\parindent=0pt{\textbf{Proof: }}}{\hfill $\blacksquare$ \\} \newcounter{ctr} \newtheorem{thm}{Theorem} \newtheorem{lemma}{Lemma} \newtheorem{lem}{Lemma} \newtheorem{cor}{Corollary} \newtheorem{con}{Conjecture} \newtheorem{prop}{Proposition} \newtheorem{defi}{Definition} \newtheorem{example}{Example} \newtheorem{rem}{Remark} \newtheorem{alg}{Algorithm} \newtheorem{ex}{Exercise} \newcommand{\adj}{{\rm adj\,}} \newcommand{\trace}{{\rm trace\,}} \newcommand{\order}{{\rm order\,}} \newcommand{\spanl}{{\rm span\,}} \newcommand{\tr}{{\rm trace\,}} \newcommand{\Rn}{\mathbb{R}^{n}} \newcommand{\St}{{\mathcal S^{n-1}\,}} \newcommand{\Sn}{{\mathcal S^n\,}} \newcommand{\Mn}{{\mathcal M^n\,}} \newcommand{\Mmn}{{\mathcal M^{n-1}\,}} \newcommand{\beq}{\begin{equation}} \newcommand{\eeq}{\end{equation}} \newcommand{\beqr}{\begin{eqnarray}} \newcommand\C{\mathbb C} \newcommand{\II}{\mathcal I} \newcommand{\DD}{\mathcal D} \newcommand{\EE}{\mathcal E} \newcommand{\NN}{\mathcal N} \newcommand{\PP}{\mathcal P} \newcommand{\LL}{\mathcal L} \newcommand{\A}{\mathcal A} \newcommand{\RR}{\mathbb R} \newcommand{\benum}{\begin{enumerate}} \newcommand{\eenum}{\end{enumerate}} \newcommand{\st}{\textnormal{s.t.}} \newcommand{\TRSe}{TRS$_=\;$} \newcommand{\disp}{\displaystyle} \newcommand{\QED}{\hfill ~\rule[-1pt] {8pt}{8pt}\par\medskip ~~} \newcommand{\bpr}{\textsc{proof: }} \newcommand{\epr}{{\red \QED}} \newcommand{\Diag}{{\rm Diag\,}} \newcommand{\diag}{{\rm diag\,}} \newcommand{\offDiag}{{\rm offDiag\,}} \newcommand{\usMat}{{\rm us2Mat\,}} \newcommand{\sMat}{{\rm sMat\,}} \newcommand{\usvec}{{\rm us2vec\,}} \newcommand{\svec}{{\rm svec\,}} \newcommand{\kvec}{{\rm vec\,}} \newcommand{\ZS}{{\mathcal Z_S} } \newcommand{\XS}{{\mathcal X} } \newcommand{\WW}{{\mathcal W} } %\newcommand{\XS}{{\mathcal X_s} } \newcommand{\XSu}{{\mathcal X_{u}} } \newcommand{\XSd}{{\mathcal X_{d}} } \newcommand{\XX}{{\mathcal X} } %\newcommand{\SO}{{\mathcal S^s} } \newcommand{\SO}{{\mathcal S} } \newcommand{\SOu}{{\mathcal S} } %\newcommand{\SOu}{{\mathcal S^{su}} } \newcommand{\Ss}{{\mathcal S} } %\newcommand{\Stn}{{{\mathcal S}^{\scriptsize{\pmatrix{n\cr 2}}}\,}} \newcommand{\Stn}{{\mathcal S^{n-1}\,}} \newcommand{\tn}{{\scriptsize{\pmatrix{n\cr 2}}}\,} \newcommand{\Rtnp}{{\R^{\scriptsize{\pmatrix{n+1\cr 2}}}\,}} \newcommand{\EDMC}{{\bf{\rm EDMC\,}}} \newcommand{\EDM}{{\bf{\rm EDM\,}}} \newcommand{\NEDM}{{\bf{\rm NEDM\,}}} \newcommand{\SDP}{{\bf{\rm SDP\,}}} % Nick's definitions. \def\Rnbyn{\mathbb{R}^{n\times n}} \def\R{\mathbb{R}} \def\normF#1{\|#1\|_F} \def\normt#1{\|#1\|_2} \def\norm#1{\|#1\|} \def\cm{correlation matrix} \def\cms{correlation matrices} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % All refs in roman. \def\eqref#1{{\normalfont(\ref{#1})}} \title{ \blue Robust Solutions of LP and SDP with Applications to the Nearest EDM Problem } \author{Henry Wolkowicz ~~~ (with: Suliman Al-Homidan, Hua Wei) } %%\subtitle{ \small{(in order of appearance)} } \institution{Department of Combinatorics and Optimization\\ University of Waterloo\\ \begin{figure} \psfig{file=UWlogori.ps,height=10mm} \end{figure} } \email{hwolkowicz@uwaterloo.ca } %%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%% \begin{document} \bibliographystyle{plain} \maketitle \begin{slide}{} \begin{center} {\bf Robust Solutions of\\ Large Sparse Semidefinite Programming\\ with Applications to the\\ Nearest Euclidean Distance Matrix Problem } \end{center} \vspace{3mm} Henry Wolkowicz \\ Department of Combinatorics \& Optimization \\ University of Waterloo \mbox{ \begin{figure} \psfig{file=UWlogori.ps,height=20mm} \end{figure} } (with: Suliman Al-Homidan, Hua Wei) \end{slide} \begin{slide}{OUTLINE} \begin{description} \item $\bullet$ Background on SDP; Notation and Motivation \item $\bullet$ Robust, (`non-interior') path-following algorithm for SDP \\ \qquad(outline of GN PCG method using LP) \item $\bullet$ Application to Nearest Euclidean Distance Matrix Problem \item $\bullet$ Numerics (Comparisons with a dual algorithm) \end{description} \end{slide} \begin{slide}{Notation and Motivation} \label{sect:unc} \[ \mbox{(SDP)}\quad \begin{array}{lcl} & \min & f(X)\\ & \mbox{subject to} & \A X=b \\ && X \succeq 0, \end{array} \] ~~\\ where: $f:\Sn \rightarrow \RR$ convex function\\ $~\qquad \Sn$ $n \times n$ real symmetric matrices\\ $~\qquad \A:\Sn \rightarrow \RR^m \mbox{ linear operator, }$\\ $X (\succeq) \succ 0$ denotes positive (semi)definite $\left(\quad (\A X)_i = \left< A_i,X\right> = \trace A_i X,\quad A_i=A_i^T, i = 1\ldots n) \quad \right)$ \end{slide} \begin{slide}{Linear Primal-Dual Pair of SDPs} (looks/behaves like Linear Program, LP) \[ \mbox{(PSDP)}\quad \begin{array}{lcl} & \min & \left=\trace CX\\ & \mbox{subject to} & \A X=b \\ && X \succeq 0, \end{array} \] \[ \mbox{(SDP)}\quad \begin{array}{lcl} & \max & b^Ty\\ & \mbox{subject to} & \A^* y+Z =C \\ && Z \succeq 0, \end{array} \] {\em adjoint operator}: $ \A^* y = \sum^m_{i=1} y_i A_i $ \end{slide} \begin{slide}{(Perturbed) Optimality Conditions} For {\em barrier parameter} $\mu > 0$: \[ \begin{array}{rcl} F_{\mu}(X,y,Z):=\pmatrix{ \A^* y + Z - C \cr \A X -b \cr ZX-\mu I } = 0 \quad \pmatrix{ \mbox{dual feasibility} \cr \mbox{primal feasibility} \cr \mbox{pert. compl. slack.} } \end{array} \] For SDP:\\ $\qquad\qquad F_\mu: \Sn \times \Re^m \times \Sn \rightarrow \Sn \times \Re^m \times \Mn$ \\ i.e. overdetermined nonlinear system \end{slide} \begin{slide}{(Non) Interior Path-Following} \subsection*{Illustration/Motivation on LP Case} \[ (\mbox{LP}) \qquad \begin{array}{ccl} p^* :=&\min & c^Tx \quad (\mbox{or }\left)\\ &\mbox{s.t.} & Ax = b \\ & & x \geq 0 \quad (\mbox{or } x \succeq 0) \end{array} \] \[ (\mbox{DLP}) \qquad \begin{array}{ccl} d^* :=&\max & b^Ty\\ &\mbox{s.t.} & A^T y +z = c \\ & & z \geq 0 \quad (\mbox{or } z \succeq 0) \end{array} \] $A \in \Re^{m \times n}$ full rank (onto); LP and DLP strictly feasible \end{slide} \begin{slide}{Dual log-barrier problem with parameter $\mu >0$} \[ \begin{array}{ccl} d_\mu^* :=&\max & b^Ty+\mu \sum_{j=1}^n \log z_j \quad (+\mu \log \det (z))\\ &\mbox{s.t.} & A^T y +z = c \\ & & z > 0 \qquad (z \succ 0). \end{array} \] stationary point of the Lagrangian / optimality conditions \[ F_\mu(x,y,z)= \left( \begin{array}{cl} A^T y +z - c \\ A x -b \\ X-\mu Z^{-1} \end{array} \right) =0, \begin{array}{cl} \quad x,z>0, \quad (\succ 0)\\ X=\Diag(x),\\Z=\Diag(z) \end{array} \] {\em central path}: set of these solutions $(x_\mu,y_\mu,z_\mu), \mu > 0$ \end{slide} \begin{slide}{Ill-Conditioning} As $\mu \rightarrow 0$, Jacobian $F_\mu^\prime(x,y,z)$ grows ill-conditioned near central path {\bf Cure/Fix:} Make nonlinear equations {\em less nonlinear}, i.e. preconditioning for Newton type methods;\\ premultiply by block-diag matrix with blocks $(I,I,Z)$: \[ \begin{array}{rcl} F_\mu(x,y,z) \leftarrow \pmatrix{I&0&0\cr 0&I&0 \cr 0&0&Z} F_\mu(x,y,z) &=& \left( \begin{array}{cl} A^T y +z - c \\ A x -b \\ ZX-\mu I \end{array} \right) \\ & =:& \left( \begin{array}{cl} R_d \\ r_p \\ R_{ZX} \end{array} \right) \end{array} \] \end{slide} \begin{slide}{Special Structure} Special structure of linearized system can be exploited;\\ linearization for the Newton direction $\Delta s=\pmatrix{\Delta x \cr \Delta y \cr \Delta z}$ is \[ F^{\prime}_\mu(x,y,z) \Delta s= \pmatrix{ 0 & A^T & I \cr A & 0 & 0 \cr Z & 0 & X } \Delta s = -F_\mu(x,y,z). \] \end{slide} \begin{slide}{Overdetermined system in SDP case} \[ \Sn \times \Re^m \times \Sn \rightarrow \Sn \times \Re^m \times \Mn \] apply symmetrization 'undoes preconditioning' \[ \pmatrix{ I & 0 & 0 \cr 0 & I & 0 \cr 0 & 0 & {\cal S} } \pmatrix{ 0 & A^T & I \cr A & 0 & 0 \cr Z & 0 & X } \] e.g. last equation is linearization of:\\ $ZX+XZ - 2 \mu I = 0$ (AHO search direction) \end{slide} \begin{slide}{Reduction/block elimination for the normal equations} \label{sect:rednormal} Step 1 (Eliminate $\Delta z$): \[ \pmatrix{ I & 0 & 0 \cr 0 & I & 0 \cr -X & 0 & I } \pmatrix{ 0 & A^T & I \cr A & 0 & 0 \cr Z & 0 & X } = \left( \begin{array}{ccc} 0 & A^T & I \cr A & 0 & 0 \cr Z & -XA^T & 0 \end{array} \right). \] We let \[ P_Z= \pmatrix{ I & 0 & 0 \cr 0 & I & 0 \cr -X & 0 & I }, \quad K= \left( \begin{array}{ccc} 0 & A^T & I \cr A & 0 & 0 \cr Z & -XA^T & 0 \end{array} \right). \] \end{slide} \begin{slide}{with right-hand side} \[ -\pmatrix{ I & 0 & 0 \cr 0 & I & 0 \cr -X & 0 & I } \left( \begin{array}{cl} R_d \\ r_p \\ R_{ZX} -\mu e \end{array} \right) = \left( \begin{array}{cl} -R_d \\ -r_p \\ XR_d -R_{ZX} \end{array} \right) \] \end{slide} \begin{slide}{Step 2 (Eliminate $\Delta x$)} \[ \begin{array}{rcl} \label{eq:elimdx} \vspace{.1in} F_n:= P_nK & := & \pmatrix{ I & 0 & 0 \cr 0 & I & -AZ^{-1} \cr 0 & 0 & Z^{-1} } \left( \begin{array}{ccc} 0 & A^T & I \cr A & 0 & 0 \cr Z & -XA^T & 0 \end{array} \right) \\ &=& \left( \begin{array}{ccc} 0 & A^T & I_n \cr 0 & AZ^{-1}XA^T & 0 \cr I_n & -Z^{-1}XA^T & 0 \end{array} \right) \end{array} \] $AZ^{-1}XA^T$ can have:\\ $\bullet$ uniformly bounded condition number, e.g. G{\"u}ler et al 1993\\ $\bullet$ structured singularity, e.g. M. Wright 1997, S. Wright 1997\\ But ${\rm cond}(F_n) \rightarrow \infty$. \end{slide} \begin{slide}{The right-hand side becomes:} \[ \begin{array}{l} \vspace{8mm} -P_n P_Z \left( \begin{array}{cl} R_d \\ r_p \\ R_{ZX} \end{array} \right) =\\ \qquad \qquad \left( \begin{array}{cl} -R_d \\ -r_p +A(x-Z^{-1}XR_d -\mu Z^{-1} e) \\ Z^{-1}XR_d -x +\mu Z^{-1} e \end{array} \right) \end{array} \] \end{slide} \begin{slide}{Condition Number} {\bf Proposition } The condition number of $F_n^TF_n$ diverges to infinity if $x(\mu)_i/z(\mu)_i $ diverges to infinity, for some $i$, as $\mu$ converges to 0. The condition number of $(F_\mu^{\prime})^TF_\mu^{\prime}$ is uniformly bounded if there exists a unique primal-dual solution.\\ \bpr Note that \[ F_n^TF_n = \pmatrix{ I_n & -Z^{-1}XA^T & 0 \cr -AXZ^{-1} & (AA^T+(AZ^{-1}XA^T)^2+AZ^{-2}X^2A^T) & A \cr 0 & A^T & I_n }. \] By interlacing of eigenvalues, ... \epr {\bf Corollary } The condition number of $F_n$ is at least $O(1/\mu)$. \epr \end{slide} \begin{slide}{Roundoff Error} {\bf EXAMPLE:} $ A=\pmatrix{1&1}, c=\pmatrix{-1 \cr 1}, b=1. $\\ $ x^*=\pmatrix{1 \cr 0}, y^*=-1, z^*=\pmatrix{0 \cr 2};$\\ initial points: $ \begin{array}{lc} x = \pmatrix{9.183012e-001\cr 1.356397e-008}, z = \pmatrix{2.193642e-008\cr 1.836603e+000},\\ y = -1.163398e+000. \end{array} $ residuals and duality gap: $ \|r_b\| = 0.081699, \quad \|R_d\| = 0.36537, \quad \mu=x^Tz/n = 2.2528e-008 $ 5 decimals rounding before/after arithmetic\\ centering with $\sigma = .1$\\ {\bf BUT:} residuals are NOT order $\mu$. \end{slide} \begin{slide}{Two Search Directions} found using full matrix $F^{\prime}_{\mu}$ and backsolve matrix $F_n$: \[ \pmatrix{\Delta x\cr \Delta y\cr \Delta z}= \pmatrix{ 8.17000e-02 \cr -1.35440e-08 \cr 1.63400e-01 \cr -2.14340e-08 \cr 1.63400e-01 }; \quad = \pmatrix{ -6.06210e+01 \cr -1.35440e-08 \cr 1.63400e-01 \cr 0.00000e+00 \cr 1.63400e-01 } \] error in $\Delta y$ is small; error after backsubstitution for $(\Delta x)_1$ is large.\\ \[ \pmatrix{AZ^{-1}XA^T \cr -Z^{-1}XA^T}= \pmatrix{ 4.18630e+07 \cr -4.18630e+07 \cr -7.38540e-09 } \] %\end{example} \end{slide} \begin{slide}{Stable Reduction} \subsection*{Alternate Second Step} {\bf Assuming!} $A=[I_m ~ E]$.\\ Partition diagonal matrix $Z,X$ using vectors $z=\pmatrix{z_m \cr z_v}, x=\pmatrix{x_m \cr x_v}$, $XA^T=\pmatrix{X_m \cr X_v E^T}$ \[ \begin{array}{rcl} \vspace{.1in} F_s:&=&P_sK=\pmatrix{ I_n & 0 & 0 & 0 \cr 0 & I_m & 0 & 0 \cr 0 & -Z_m & I_m & 0 \cr 0 & 0 & 0 & I_v } \pmatrix{ 0 & 0 & A^T & I_n \cr I_m & E & 0 & 0 \cr Z_m & 0 & -X_m & 0 \cr 0 & Z_v & -X_vE^T & 0 }\\ &=& \pmatrix{ 0 & 0 & A^T & I_n \cr I_m & E & 0 & 0 \cr 0 & -Z_mE & -X_m & 0 \cr 0 & Z_v & -X_vE^T & 0 }. \end{array} \] \end{slide} \begin{slide}{} The right-hand side becomes \[ \begin{array}{rcl} -P_s P_Z \left( \begin{array}{cl} A^T y +z - c \\ A x -b \\ ZXe -\mu e \end{array} \right) = \vspace{.5mm} -P_s \left( \begin{array}{cl} R_d \\ r_p \\ -XR_d +ZXe -\mu e \end{array} \right)\\ \qquad= \left( \begin{array}{cl} -R_d \\ -r_p\\ -Z_mr_p -X_m(R_d)_m +Z_mX_me -\mu e\\ -X_v(R_d)_v +Z_vX_ve -\mu e \end{array} \right) \end{array} \] \end{slide} \begin{slide}{Summary} \subsection*{ Path-following and NOT Interior-point} $\bullet$ staying interior is a heuristic for staying within a neighbourhood of the central path\\ $\bullet$ staying interior is required for numerical accuracy when solving the {\em current} ill-conditioned reduced systems \end{slide} \begin{slide}{Iterations; Backtracking Strategies} \begin{figure}[htb] \epsfxsize=280pt \centerline{\epsfbox{backtrack.eps}} \caption{Iterations for Backtracking Strategies} \label{fig:backtracking} \end{figure} \end{slide} \begin{slide}{(Nearest) Euclidean Distance Matrix Completion using SDP} {\bf Given}:\\ {\em pre-distance matrix} $A\in\Sn$ (nonnegative with zero diagonal)\\ {\em weight matrix} $H\in \Sn$: \[ \mbox{(NEDM)}\quad \mu^*= \min \frac12 \|{H\circ (A-D)}\|_F^2 ~ \mbox{subject to: } D \in \EDM \] \begin{tiny} $\EDM = \{ D=(d_{ij})\in \Sn : d_{ij} = \|x_i-x_j\|^2, \mbox{ for some } x_i \in \Re^k \}$, $k$ is {\em embedding dimension} \end{tiny} \end{slide} \begin{slide}{Applications} e.g. molecular conformation problems in chemistry; multidimensional scaling and multivariate analysis problems in statistics; genetics, geography, .... \end{slide} \begin{slide}{Mixed-Cone Formulation} direct approach using a mixed SDP and second-order (or Lorentz) cone problem: \[ \begin{array}{rcl} &\min &\alpha \\ &\mbox{s.t.}& Y=H\circ(\LL(X)-A),~ \|F{Y}\|_F \le \alpha \\ && X \in \Stn, Y \in \Sn, X \in \SDP \end{array} \] where $X \in SDP \Rightarrow \LL(X) \in \EDM$\\ (Public domain software packages are available) \end{slide} \begin{slide}{Connection between \SDP and \EDM} \[ B=[x_1~x_2~\ldots ~ x_n], \quad k \times n \] \[ D_{ij}=\|x_i-x_j\|^2 = -2x_i^Tx_j+ \|x_i\|^2 + \|x_j\|^2 \] \[ D= -2B^TB+ e \left(\diag (B^TB)\right)^T+ \left(\diag (B^TB)\right) e^T \] ~~\\ With $X=B^TB \succeq 0$\\ (relaxations with $D,X$ lose embedding dimension $k$ information) \end{slide} \begin{slide}{Operator Notation} $\usvec, \usMat, \svec, \sMat $ \[x=\svec(X) \in \R^{n+1 \choose 2}, \quad X= \sMat(x) \] $\sqrt{2}$ times vector (columnwise) from upper-triang of $X$. ${n+1\choose 2}=n(n+1)/2$; $\sqrt{2}$ guarantees isometry. $\sMat:=\svec^{-1}$ mapping into $\Sn$\\ adjoint transformation $\sMat^*= \svec$: \begin{eqnarray*} \left<\sMat(v),S\right> &=& \trace \sMat(v)S \\ &=& v^T \svec(S) = \left \end{eqnarray*} \end{slide} \begin{slide}{PSD Equivalence} \begin{center} $D$ is \EDM \quad ($\subset \Sn$)\\ ~\\ \fbox{ {\em iff} } ~\\ \end{center} \[ D=\LL( X):=\pmatrix{0 & \diag ( X)^T \cr \diag ( X) & \diag(X)e^T+e\diag(X)^T -2X }, \] \begin{center} for some $X \succeq 0, X \in {\mathcal S}^{n-1}$ \end{center} ($e$ is vector of ones) \[ \LL: \Stn \rightarrow \Sn, \quad \LL({\mathcal S}^{n-1}_+) = \EDM \] \end{slide} \begin{slide}{} with partition: \[ D=\left[\matrix{\alpha & d^T\cr d& {\bar D} \cr}\right], \] where $\alpha \in \R$ \[ \LL^*(D)=2\left(\Diag(d)+\Diag(\bar De) -\bar D\right) \] \[ \LL^\dagger (D)= \frac 12 \left( de^T+ed^T - \bar D \right) \] \[ \LL^*,\LL^\dagger: \Sn \rightarrow \Stn, \quad \LL^\dagger (EDM)={\mathcal S}^{n-1}_+ \] \end{slide} \begin{slide}{Duality and Optimality Conditions} (using $X = \sMat(x)+I$) an equivalent problem is: \[ \mu^* := \min \frac 12 \|H \circ (A-\LL(X)) \|_F^2 \quad\mbox{subject to \quad $X \succeq 0$} \] strong (Lagrangian) duality holds (Slater's holds for primal and holds for dual if the graph is complete) \[ \mu^*=\nu^*:= \max_{\Lambda \succeq 0} \min_{X} \frac 12 \|H \circ (A-\LL(X)) \|_F^2 - \trace \Lambda X \] \end{slide} \begin{slide}{{\bf Wolfe dual} and optimality conditions} With \[ C:=\LL^*(H^{(2)} \circ A), \] optimality conditions are: \[ \begin{array}{rcll} X&:=&\sMat(x) \succeq 0 & \mbox{(primal feasibility)}\\ \Lambda &:=& \LL^* \left\{ H^{(2)} \circ( \LL(X) ) \right\}-C , ~~\Lambda \succeq 0 & \mbox{(dual feasibility)}\\ \Lambda X &:=& 0 & \mbox{(compl. slack.)}\\ \end{array} \] equivalent dual problem: \begin{equation} \label{duals} \begin{array}{rcl} &\max &\frac 12 \|H \circ (A-\LL(X)) \|_F^2 - \trace \Lambda X \\ & \mbox{subject to}& \Lambda =\LL^* \left\{ H^{(2)} \circ( \LL(X) ) \right\}-C \\ && \Lambda \succeq 0. \end{array} \end{equation} \end{slide} \begin{slide}{Bilinear System} eliminate $\Lambda$\\ exact primal-dual feasibility during iterations\\ full rank Jacobian at optimality. ~~\\ {\em single bilinear (perturbed) equation} in $x$; $\quad F_{\mu}(x): \R^{n \choose 2} \rightarrow \Mmn$ \fbox{ \begin{minipage}{4.05in} \[ \begin{array}{c} \vspace{.5mm} F_{\mu}(x):= \left[ \LL^* \left\{ H^{(2)} \circ( \LL(X) ) \right\}-C \right]X - \mu I = 0 \end{array} \] \end{minipage} } ~\\ typical SDP - overdetermined system of bilinear equations\\ current approach is to symmetrize - which results in ill-conditioning! from rank deficient Jacobian at optimality.\\ BUT, here, no symmetrization used;\\ solve using (an inexact) Gauss-Newton method - with PCG \end{slide} \begin{slide}{Linearization} Let $\WW(x):= \LL^* \left\{ H^{(2)} \circ( \LL(x) ) \right\}$ ~\\ Linearization for search direction $\Delta x$ at current $x=\svec (X)$: \fbox{ \begin{minipage}{4.05in} \[ \begin{array}{lcl} F^{\prime}_{\mu}(x) \Delta x = \left[ \WW(x) -C \right] \Delta x + \left[ \WW( \Delta x) \right]X \end{array} \] \end{minipage} } ~\\ This is a linear, full rank, overdetermined system.\\ Our search direction $\Delta x$ is its (approx.) least squares solution. \end{slide} \begin{slide}{Algorithm: p-d i-e-p framework} {\bf $\bullet$ Initialization:} \begin{quote} {\bf $\bullet\bullet$ Input data:} a pre-distance $n\times n$ matrix $A$ \\ {\bf $\bullet\bullet$ Positive tolerances:} \[\epsilon_1 \mbox{ (stopping), } \epsilon_2 \mbox{ (lss accuracy), } \epsilon_3 \mbox{ (crossover), } \] \\ {\bf $\bullet\bullet$ Find initial strictly feasible points:} both $X^0,\Lambda^0:=\left(\WW (X) -C \right)~\succ~0$; $\mu>0$\\ {\bf $\bullet\bullet$ Set initial parameters:} \[ {\rm gap}=\trace \Lambda^0X^0;~~ \mu={\rm gap}/n;~~{\rm objval}=f(X^0);~~k=0. \] \end{quote} {\bf $\bullet$ while} $\min \{ \frac {\rm gap}{{\rm objval}+1} , {\rm objval} \} > \epsilon_1$ \begin{quote} {\bf $\bullet$$\bullet$ solve lss for search direction} (accuracy $\epsilon_2\min\{\mu,1\}$) \[ ~~~~F^{\prime}_{\sigma \mu}(x^k) \left( \Delta x^k \right) = -F_{\sigma \mu}(x^k), \] where $\sigma_k$ centering, $\mu_k=\frac{1}{n}\trace (\WW(X^{k})-C)X^k$ \[ X^{k+1} = X^k + \alpha_k \Delta X^k, ~~ \alpha_k > 0, \] so that both $X^{k+1},(\WW(X^{k+1})-C) \succeq 0$\\ ($\alpha_k=1$ after the crossover.)\\ {\bf $\bullet$$\bullet$ update} \[ k \leftarrow k+1 \quad \mbox{ and then} \] \[ \sigma_k \quad \left(\mbox{set } \sigma_k=0 \mbox{ if } \min \{ \frac {\rm gap}{{\rm objval}+1} , {\rm objval} \} < \epsilon_3\right) \] \end{quote} {\bf $\bullet$ end (while)}.\\ {\bf $\bullet$ Conclusion: } $D=\LL(X) \in \EDM$ is approx. to $A$ \end{slide} \begin{slide}{Crossover} After the {\bf crossover}, centering $\sigma=0$ and steplength $\alpha=1$, we get q-quadratic convergence; allows for {\em warm starts}. {\bf Long steps} can be taken {\em beyond} the positivity boundary. (tests show improved convergence rates) \end{slide} \begin{slide}{Preconditioning} \[ \left(\Lambda+\XS \WW \right) P^{-1}(\widehat{\Delta x}) = - F_{\mu}(x), \] where \[ \widehat{\Delta x}=P(\Delta x) \] \end{slide} \begin{slide}{Diagonal Preconditioning} Optimal scaling Dennis and W. (1993) full rank matrix $A\in\R^{m \times n}$, $m \ge n$, with condition number $\omega (K) := n^{-1}\trace (K)/\det(K)^{1/n}$, the optimal scaling \[ \min \omega ( (AD)^T(AD)) \quad \mbox{subject to: $D$ positive and diagonal} \] solution: $d_{ii}= 1/\normt{A_{:i}}, i=1,\ldots,n$ {\em explicit} expressions for preconditioner\\ inexpensive \end{slide} \begin{slide}{} {\em diagonal} operator $P$; evaluate using columns of $F_\mu^\prime(v)$. \\ $k\cong (i,j),~1\leq i < j\leq n$, strictly upper triangular part\\ \[ \begin{array}{rcl} \| (\Lambda+\XS\WW)(e_k) \|^2_F &=& \| \Lambda(e_k) \|^2_F + \| (\WW(e_k))X \|^2_F\\ && \qquad +\left<\Lambda(E_{ij}),(\WW(E_{ij}))X\right>, \end{array} \] where \[ \begin{array}{rcl} \Lambda (e_k) &=& \left\{ \begin{array}{cc} \frac 1{\sqrt{2}} \left( \Lambda_{:i}e_j^T+\Lambda_{:j}e_i^T \right), & \mbox{if } i < j\\ \left( \Lambda_{:i}e_i^T \right), & \mbox{if } i = j. \end{array} \right. \end{array} \] and $\XS \WW $ ..... inexpensive - 50\% reduction in LSQR iterations \end{slide} \begin{slide}{Numerical Tests} Pentium 4; MATLAB 6.5; 1 GIG RAM. crossover heuristic: relative duality gap $< .1$. Stopping criteria (relative duality gap) $< \epsilon_1 = 1e-10$.\\ (But - average accuracy attained $1e-13$, q-quadratic convergence.) \end{slide} \begin{slide}{} \begin{figure}[htb] \epsfxsize=220pt \centerline{\epsfbox{dens0050103.ps}} \caption{density .0005:.001:.003, CPU times and nnz($\Lambda$), n=200} \label{fig:densityn200} \end{figure} \end{slide} \begin{slide}{Conclusion} Gauss-Newton direction: \begin{description} \item Advantages/Disadvantages: \begin{description} \item Robust, warm starts are simple, longer steps \item exact primal and dual feasibility at each iteration \item Can apply CG-type approaches \item q-quadratic convergence \item scale-invariant on the right \end{description} \item Future: \begin{description} \item Need large sparse QR efficient as Cholesky \item predictor-corrector \end{description} \end{description} \end{slide}{} \begin{slide}{\EDM\ Completion Problem, \EDMC} \label{sect:completion} $\bullet$ given certain fixed elements of a \EDM matrix $A$\\ $\bullet$the other elements are unknown (free)\\ $\bullet$complete this matrix to an \EDM \[ \Ss =\{ (i,j): A_{i,j}=\frac 1{\sqrt 2}b_k ~\mbox{is known, fixed},~i < j\}, ~|S|=m, \] \[ (\EDMC)\qquad \begin{array}{ccc} \mu^* := & \min & f(X):=\frac 12\|X\|^2_F \\ & \mbox{subject to} & \A(X) =b \\ & & X \succeq 0, \end{array} \] constraint $\A=\II\cdot \LL:\St \rightarrow \R^{|S|}$ yields interpolation conditions \[ \A(X)_{ij}=\trace E_{ij}\LL(X)=b_{k}, \quad \forall k\cong (ij)\in S, \] \end{slide}{} \begin{slide}{Duality/Optimality for EDMC} $\bullet$strict convexity, coercivity {\blue implies} compact level sets\\ $\bullet$\EDMC attained and no duality gap (actually primal and dual attainment) Lagrangian dual \[ \mu^*=\nu^*:= \max_{\Lambda \succeq 0,y\in \R^{|S|}} \min_{X} \frac 12 \|X \|_F^2 +y^T (b-\A(X)) - \trace \Lambda X \] \end{slide}{} \begin{slide}{characterization of optimality} \begin{thm} \label{eq:optcharctexpl} Suppose that the feasible set of \EDMC is not the empty set. Then the optimal solution of \EDMC is $D = \LL \left([\A^*(y)]_+\right)$, where $y$ is the unique solution of the single equation \[ \A \left([\A^*(y)]_+\right) = b, \] and $B_+$ denotes the projection of the symmetric matrix $B\in \St$ onto the cone $\PP_{n-1}$. \end{thm} \end{slide}{} \begin{slide}{Proof} optimality conditions after differentiation {\blue \[ \begin{array}{cc} X=\A^*(y)+\Lambda \succeq 0, \quad \Lambda \succeq 0, & \mbox{dual feasibility} \\ \A(X)=b & \mbox{primal feasibility} \\ \quad \Lambda X=0 & \mbox{complementary slackness} \end{array} \] } This means that $\A^*(y)=X-\Lambda$, where both $X\succeq 0, \Lambda \succeq 0$, and $\Lambda X = 0$. Therefore the three symmetric matrices $W=\A^*(y),X,\Lambda$ are mutually diagonalizable. We write $X=PD_XP^T,~\Lambda=PD_\Lambda P^T$, i.e. we conclude that $W=\A^*(y)=P\left(D_X-D_\Lambda\right)P^T$, $D_XD_\Lambda = 0$. Therefore $[\A^*(y)]_+ = PD_XP^T = X$. \epr \end{slide}{} \begin{slide}{Efficient/Explicit Solution if $y \geq 0$} large class ({\red generic?}) can be solved in polytime. \begin{cor} \label{cor:explicitopt} {\blue The linear operator $\A$ is onto and $\A \A^*$ is nonsingular. Suppose that $y=(A\A^*)^{-1}b \in \R^m_+$. Then \[ D=\LL \left(\A^*(y)\right) \] is the unique solution of \EDMC\@. } \end{cor} \bpr That $\A$ is onto follows from the definitions. If $y \geq 0$, then the matrix $\II(y) \geq 0$ with 0 diagonal. Therefore, $X=\LL^*(\II(y))$ is diagonally dominant with nonnegative diagonal, i.e. $X \succeq 0$ by Gersgorin's disk theorem. This implies that $D$ is a distance matrix and it satisfies the interpolation conditions, i.e. it satisfies the optimality conditions in Theorem \ref{eq:optcharctexpl}. \epr \end{slide} \begin{slide}{Numerics} though not $y = \A^\dagger b\geq 0$ in general, still get\\ {\blue a distance matrix $D$, i.e. $\A^* (y) \succeq 0$}. \\ $n=10:10:100$; density $.1:.1:.8$.\\ dimension versus density with number of failures in $100$ tests: \[ \pmatrix{ \qquad ~~~ n\backslash\mbox{density} & .1 &.2&.3&.4&.5&.6&.7&.8\cr 10& 19 & 27 & 29 & 25 & 32 & 27 & 20 & 38\cr 20& 6 & 20 & 23 & 22 & 27 & 21 & 28 & 28\cr 30& 8 & 8 & 9 & 9 & 11 & 16 & 17 & 24\cr 40& 2 & 2 & 6 & 5 & 14 & 17 & 20 & 17 \cr 50& 2 & 0 & 2 & 8 & 7 & 8 & 15 & 12 \cr 60& 1 & 1 & 1 & 1 & 3 & 8 & 15 & 11 \cr 70& 2 & 0 & 3 & 1 & 5 & 7 & 6 & 15 \cr 80& 1 & 0 & 0 & 4 & 2 & 4 & 9 & 9 \cr 90& 1 & 0 & 0 & 1 & 3 & 2 & 5 & 6 \cr 100& 0 & 0 & 0 & 0 & 1 & 6 & 5 & 5 }. \] \end{slide} %\end{multicols} \end{document}