The finite positive square root preserves operator order

First suppose $P,Q$ are positive definite and $P⪯Q$. By [F2], their eigenvalues are strictly positive: applying the positive quadratic form to an eigenvector gives $\lambda \mathrm{\parallel }v{\mathrm{\parallel }}^2>0$. Define $P^{-1\mathrm{/}2}$ spectrally by replacing each eigenvalue $\lambda$ with ${\lambda }^{-1\mathrm{/}2}$, and set $P^{-1}:=(P^{-1\mathrm{/}2}{)}^2$; multiplication in the eigenbasis shows these are the inverse square root and inverse of $P$. Put $C:=P^{-1\mathrm{/}2}Q{P}^{-1\mathrm{/}2}$. By the adjoint product rule in [F7], $C$ is self-adjoint, and [F7] gives $⟨v,Cv⟩=⟨P^{-1\mathrm{/}2}v,Q{P}^{-1\mathrm{/}2}v⟩>0$ for $v\mathrm{\ne }0$, so $C$ is positive definite. Also $C-I=P^{-1\mathrm{/}2}(Q-P)P^{-1\mathrm{/}2}⪰0$ by the same quadratic-form identity and [F1]. Hence every eigenvalue $\mu$ of $C$ is at least one by [F2]; define $C^{-1}$ spectrally by reciprocating its positive eigenvalues, giving $C^{-1}⪯I$. Since $Q=P^{1\mathrm{/}2}C{P}^{1\mathrm{/}2}$, multiplying this factorization by $P^{-1\mathrm{/}2}{C}^{-1}{P}^{-1\mathrm{/}2}$ on either side gives the identity, so $Q^{-1}=P^{-1\mathrm{/}2}{C}^{-1}{P}^{-1\mathrm{/}2}$. Congruence of $C^{-1}⪯I$ by $P^{-1\mathrm{/}2}$ now yields $Q^{-1}⪯P^{-1}$.

For $t>0$ and nonzero $v$, $⟨v,(A+tI)v⟩\ge t\mathrm{\parallel }v{\mathrm{\parallel }}^2>0$, and likewise for $B$, so $A+tI$ and $B+tI$ are positive definite and ordered. Step 1.1 gives $(A+tI{)}^{-1}⪰(B+tI{)}^{-1}$. Multiplying the order by $t$ and subtracting from $I$ yields $A(A+tI{)}^{-1}=I-t(A+tI{)}^{-1}⪯I-t(B+tI{)}^{-1}=B(B+tI{)}^{-1}\mathrm{.}$

Set $f(s)=s^{-1\mathrm{/}2}\mathrm{/}(1+s)$ for $s>0$. It is continuous, hence integrable on every compact interval in its domain by [F10]. For $0<s\le 1$, $f(s)\le s^{-1\mathrm{/}2}$; for $s\ge 1$, $f(s)\le s^{-3\mathrm{/}2}$. For each integer $k\ge 1$, the sequence $a_k={\int }_{1\mathrm{/}k}^{1}f(s)ds$ is nondecreasing by positivity and additivity, and is bounded above by ${\int }_{1\mathrm{/}k}^1{s}^{-1\mathrm{/}2}ds=2(1-k^{-1\mathrm{/}2})\le 2$. For $k>1$ this evaluation follows from the primitive $2\sqrt{s}$, whose derivative is $s^{-1\mathrm{/}2}$ by [F4], and the fundamental theorem in [F4]; for $k=1$ both sides are zero by the equal-endpoint convention in [F5]. The sequence $b_k={\int }_1^{k}f(s)ds$ is likewise nondecreasing and bounded above by ${\int }_1^k{s}^{-3\mathrm{/}2}ds=2(1-k^{-1\mathrm{/}2})\le 2$, using the primitive $-2s^{-1\mathrm{/}2}$ and [F4] when $k>1$, and [F5] when $k=1$. Both sequences therefore converge by [F9]. By additivity [F5], $c_k:={\int }_{1\mathrm{/}k}^{k}f(s)ds=a_k+b_k$ converges to a finite $c$. For every $k\ge 2$, monotonicity [F8] and $f(s)\ge 1\mathrm{/}(3\sqrt{2})$ on $[1,2]$ give $c_k\ge {\int }_1^{2}f(s)ds\ge 1\mathrm{/}(3\sqrt{2})$, so $c>0$. For each fixed $\lambda >0$, define $c_{k,\lambda }:={\int }_{1\mathrm{/}(k\lambda )}^{k\mathrm{/}\lambda }f(s)ds$. If $m\ge 2$ is fixed and $k$ is large enough that $1\mathrm{/}(k\lambda )\le 1\mathrm{/}m$ and $k\mathrm{/}\lambda \ge m$, additivity and monotonicity of the integral with the two displayed tail bounds give $c_m\le c_{k,\lambda }=c_m+{\int }_{1\mathrm{/}(k\lambda )}^{1\mathrm{/}m}f(s)ds+{\int }_m^{k\mathrm{/}\lambda }f(s)ds\le c_m+2(m^{-1\mathrm{/}2}-(k\lambda {)}^{-1\mathrm{/}2})+2(m^{-1\mathrm{/}2}-(\lambda \mathrm{/}k{)}^{1\mathrm{/}2})\le c_m+\frac{4}{\sqrt{m}}\mathrm{.}$ Since $c_m\to c$ and $4\mathrm{/}\sqrt{m}\to 0$, this proves $c_{k,\lambda }\to c$ for each fixed $\lambda >0$. Now [F6] with $t=\lambda s$ in the coefficient integral over the fixed interval $[1\mathrm{/}k,k]$ gives ${\int }_{1\mathrm{/}k}^k{t}^{-1\mathrm{/}2}\frac{\lambda }{\lambda +t}dt=\sqrt{\lambda }{c}_{k,\lambda }⟶c\sqrt{\lambda };$ for $\lambda =0$ the coefficient is identically zero. For a nonnegative $X$, let $P_{\lambda }$ be its spectral projection at eigenvalue $\lambda$ and define $J_k(X):={\sum }_{\lambda \in \mathrm{spec}(X)}\left({\int }_{1\mathrm{/}k}^k{t}^{-1\mathrm{/}2}\frac{\lambda }{\lambda +t}dt\right)P_{\lambda }\mathrm{.}$ The coefficient depends only on $\lambda$, so this definition is independent of choices of eigenvectors within eigenspaces. The coefficient convergence and finite spectral sum imply $⟨v,J_k(X)v⟩\to c⟨v,X^{1\mathrm{/}2}v⟩$ for every vector $v$. For $k\ge 2$, this quadratic form is also the integral of the continuous scalar function $t\mapsto ⟨v,t^{-1\mathrm{/}2}X(X+tI{)}^{-1}v⟩$ on $[1\mathrm{/}k,k]$, as follows from the finite spectral sum. By step 2.1 these functions for $X=A$ and $X=B$ are pointwise ordered; [F8] therefore gives $⟨v,J_k(A)v⟩\le ⟨v,J_k(B)v⟩$. Taking limits and using $c>0$ yields $⟨v,A^{1\mathrm{/}2}v⟩\le ⟨v,B^{1\mathrm{/}2}v⟩$ for every $v$, which is $A^{1\mathrm{/}2}⪯B^{1\mathrm{/}2}$.