# Machine Assistance

• Here is Mathematica code written over the last few years and optimized in 2006 to fit in 6 lines which computes the Laplacian L = d d* + d* d for an arbitrary discrete geometry G. The frame work used here is the language of simplicial complexes. They are even simpler than graphs. A simplicial complex is just a set of non-empty finite sets which have the property that it is closed under the operation of taking non-empty finite subsets. Graphs define simplicial complexes, the so called Whitney complex: look at all complete subgraphs and take the set of vertices of each of this subgraphs. This is a simplicial complex. For a triangular graph for example, the complex is G={{1,2,3},{1,2},{1,3},{2,3},{1},{2},{3}}. There are 3 sets with one elements which are the vertices. Then there are three sets with two elements which are the edges. Then there is one set with three elements which are the triangles. If we write down G as above, we also have already fixed an orientation of each of its elements. A discrete differential form is just a function on G. When restricted to sets of size 1, it is a 0-form, if restricted to sets of size 2, it is a 1-form, if restricted to sets of size 3, it is a 2-form, etc. The Laplacian L leaves these k-forms invariant. On 0-forms, it is a discrete analogue of the classical Laplacian L0 = div(grad) which appears in Gauss's reformulation L0 V = 4π σ of the laws of gravity which leads to the Newton laws (as we have seen in the divergence lecture). Then there is the Laplacian L1 which acts on 1-forms. The equation L1 F = j is completely equivalent to the Maxwell equations if things are done in space-time R4. But we can do the same physics in any discrete network like our small world Gaia from Lecture 37. Here is the code. The procedure R[n,m] generates a random complex. I added the complex of Gaia, the discrete world [PDF] we have looked at in Lecture 37. You have computed three eigenvalues in that lecture (they were all 2). There is no accident that 2 appeared in each sector. There is a general McKean-Singer super symmetry: the non-zero eigenvalues on the even forms are the same than the non-zero eigenvalues on the odd forms. See this paper, where it was noticed to hold in the discrete. The code below computes also the dimension of the kernels of the Li which are called the Betti numbers. The first one b0 tells how many connectivity components the world has. The second one b1 is zero if and only if the world is simply connected. In that case, any vector field F which has zero curl is a gradient field. The story is not only similar, it is exactly the same than the continuum.
```(* "Cohomology in 6 lines", February 16, 2018, Oliver Knill                *)
Generate[A_]:=Delete[Union[Sort[Flatten[Map[Subsets,A],1]]],1]
R[n_,m_]:=Module[{A={},X=Range[n],k},Do[k:=1+Random[Integer,n-1];
A=Append[A,Union[RandomChoice[X,k]]],{m}];Generate[A]]; G=R[5,7];
G={{1,2,3},{2,3,4},{1,2},{1,3},{2,3},{2,4},{3,4},{1},{2},{3},{4}}; (*Gaia  *)
G=Sort[G]; n=Length[G]; Dim=Map[Length,G]-1;f=Delete[BinCounts[Dim],1];
Orient[a_,b_]:=Module[{z,c,k=Length[a],l=Length[b]}, If[SubsetQ[a,b] &&
(k==l+1),z=Complement[a,b][];c=Prepend[b,z];Signature[a]*Signature[c],0]];
d=Table[0,{n},{n}]; d=Table[Orient[G[[i]],G[[j]]],{i,n},{j,n}];
Dirac=d+Transpose[d]; L=Dirac.Dirac; f=Prepend[f,0]; m=Length[f]-1;
U=Table[v=f[[k+1]];Table[u=Sum[f[[l]],{l,k}];L[[u+i,u+j]],{i,v},{j,v}],{k,m}];
Cohomology=Map[NullSpace,U]; Betti=Map[Length,Cohomology]
MatrixForm[L]            (* The full Laplacian                             *)
Map[Eigenvalues,U]       (* The eigenvalues of L_0, L_1 and L_2            *)
```
• Here is one possible computation of the volume π2/2 of the 3-ball x2 + y2 + z2 + w2 ≤ 1 in R4:
```Integrate[1,{x,-1,1},{y,-Sqrt[1-x^2],Sqrt[1-x^2]},
{z,-Sqrt[1-x^2-y^2],Sqrt[1-x^2-y^2]},
{w,-Sqrt[1-x^2-y^2-z^2],Sqrt[1-x^2-y^2-z^2]}]
```
This is the sledge hammer method using Euclidean coordinates which Mathematica can do but which we would have a bit of difficulty with without doing a couple of trig substitutions. Better is to use spherical or Hopf coordinates. (see section 26). The Jacobean determinant is r3 sin(2 φ)/2. Now just integrate over the toral variables θ1 and θ1 as well as φ from 0 to π/2. Even more elegant is to use that 2π2 is the surface area (you have computed that in Homework Problem 26.5), then note that this is A(r) = 2π^2 r3 if the radius parameter r is added. Now integrate ∫0R A(r) dr = 2π2 R4/4 = 2π^2 R4/4 which is for R=1 the expression π^2/2.
• Here is the proof of the ``important formula" in the proof of Stokes theorem:
```Simplify[{Ry-Qz,Pz-Rx,Qx-Py}.Cross[{xu,yu,zu},{xv,yv,zv}]==
{Px*xu+Py*yu+Pz*zu,Qx*xu+Qy*yu+Qz*zu,Rx*xu+Ry*yu+Rz*zu}.{xv,yv,zv}-
{Px*xv+Py*yv+Pz*zv,Qx*xv+Qy*yv+Qz*zv,Rx*xv+Ry*yv+Rz*zv}.{xu,yu,zu}]
```
• Here is the proof that the planimeter vector field has curl 1:
```s = Solve[{(x - a)^2 + (y - b)^2 == 1, a^2 + b^2 == 1}, {a, b}];
{A, B} = First[{a, b} /. s]; F = {-(y - B), x - A};
Simplify[ Curl[F, {x, y}]]
```
• In unit 15, we want to verify that functions solve a PDE. Here is an example showing that the gaussian = normal distribution solves the heat equation:
```f=PDF[NormalDistribution[m,Sqrt[2 t]],x];
Simplify[D[f,t]==D[f,{x, 2}]]
```
This actually also works In Wolfram alpha, but you have to rerun the computation. • For the heat equation, here is a Wolfram Alpha verification The direct code
```f=E^(-(x-m)^2/(t))/(Sqrt[Pi*t]);Simplify[D[f,{x,2}]==4 D[f,t]]
```
does currently not work in Wolfram alpha.
• In unit 10, we have computed some distortion factors. Here is how you can check it with Mathematica. First the case of polar coordinates
```{f1,f2}={r*Cos[t], r*Sin[t]};
df = {{D[f1,r], D[f1,t]},
{D[f2,r], D[f2,t]}};
Simplify[Det[df]]
```
And here are the spherical coordinates:
```{f1, f2, f3} = {r Sin[s] Cos[t], r Sin[s] Sin[t], r Cos[s]};
df = {{D[f1, r], D[f1, s], D[f1, t]},
{D[f2, r], D[f2, s], D[f2, t]},
{D[f3, r], D[f3, s], D[f3, t]}};
Simplify[Det[df]]
```
We urge you still to do this calculation once by hand! You only appreciate the machine assistance if you have gone through the simplifications yourself.
• A special page About HW 12.
• In Units 8 and 9, it is good to be able to plot a curve. Here is an example on how to plot a spiral
```ParametricPlot[ {t Cos[t],t Sin[t]},{t,0,10Pi}]
```
And here is the Wolfram Alpha. Wolfram alpha also right away computes the arc length. In Mathematica, you would do it
```Integrate[ Norm[D[{t Cos[t], t Sin[t]}, t]], {t, 0, 10 Pi}]
```
In this case the integral was known. Numerically you get an answer with
```NIntegrate[ Norm[D[{t Cos[t], t Sin[t]}, t]], {t, 0, 10 Pi}]
```
• In Unit 5, you deal with surfaces having holes. Here is how you draw a bagle
```ContourPlot3D[(3+x^2+y^2+z^2)^2-16(x^2+y^2)==0,{x,-3,3},{y,-3,3},{z,-3,3}]
```
In Wolfram Alpha. or two bagles
```ContourPlot3D[((3+x^2+y^2+z^2)^2-16(x^2+y^2))((3+(x-4)^2+y^2+z^2)^2-16((x-4)^2+y^2))==0,{x,-3,7},{y,-3,3},{z,-3,3}]
```
In Wolfram Alpha.
• In Unit 2, we row reduced. You can get to the unique row reduced echelon form as follows. This is an example with a 3x4 matrix:
``` RowReduce[{{2,2,2,2},{1,1,1,1},{3,4,6,1}}]
```
And here is the implementation in Wolfram Alpha
• In Unit 5, you are asked to glue together some manifolds. Here is an example how you can glue together two spheres. Here is how you plot in Mathematica. If we have a product of two expressions, then one of them has to be zero. We see therefore the union of two spheres. Now change the constant 0 on the right with say 1 to get the two spheres merged.
```ContourPlot3D[(x^2+y^2+z^2-1)((x-1)^2+y^2+z^2-1)==0,{x,-3,4},{y,0,3},{z,-2,2}]
```
```  v={v1,v2,v3}; w={w1,w2,w3};