WOLFRAM NOTEBOOK

In[]:=
deploy
Tue 8 Oct 2019 14:14:46

Cheap version of B3

In[]:=
SeedRandom[1];d=5;normalize[x_]:=x/Total[x,10];u=normalize[RandomReal[{0,1},d]];Unprotect[D];D=DiagonalMatrix[u];A=D-Outer[Times,u,u];cheap=DiagonalMatrix[Sqrt[u]]-Outer[Times,u,
u
];Assert[cheap.cheapA]

Compare against symmetric decomposition

In[]:=
SeedRandom[1];Unprotect[D];d=5;normalize[x_]:=x/Total[x,10];u=normalize[RandomReal[{0,1},d]];D=DiagonalMatrix[u];A=D+uu;uu=Outer[Times,u,u];(*uu'*)iD=DiagonalMatrix[1/u];(*
-1
D
*)x=
u.iD.u+1
-1
u.iD.u
;B2=
1/2
D
+x
1/4
iD
.uu.
1/4
iD
;B3=
1/2
D
+xuu.
1/2
iD
;Print["error exact=",Norm[MatrixPower[symsqrt[A],2]-A]]Print["error rank1=",Norm[B2-symsqrt[A]]]Print["error rank1 fixed=",Norm[B3-symsqrt[A]]]Print["error diag=",Norm[MatrixPower[D,1/2]-symsqrt[A]]]
error exact=3.99455×
-16
10
error rank1=0.00876574
error rank1 fixed=0.0414601
error diag=0.219797

Errors of two approximations

In[]:=
Clear[doit];doit[idx_]:=Module{},d=5;u=normalize[RandomReal[{0,1},d]];D=DiagonalMatrix[u];A=D+uu;uu=Outer[Times,u,u];(*uu'*)iD=DiagonalMatrix[1/u];(*
-1
D
*)x=
u.iD.u+1
-1
u.iD.u
;B2=
1/2
D
+x
1/4
iD
.uu.
1/4
iD
;B3=
1/2
D
+xuu.
1/2
iD
;truth=symsqrt[A];{Norm[B2-truth],Norm[B3-truth]};vals=doit/@Range[1000000];Count[Thread[vals[[All,1]]<vals[[All,2]]],True]
Out[]=
897509
Wolfram Cloud

You are using a browser not supported by the Wolfram Cloud

Supported browsers include recent versions of Chrome, Edge, Firefox and Safari.


I understand and wish to continue anyway »

You are using a browser not supported by the Wolfram Cloud. Supported browsers include recent versions of Chrome, Edge, Firefox and Safari.