share/tensor/schwarz.dem

   1 /*
   2  * This program is free software; you can redistribute it and/or
   3  * modify it under the terms of the GNU General Public License as
   4  * published by the Free Software Foundation; either version 2 of
   5  * the License, or (at your option) any later version.
   6  *
   7  * This program is distributed in the hope that it will be
   8  * useful, but WITHOUT ANY WARRANTY; without even the implied
   9  * warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  10  * PURPOSE.  See the GNU General Public License for more details.
  11  *
  12  * Obtaining the Schwarzschild metric with MAXIMA's CTENSOR
  13  *
  14  */
  15 kill(all);
  16
  17 ("Finding the Schwarzschild solution of the Einstein vacuum equations" )$
  18
  19 showtime:all;
  20 if get('ctensor,'version)=false then load(ctensor);
  21 (ratwtlvl:false,ratfac:true);
  22 ("Specify the dimension of the manifold and the coordinate labels.")$
  23 (dim:4,ct_coords:[r,th,ph,t]);
  24 ("Enter the general static spherically symmetric metric.")$
  25 lg:matrix([a,0,0,0],[0,r^2,0,0],[0,0,r^2*sin(th)^2,0],[0,0,0,-d]);
  26 ug:invert(lg)$
  27 ("Specify functional dependencies")$
  28 depends([a,d],r);
  29 ("computes inverse metric and specifies diagonality")$
  30 ug:invert(lg);
  31 ("computes the mixed Christoffel symbols but not display them")$
  32 christof(false);
  33 ("computes and ratsimps Ricci tensor")$
  34 uricci(false);
  35 ("computes and displays the Einstein tensor")$
  36 einstein(true);
  37
  38 ("List of the non-zero components of the rank 2 Einstein tensor")$
  39 exp:findde(ein,2);
  40
  41 ("Now begin to solve the field equations")$
  42 exp1:ode2(last(exp),a,r);
  43
  44 ("A kludge to get the solution (the 1,1 component) explicitly")$
  45 solve(exp1,r);
  46 resultlist:solve(%,a)$
  47 h:ev(part(resultlist,1),eval);
  48 ("To cast the solution into standard form")$
  49 h1:h,exp(%c) = 1/(2*m),factor;
  50 ("Now to find the 4,4 component")$
  51 ev(first(exp),h1,diff,factor);
  52 ode2(num(%),d,r);
  53 expand(radcan(%));
  54 h2:ev(%,%c = 1);
  55 ("h1 and h2 should be the solution and to check")$
  56 sol:[h1,h2];
  57 exp,sol,diff,ratsimp;
  58
  59 kill(all);
  60 if get('ctensor,'version)=false then load(ctensor);
  61 (dim:4,ct_coords:[r,th,ph,t]);
  62
  63 ("Enter the Schwarzschild metric in standard coordinates.")$
  64 lg:matrix([1/(1-2*m/r),0,0,0],[0,r^2,0,0]
  65   ,[0,0,r^2*sin(th)^2,0],[0,0,0,(2*m/r-1)]);
  66 ug:invert(lg)$
  67
  68 ("Compute and display mixed Christoffel symbols")$
  69 christof(all);
  70 uricci(true);
  71 ("computes scalar curvature")$
  72 scurvature();
  73 ("computes Riemann tensor")$
  74 lriemann(true);
  75 ("computes contravariant Riemann tensor")$
  76 uriemann(false);
  77 ("computes the Kretchmann invariant Rijkl^2")$
  78 rinvariant();
  79 diagmetric:true;
  80 ("Compute the covariant form of geodesic equations")$
  81 cgeodesic(true);
  82
  83 block([title: "Schwarzschild Potential for Mass M=2",m:2.],
  84 plot3d([r*cos(th),r*sin(th),ug[1,1]],[r,.4,4.],[th,-%pi,%pi],['grid,50,15]));
  85 showtime:false;
  86 /* End of demo -- comment line needed by MAXIMA to resume demo menu */