\begin{filecontents}{seceqn.sty}
% seceqn.sty
% From TeXhax 89.37
%
% Substyle file for use with "article" to cause equations to be numbered
% within sections.
%
% 13 Apr 89	Jerry Leichter
\typeout{Document Option `seceqn':  13 Apr 89}
\@addtoreset{equation}{section}		%Make equation=0 when section steps
\def\theequation{\thesection.\arabic{equation}}
					%How an equation number looks

% From NUMINSEC.STY (with SIAM stuff).
\@addtoreset{figure}{section}
\def\thefigure{\thesection.\@arabic\c@figure}

\@addtoreset{table}{section}
\def\thetable{\thesection.\@arabic\c@table}
\end{filecontents}

%\documentclass[landscape,slideColor,colorBG]{seminar}

\documentclass[autumn_m,distiller,nototal,slideColor,colorBG]{prosper}

%\documentclass[a4,portrait]{seminar}
%\usepackage{semcolor}
%\input{seminar.bug}
%\usepackage{theorem,exscale}
%\usepackage{fullpage}        % was not available on PC?
%\usepackage{float}
%\usepackage{psfrag}
%\usepackage{algorithmic}     % was not available on PC?
%\usepackage{hyperref}
%\usepackage{seceqn}
%\usepackage{epsfig}
%\usepackage{amssymb}
%\usepackage{amssymb,epsfig}
%\usepackage{amsmath,amsthm,amsfonts,latexsym,amssymb,epsfig}
%\usepackage{fullpage}        % was not available on PC?
%\usepackage{float}
%\usepackage{algorithmic}        % was not available on PC?
%\usepackage{epsf}
%\usepackage{hyperref}
%\usepackage{seceqn}
%\usepackage[dvips]{graphicx}

\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}[section]
\newtheorem{lemma}{Lemma}[section]
\newtheorem{cor}{Corollary}[section]
\newtheorem{con}{Conjecture}[section]
\newtheorem{prop}{Proposition}[section]
\newtheorem{defi}{Definition}[section]
\newtheorem{example}{Example}[section]
\newtheorem{rem}{Remark}[section]
\newtheorem{alg}{Algorithm}[section]
\newtheorem{ex}{Exercise}[section]
\newcommand{\adj}{{\rm adj\,}}
\newcommand{\trace}{{\rm trace\,}}
\newcommand{\spanl}{{\rm span\,}}
\newcommand{\tr}{{\rm trace\,}}
\newcommand{\Rn}{\mathbb{R}^{n}}
\newcommand{\Sn}{{\mathcal S^n\,}}
\newcommand{\St}{{\mathcal S^{n-1}\,}}
\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{\NN}{\mathcal N}
\newcommand{\LL}{\mathcal L}
\newcommand{\A}{\mathcal A}
\newcommand{\PP}{\mathcal P}
\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}{\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{\II}{{\mathcal I} }
%\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{\EDM}{{\bf{\rm EDM\,}}}
\newcommand{\EDMC}{{\bf{\rm EDMC\,}}}
\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{
%\cyan
Robust Algorithms for
Large Sparse Semidefinite Programming (SDP)
}
\subtitle{ \small{
with Applications to the Nearest
Euclidean Distance Matrix Problem
} }
\author{Henry Wolkowicz}
\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}{
 Toulouse 2004}
\begin{center}
{\bf First Joint Canada-France meeting of the Mathematical Sciences }
\end{center}

\begin{center}
\mbox{
\begin{figure}
\psfig{file=leftcornertop3.ps,width=32mm,height=23mm}
\psfig{file=headerleft.ps,width=50mm,height=15mm}
\end{figure}
}
\end{center}
\begin{center}
{\bf Robust Algorithms\\ for\\ Large Sparse Semidefinite Programming (SDP)}
\end{center}
\vspace{.5in}
Friday July 15, 2004.




\end{slide}

\begin{slide}{OUTLINE}
\begin{description}
\item $\bullet$
Background on SDP; Notation and Motivation
\item 
\item $\bullet$
Robust, (`non-interior') path-following algorithm for SDP \\
\qquad(outline  of GN PCG method using LP)
\item 
\item $\bullet$
Application to Nearest Euclidean Distance Matrix Problem
\item 
\item $\bullet$
Numerics (Comparisons with a dual algorithm)
\end{description}




\end{slide}
\begin{slide}{Notation and Motivation}
\[
\mbox{(SDP)}\quad
\begin{array}{lcl}
& \min  & f(X)\\
&  \mbox{subject to} &  \A X=b \\
 &&   X \succeq 0,
\end{array}
\]
where:
\[
\begin{array}{rlc}
f:\Sn \rightarrow \RR &  \mbox{convex function}\\
 \Sn                  & n \times n \mbox{ real symmetric matrices}\\
X (\succeq) \succ 0  &  \mbox{denotes positive (semi)definite} \\
  \A:\Sn \rightarrow \RR^m &  \mbox{ linear transformation }\\
\end{array}
\]
$\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<C,X\right>=\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}{(some of the) APPLICATIONS}

\begin{description}
\item[$\bullet$] 
Relaxations of
\underline{\cyan \bf hard combinatorial problems}: e.g.
max-cut; graph partitioning; quadratic assignment problem;
max-clique.

\item[$\bullet$] 
\underline{\cyan \bf NLP} e.g.:
quasi-Newton updates that preserve positive definiteness;
Trust region algorithms for large scale minimization;
Extended SQP techniques for constrained minimization.

\item[$\bullet$] 
Partial Hermitian \underline{\cyan \bf matrix completion problems} and
Euclidean distance matrix completion problems.

\item[$\bullet$] 
\underline{\cyan \bf Engineering problems} such as:
Ricatti equations; 
min-max eigenvalue problems; matrix norm minimization; eigenvalue
localization.
\end{description}

\end{slide}
\begin{slide}{ SIMILARITIES to LP: (i) Duality}
payoff function, player $Y$ to player $X$ (Lagrangian)
\[ L(X,y) :=  \tr (CX) +y^t(b-{\cal A}X)
\]

Optimal (worst case) strategy for player $X$:
\[p^* = 
      \max_{X \succeq 0 }  \min_{y }  L(X,y) 
\]
Using the {\red \em hidden constraint} \quad  $b-{\cal A}X=0$,\\
 recovers primal problem.

\end{slide}
\begin{slide}{apply adjoint}
\[ 
\begin{array}{rcl}
L(X,y) &=&  \tr (CX) +y^t(b-{\cal A}X)\\
       &=&  b^ty + \tr \left(C -{\cal A}^*y \right) X
\end{array}
\]

~\\
~\\
adjoint operator,~~ ${\cal A}^*y= \sum_i y_i A_i $
\[
\left<{\cal A}^*y,X \right>= \left<y,{\cal A}X \right>, ~~~ \forall X,y
\]

{\red \em Hidden Constraint:} $C-{\cal A}^*y \preceq 0$

\end{slide}
\begin{slide}{exploit {\em Hidden Constraint}}
\vspace{-.1in}
\[p^* = 
       \max_{X \succeq 0 } \min_{y }  L(X,y) 
\leq d^*:=\min_y \max_{X \succeq 0} L(X,y) 
\]
dual obtained from optimal strategy of competing player, Y. \\
{\red \em Hidden Constraint:} $C-{\cal A}^*y \preceq 0$
 yields the dual
\[ {\bf (DSDP)}
\begin{array}{ccc}
    d^*=& \min &b^ty \\
 &  \mbox{s.t.} & {\cal A}^*y \succeq  C \\
    \end{array}
\]
\vspace{-.1in}
for the primal
\vspace{-.1in}
\[ {\bf (PSDP)}
\begin{array}{ccc}
    p^*=  & \max &\tr CX \\
 &  \mbox{s.t.} & {\cal A}X = b\\
  && X \succeq 0
    \end{array}
\]
\end{slide}
\begin{slide}{Characterization of Optimality}
   primal-dual pair $X,y$~~(slack $Z\succeq 0$)
 \[ 
\begin{array}{cc} 
    {\cal A}^*y -Z = C  & \mbox{dual feasibility}\\
~\\
    AX = b & \mbox{primal feasibility}\\
~\\
    ZX  = 0 & \mbox{complementary slackness}
\end{array}
\]
\[
    ZX = \mu I ~~~~~ {\cyan \mbox{perturbed C.S., } \mu > 0}
\]

Basis for methods:
\begin{description}
\item[$\bullet$] 
primal simplex \quad (maintain: primal feas. \& compl. slack.)
\item[$\bullet$] 
dual simplex \quad (maintain: dual feas. \& compl. slack.)
\item[$\bullet$] 
interior point \quad (maintain: primal feas. \& dual feas.)
\end{description}


\end{slide}
\begin{slide}{SDP Application: (Direct) Max-Cut Relaxation}
Graph $G=(E,V)$; $|V|=n$ (nodes); $w_{ij}$ weights on edges;
\[
 \begin{array}{c}
    \max ~ \frac 12 \sum_{i<j} w_{ij}(1-x_ix_j),~~~x \in  \{ \pm 1 \}^n.
\end{array}
\]
Equate $x_i=1$ with $i$ in set  $\cal I$ and $x_i=-1$ otherwise.\\
Equivalent problem: homogeneous {\em $(\pm 1)$-QQP}
\[
\mu^*:=\max ~ q(x):=x^tQx=\trace Qxx^T, \quad x \in \{ \pm 1 \}^n.
\]
{\red REPLACE $x \in \{ \pm 1 \}^n$ WITH CONSTRAINTS $x_i^2=1$ ??!! }
\vspace{-.1in}
\[
\begin{array}{lc}
\mbox{\cyan LIFTING}:\qquad\qquad &
X=xx^t 
\end{array}
\]
\vspace{-.1in}
Relax the rank-1 condition on $X$ to get linear SDP.
\[ 
\mu^* \leq \max \{ \tr QX:   \diag (X) = e, X\succeq 0\}  
\] 
\end{slide}
\begin{slide}{SDP from general quadratic approx? {\cyan (Lagr. Relax.!)}}
\[q_i(y):=\frac 12 y^tQ_iy+y^tb_i + c_i,~y\in \Re^n \]
\[ {\bf (QQP)}
\begin{array}{ccc}
    q^*=  & \min &q_0(y) \\
 &  \mbox{s.t.} & q_i(y)\leq 0\\
     &&  i=1,\ldots m
    \end{array}
\]
\[ 
\begin{array}{ccc}
\mbox{Lagrangian}: \qquad
  L(y,x) &=& q_0(y) + \sum_{i=1}^m x_iq_i(y)
\end{array}
\]
or equivalently
\vspace{-.2in}
\[ 
\begin{array}{ccc}
  L(y,x) &=& \frac 12 y^t (Q_0 +\sum_{i=1}^m x_iQ_i)y 
             \quad \mbox{ (quadratic in } y) \\
 &&   +y^t(b_0 +\sum_{i=1}^m x_ib_i)
             \quad \mbox{ (linear in } y)\\
   &&+ (c_0 +\sum_{i=1}^m x_ic_i) 
             \quad \mbox{ (constant in } y)
\end{array}
\]
\end{slide}
\begin{slide}{Weak Duality}
follows from definition of dual program and hidden constraints:
\[d^*=\max_{x\geq 0} \min_y L(y,x) \leq q^* = \min_y \max_{x\geq 0} L(y,x).
\]

Now {\cyan homogenize}; multiply linear term by new variable $y_0$
\[ y_0y^t(b_0 +\sum_{i=1}^m x_ib_i), ~~ y_0^2=1.  \]
and add new constraint to Lagrangian (Lagrange multiplier $t$)
\[
t(y_0^2-1)
\]

\end{slide}
\begin{slide}{Homogenization}
\[
\begin{array}{lclc}
  d^*
&=&
  \max\limits_{x \geq 0} \min\limits_y &L(y,x)
\\&=& 
\max\limits_{x \geq 0} 
         \min\limits_{y_0^2=1}& \frac 12 y^t (Q_0 +\sum_{i=1}^m x_iQ_i)y 
                            ~~~~~+ty_0^2 \\
 &&&   +y_0y^t(b_0 +\sum_{i=1}^m x_ib_i)\\
 &&&+ (c_0 +\sum_{i=1}^m x_ic_i)
                            ~~~~~-t\\
  &=& \max\limits_{x \geq 0,t} 
         \min\limits_{y}& \frac 12 y^t (Q_0 +\sum_{i=1}^m x_iQ_i)y 
                            ~~~~~+ty_0^2 \\
 &&&   +y_0y^t(b_0 +\sum_{i=1}^m x_ib_i)\\
 &&&+ (c_0 +\sum_{i=1}^m x_ic_i)
                            ~~~~~-t\\
\end{array}
\]
{\cyan hidden semidefinite constraint yields SDP}
\end{slide}
\begin{slide}{Apply Hidden SDP Constraint\\ ({\cyan Hessian psd})}
\[ 
B:=\left( \begin{array}{cc}
      0 & b_0^t \\ b_0 &Q_0
   \end{array}  \right)
 \mbox{ and } A: \Re^{m+1} \rightarrow {\cal S}_{n+1}
\] 
\[ 
A \left( \begin{array}{c}
      t \\ x
   \end{array}  \right)
   := -\left[ \begin{array}{cc}
        t &  \sum_{i=1}^m x_ib_i^t \\
      \sum_{i=1}^m x_ib_i  & \sum_{i=1}^m x_i Q_i
        \end{array}   \right]
\]
\[
\mbox{ Lagrangian psd }:\quad
B-A \left( \begin{array}{c}
      t \\ x
   \end{array}  \right)
 \succeq 0.
\]

{\bf NOTE There is NO hidden constraint on the $Q_i$ 
if all $q_i$ are convex.}
Better algorithms exist for the convex case, 
e.g. proximal methods, using quadratic cones, ...

\end{slide}
\begin{slide}{Dual of Dual $\rightarrow$ SDP Relaxation}
dual program is equivalent to SDP (with $c_0=0$)
\[ {\bf (D)}
\begin{array}{ccc}
    d^*=  & \mbox{sup} & -t +\sum_{i=1}^m x_ic_i \\
 &  \mbox{s.t.} & A\left( \begin{array}{c}
      t \\ x
   \end{array}  \right)
 \preceq B\\
  && x \in \Re^m, t \in \Re
    \end{array}
\]
As in LP, dual of dual is obtained from optimal
strategy of the competing player:
\vspace{-.2in}
\[ {\bf (DD)}
\begin{array}{ccc}
    d^*=& \inf &\tr BU \\
 &  \mbox{s.t.} & A^*U = \left( \begin{array}{c}
      -1 \\ c
   \end{array}  \right) \\
  && U \succeq 0.
    \end{array}
\]
\end{slide}
\begin{slide}{Tractable Relaxations}
In some sense, Lagrangian relaxation is {\bf best tractable
relaxation}. 

There are {\em higher order} relaxations:\\
 e.g. from $X=xx^T$ from
max-cut relaxation (from $x_j^2=1$)
\[
\mbox{ 2nd LIFTING}: \quad
x_i x_j^2x_k =x_i x_k, \quad   Y=\pmatrix{1\cr \svec X}\pmatrix{1& \svec X}
\]

Public domain software: e.g. {\cyan NEOS}\\
URL: www-neos.mcs.anl.gov


\end{slide}


\begin{slide}{({\cyan Perturbed}) Optimality Conditions}
For {\em barrier parameter} ${\cyan \mu} > 0$:
\[
\begin{array}{rcl}
F_{\mu}(X,y,Z):=\pmatrix{
    \A^* y + Z - C \cr
    \A X -b \cr
    ZX-{\cyan \mu} I 
} = 0 \qquad
\pmatrix{
    \mbox{dual feasibility} \cr
    \mbox{primal feasibility} \cr
    \mbox{{\cyan 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<c,x\right>)\\
        &\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}
\]
Assume:
{\yellow
$A \in \Re^{m \times n}$ full rank (onto); LP, DLP strictly feasible
}


\end{slide}
\begin{slide}{dual log-barrier problem;\\ parameter $\mu >0$}
\[
        \begin{array}{cclc}
        d_\mu^* :=&\max   & b^Ty+\mu \sum_{j=1}^n \log z_j &
                                    ({\yellow+\mu \log \det (z)})\\
        &\mbox{s.t.} & A^T y +z = c & ({\yellow A^T \cong A^*})\\
        &     &  z >  0 & ({\yellow 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 ({\yellow \succ 0})\\ {\green X=\Diag(x)}\\
                                  {\green Z=\Diag(z)}
\end{array}
 \]

{\cyan \em central path}:=
set of solutions $(x_\mu,y_\mu,z_\mu), \mu > 0$



\end{slide}
\begin{slide}{Jacobian Ill-conditioning}


As $\mu \rightarrow 0$,  Jacobian $F_\mu^\prime(x,y,z)$
ill-conditioned near central path

{\cyan \bf Cure/Fix:} Make nonlinear equations {\red \em less nonlinear},
i.e. {\red 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}
 \]
{\yellow recovers modern primal-dual optimality paradigm}



\end{slide}
\begin{slide}{Exploited Special Structure}
linearization for the Newton
direction 
\[
\Delta s=\pmatrix{\Delta x \cr \Delta y \cr \Delta z}
\]
\[ 
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}{{\cyan Overdetermined} system in SDP case}
\[
\Sn \times \Re^m \times \Sn \rightarrow  \Sn \times \Re^m \times \Mn
\]
 apply symmetrization; {\red 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 after symmetrization:\\
$ZX+XZ - 2 \mu I = 0$ (AHO search direction)



\end{slide}
\begin{slide}{Reduction/{Block-Elimination}\\ $\rightarrow$ Normal Equations}
\label{sect:rednormal}
Step 1 ({\green Eliminate $\Delta z$ from row 3}):
\[
 \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).
\]
Define:
 \[ 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$ from row 2}
({\green and scale row 3})
\[
\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}
\]
${\cyan AZ^{-1}XA^T}$ can have:\\
 $\bullet$ {\yellow uniformly bounded condition number}, 
         e.g. G{\"u}ler et al 1993\\
 $\bullet$ {\green structured singularity}, e.g. S. Wright 95,97/ M. Wright 1999\\
But {\red ${\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)
=
\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}{Ill-conditioning}
{\green \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, ...
{\black \epr}
{\green \bf Corollary }
The condition number of $F_n$ is at least $O(1/\mu)$.
{\black \epr}



\end{slide}
\begin{slide}{EXAMPLE} 
$ A=\pmatrix{1&1}, c=\pmatrix{-1 \cr 1}, b=1,  $\\
$ x^*=\pmatrix{1 \cr 0}, y^*=-1, z^*=\pmatrix{0 \cr 2};$\\
{\yellow 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}
$

{\cyan residuals and duality gap}:
$
\|r_b\| = 0.081699, 
\|R_d\| = 0.36537,  \mu=x^Tz/n = 2.2528e-008
$

5 decimals rounding before/after arithmetic\\
centering with $\sigma = .1$\\
{\bf BUT:} {\red residuals are NOT order $\mu$}.

\end{slide}
\begin{slide}{search directions found}
\[
\begin{array}{cccc}
\mbox{using}: &  \qquad \mbox{full matrix } F^{\prime}_{\mu} &\mbox{and}& 
                  \mbox{backsolve matrix } F_n \\
&\pmatrix{\Delta x\cr \Delta y\cr \Delta z}= \pmatrix{ 
   {\yellow 8.17000e-02}  \cr
   -1.35440e-08   \cr
    1.63400e-01   \cr
   -2.14340e-08   \cr
    1.63400e-01   
           };  &&
   \pmatrix{
       {\green -6.06210e+01} \cr
       -1.35440e-08 \cr
         1.63400e-01 \cr
        0.00000e+00 \cr
         1.63400e-01 
             }
\end{array}
\]
error in $\Delta y$ is small;\\
 but error after backsubstitution for $(\Delta x)_1$ is {\green 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}{Alternate Second Step;\\ Stable Reduction}
{\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 & {\cyan -Z_mE} &  {\cyan -X_m}  & 0
\cr 0  & {\cyan Z_v}  & {\cyan -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}{}
\subsection*{Summary: Path-following;\\ {\red 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}{(Nearest) Euclidean Distance Matrix Completion using SDP}

{\bf Given}:\\
 {\cyan \em pre-distance matrix} $A\in\Sn$ (nonnegative with zero
diagonal)\\
 {\green \em weight matrix} $H\in \Sn$:

\[
\mbox{({\red NEDM})}\quad  \mu^*= \min  \frac12  \normF{H\circ (A-D)}^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 \}$, 
{\yellow $k$} is {\yellow \em embedding dimension}

{\green $\circ$} denotes {\green Hadamard (elementwise) matrix product}
\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),~ \normF{Y} \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, but problem size becomes
large)


\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$\\ 

\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<v,\svec(S)\right>
\end{eqnarray*}


\end{slide}
\begin{slide}{Characterization of EDM using SDP}


\begin{center}
$D$ is \EDM \quad ($\subset \Sn$)\\
~\\
{\red
\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}{adjoint/generalized inverse}


with partition:
\[
   D=\left[\matrix{\alpha & d^T\cr d& {\bar D} \cr}\right],
\]
where $\alpha \in \R$
\[
{\cyan \LL^*}(D)=2\left(\Diag(d)+\Diag(\bar De) -\bar D\right)
\]
\[
{\red \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}
\end{slide}
\begin{slide}{Algorithm continued 1}
{\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}{Algorithm continued 2}


{\bf $\bullet$ while} 
$\min \{ \frac {\rm gap}{{\rm objval}+1} , {\rm objval} \} >
              \epsilon_1$
\begin{quote}
{\bf $\bullet$$\bullet$  solve lss for search direction}\\
$\ldots$\\
{\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}{Explicit Preconditioning}

{\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}{
density .0005:.001:.003, CPU times and nnz($\Lambda$), n=200
}




\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
{\cyan 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}
THEOREM
{\cyan
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{slide}{}
\begin{slide}{Proof}
optimality conditions after differentiation
{\cyan
\[
\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$.
{\black \epr}

\end{slide}{}
\begin{slide}{Efficient/Explicit Solution if $y \geq 0$}
large class ({\red generic?}) can be solved in polytime.
\vspace{.2in}

COROLLARY
{\cyan 
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\@.
}

\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 the Theorem.
{\black \epr}

\end{slide}
\begin{slide}{Numerics: 
dim vs dens with\\ $\#$ of failures in $100$ tests
}
though $y = \A^\dagger b\geq 0$  does {\red not} hold in general,
we still get\\
{\cyan  a distance matrix $D$, i.e.  $\A^* (y) \succeq 0$}. \\
$n=10:10:100$; density $.1:.1:.8$.
\[
\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}

%\bibliography{.master,.psd,.publs,.bjorBOOK}
\end{document}