The Curvature of the Sasaki Metric
The Curvature of the Sasaki Metric
This is a modification of Mathematica notebook written to compute the MTW tensor for various cost functions. We have altered it to instead compute the bisectional curvature of a Kahler Sasaki metric. By default, this computes an example we discuss in the paper "The Kähler geometry of certain optimal transport problems". However, it can be easily edited for other examples as well.
First clear any values that may already have been assigned to the names of the various objects to be calculated. The names of the coordinates that you use are also cleared.
In[]:=
Clear[xcoord,potential,Hess,InverseHess,costxxy,costxyy,costxxyy,TwoHBC,TwoQuadHBC,n,TwoHolomorphicSectional]
First we set the dimension.
In[]:=
n=2
Out[]=
2
Note that it is important not to use the symbols, i, j, k, l, p, q, r, s or n as constants or coordinates in the costs and coordinates. The reason is that the first five of those symbols are used as summation or table indices in the calculations done below, and n is the dimension of the space.
In[]:=
xcoord={x1,x2}
Out[]=
{x1,x2}
Here we define the potential function. If you would like to experiment with other potentials, you can just define a new potential, and then run all the of remaining lines of code again.
Here we define the potential function. If you would like to experiment with other potentials, you can just define a new potential, and then run all the of remaining lines of code again.
In[]:=
potential=-x1^2/(4*x2)-Log[-2x2]
Out[]=
--Log[-2x2]
2
x1
4x2
In[]:=
In[]:=
Plot3D[potential,{x1,-1,1},{x2,-2,2}]
Out[]=
Computing the Hessian and its inverse
Computing the Hessian and its inverse
This is where we define the Hessian and its inverse.
In[]:=
Hess:=Hess=Simplify[Table[D[potential,xcoord[[i]],xcoord[[j]]],{i,1,n},{j,1,n}]]
In[]:=
InverseHess:=InverseHess=Simplify[Inverse[Hess]]
If you want to see the Hess and InverseHess explicitly, you can run the following lines.
In[]:=
MatrixForm[Hess]
Out[]//MatrixForm=
- 1 2x2 | x1 2 2 x2 |
x1 2 2 x2 | - 2 x1 2 3 x2 |
In[]:=
MatrixForm[InverseHess]
Out[]//MatrixForm=
2 x1 | x1x2 |
x1x2 | 2 x2 |
It is also worthwhile to make sure that the Hessian is actually non-singular, as this is important for to ensure that the metric is positive definite.
In[]:=
Simplify[Det[Hess]]
Out[]=
-
1
2
3
x2
In[]:=
Simplify[Log[Det[Hess]]]
Out[]=
Log-
1
2
3
x2
Computing the third and fourth derivatives
Computing the third and fourth derivatives
In order to calculate the relevant curvature tensor, we must also define two tensors of third derivatives. For simplicity, we have not changed the names of these tensors from the original MTW notebook
In[]:=
costxxy:=costxxy=Simplify[Table[D[potential,xcoord[[i]],xcoord[[j]],xcoord[[k]]],{i,1,n},{j,1,n},{k,1,n}]]
In[]:=
costxyy:=costxyy=Simplify[Table[D[potential,xcoord[[i]],xcoord[[j]],xcoord[[k]]],{i,1,n},{j,1,n},{k,1,n}]]
We must also define a tensor of fourth order derivatives.
In[]:=
costxxyy:=costxxyy=Simplify[Table[D[potential,xcoord[[i]],xcoord[[j]],xcoord[[k]],xcoord[[l]]],{i,1,n},{j,1,n},{k,1,n},{l,1,n}]]
It is possible to display cxxy,cxyy, and cxxyy, but for most practical purposes, they will be completely unwieldy.
Computing the anti-bisectional curvature
Computing the anti-bisectional curvature
With all of these, we can define the anti-bisectional curvature tensor of the Sasaki metric, which is a fourth order tensor. Note that this calculation is off from the actual anti-bisectional curvature by a factor of 2, which is why we have named it TwoHBC.
In[]:=
TwoHBC=Simplify[Table[Sum[(Sum[costxxy[[i,j,p]]*InverseHess[[p,q]]*costxyy[[q,r,s]],{p,1,n},{q,1,n}]-costxxyy[[i,j,r,s]])*InverseHess[[r,k]]*InverseHess[[s,l]],{r,1,n},{s,1,n}],{i,1,n},{j,1,n},{k,1,n},{l,1,n}]]
Out[]=
{{1,0},{0,0}},-,-1,{-1,0},-,-1,{-1,0},,0,{0,-2}
3x1
x2
3x1
x2
3
2
x1
2
x2
In this form, the tensor is not particularly useful, and to calculate whether it is non-negative definite, it is instead useful to calculate the anti-bisectional curvature for a pair of vectors v and w.
If you would like to specify what v and w are exactly, you can use the following two lines to specify them exactly. Otherwise, you can leave them as is to obtain the full quadratic curvature tensor.
In[]:=
v={Cos[θ],Sin[θ]}
Out[]=
{Cos[θ],Sin[θ]}
In[]:=
w={Sin[ϕ],-Cos[ϕ]}
Out[]=
{Sin[ϕ],-Cos[ϕ]}
Notice here that v and w are vectors in 2 dimensional Euclidean space. However, it is geometrically more natural to consider v as a real vector and w as a real covector, which we extend to their corresponding (1,0) types. Doing so, we can compute the anti bisectional curvature.
In[]:=
TwoQuadHBC=Simplify[Sum[TwoHBC[[i,j,k,l]]*v[[i]]*v[[j]]*w[[k]]*w[[l]],{i,1,n},{j,1,n},{k,1,n},{l,1,n}]]
Out[]=
-2+-++Sin[2θ]Sin[2ϕ]
2
Cos[ϕ]
2
Sin[θ]
2
Cos[θ]
2
Sin[ϕ]
6x1Cos[θ]Sin[θ]
2
Sin[ϕ]
x2
3
2
x1
2
Sin[θ]
2
Sin[ϕ]
2
x2
This does not have a sign.
In[]:=
TwoHolomorphicSectional=Simplify[Sum[TwoHBC[[i,j,k,l]]*v[[i]]*v[[j]]*v[[k]]*v[[l]],{i,1,n},{j,1,n},{k,1,n},{l,1,n}]]
Out[]=
4
Cos[θ]
6x1Sin[θ]
3
Cos[θ]
x2
3
2
x1
2
Cos[θ]
2
Sin[θ]
2
x2
2
Sin[θ]
For optimal transport, the relevant curvature tensor is R( u, OverBar[v],u,Overbar[v]) for an orthogonal (1,0) vector-covector pair. Doing so, we find the following.
In[]:=
vee={1,a}
Out[]=
{1,a}
In[]:=
vee2={a,-1}
Out[]=
{a,-1}
In[]:=
MTW=Simplify[Sum[TwoHBC[[i,j,k,l]]*vee[[i]]*vee[[j]]*vee2[[k]]*vee2[[l]],{i,1,n},{j,1,n},{k,1,n},{l,1,n}]]
Out[]=
3
2
a
2
(-ax1+x2)
2
x2
This is non-negative definite, so the associated cost is weakly regular.
It is also of interest to compute the Ricci curvature. Doing so, we find that it is non-positive.
In[]:=
Ricci=Simplify[Table[-D[Log[Det[Hess]],xcoord[[i]],xcoord[[j]]],{i,1,n},{j,1,n}]]
Acknowledgements
Acknowledgements
This notebook was written by Gabriel Khan. It is adapted from James B. Hartle’s notebook for computing curvature of a Riemannian metric.
This notebook was written by Gabriel Khan. It is adapted from James B. Hartle’s notebook for computing curvature of a Riemannian metric.