<< EDCRGTCcode.m

Bondi metric @ Newman-Penrose null frame

Off[General::"spell", General::"spell1"];

xCoord = {u, r, θ, φ};

nullFrame = {d[u], 
          Exp[2 βf[u, r, θ]] (d[r] + 
                V[u, r, θ]/2 
                   d[u]), (r Exp[γf[u, 
                      r, θ]] (d[θ] - 
                      U[u, r, θ] d[u]) + 
                I r Exp[-γf[u, r, θ]] Sin[θ] 
                  d[φ])/(Sqrt[2]), (r Exp[γf[u, 
                      r, θ]] (d[θ] - 
                      U[u, r, θ] d[u]) - 
                I r Exp[-γf[u, r, θ]] Sin[θ] 
                  d[φ])/(Sqrt[2])};

(**  V[u,r,θ] here equals V[u,r,θ]/r in the standard notation. βf and γf are used in place of β, γ to distinguish them from the spin coefficients! **)

gCoordFrame = {{Exp[2 βf[u, r, θ]] V[u, r, θ] - 
            r^2  Exp[2  γf[u, r, θ]] U[u, r, θ]^2, 
          Exp[2 βf[u, r, θ]], 
          r^2  Exp[2  γf[u, r, θ]] U[u, r, θ], 
          0}, {Exp[2 βf[u, r, θ]], 0, 0, 
          0}, {r^2  Exp[2  γf[u, r, θ]] U[u, r, θ], 
          0, -r^2 Exp[2  γf[u, r, θ]], 0}, {0, 0, 
          0, -r^2 Exp[-2 γf[u, r, θ]] 
            Sin[θ]^2}}; MatrixForm[gCoordFrame]

(**  Check line element **)

FacSimp[
      2  nullFrame[[1]] nullFrame[[2]] - 
        2  nullFrame[[3]] nullFrame[[4]] - 
        d[xCoord] . gCoordFrame . d[xCoord]]

gNPnull = {{0, 1, 0, 0}, {1, 0, 0, 0}, {0, 0, 0, -1}, {0, 
            0, -1, 0}};

simpRules = TrigRules;

RGtensors[gNPnull, xCoord, nullFrame, {0, 1}]

expRdd = Expand[Rdd]; // Timing

(**  Number of terms in each component of Rdd **)

Map[Length, expRdd, {2}]

Φ[0, 0]

Φ[0, 1]

FacSimp[Φ[1, 1] + 3  Λ]

Φ[0, 2]

\(
\)

(**  Define Main Eqs - Rest Eqs **)

MainEqs = 
        Expand[{Φ[0, 0], Φ[0, 1], 
            FacSimp[Φ[1, 1] + 3  Λ], Φ[
              0, 2]}];

RestEqs = 
        Expand[{FacSimp[Φ[1, 1] - 
                3  Λ], Φ[2, 1], Φ[2, 
              2]}];

AllEqs = Join[MainEqs, RestEqs];

Length /@ AllEqs

\(
\)

(**  Introduce series expansions **)

gamRul = 
        FuncRepRules[γf[u, r, θ], 
          Series[c1[u, θ]/r + c2[u, θ]/r^2 + 
              c3[u, θ]/r^3 + c4[u, θ]/r^4, {r, Infinity, 
              4}]];

betRul = 
        FuncRepRules[βf[u, r, θ], 
          Series[b1[u, θ]/r + b2[u, θ]/r^2 + 
              b3[u, θ]/r^3 + b4[u, θ]/r^4, {r, Infinity, 
              4}]];

UuRul = 
        FuncRepRules[U[u, r, θ], 
          Series[u1[u, θ]/r + u2[u, θ]/r^2 + 
              u3[u, θ]/r^3 + u4[u, θ]/r^4, {r, Infinity, 
              4}]];

VvRul = 
        FuncRepRules[V[u, r, θ], 
          Series[1 - 2  M[u, θ]/r + v2[u, θ]/r^2 + 
              v3[u, θ]/r^3 + v4[u, θ]/r^4, {r, Infinity, 
              4}]];

approxRules = Join[gamRul, betRul, UuRul, VvRul];

approxMain = FacSimp[MainEqs /. approxRules]; // Timing

approxMain[[1]]

approxMain[[2]]

\(
\)

(**  Redefine series expansions to satisfy leading terms - introduce nN[u,θ] in place of u3[u,θ] **)

betRul=FuncRepRules[βf[u, r, θ],Series[-c1[u, θ]^2/4/r^2-2 c1[u, θ] c2[u, θ]/3/r^3+ c4[u, θ]/r^4,{r, Infinity, 4}]];

UuRul=FuncRepRules[U[u, r, θ],Series[-(2 c1[u, θ] Cot[θ]+Derivative[0, 1][c1][u, θ])/r^2+(2  nN[u, θ]+3c1[u, θ]Derivative[0, 1][c1][u, θ]+4  c1[u, θ]^2  Cot[θ])/r^3+u4[u, θ]/r^4,{r, Infinity, 4}]];

approxRules = Join[gamRul, betRul, UuRul, VvRul];

approxMain = FacSimp[MainEqs /. approxRules]; // Timing

approxMain[[1]]

approxMain[[2]]

(**  Redefine series expansions again to satisfy leading terms: c2[u,θ]=q1[u]/Sin[θ]^2 **)

betRul=FuncRepRules[βf[u, r, θ],Series[-c1[u, θ]^2/4/r^2-2 c1[u, θ] c2[u, θ]/3/r^3-(2  c2[u, θ]^2 + 
                            3 c1[u, θ] c3[u, θ])/4/r^4,{r, Infinity, 4}]];

approxRules=Join[FuncRepRules[c2[u, θ], q1[u]/Sin[θ]^2],gamRul,betRul,UuRul,VvRul];

approxMain = FacSimp[MainEqs //. approxRules]; // Timing

approxMain[[1]]

approxMain[[2]]

(**  Solve for u4[u,θ] and generate list of replacement rules **)

u4sol=Solve[Normal[%] == 0, u4[u, θ]][[1]]

u4rules=FuncRepRules[
          u4[u, θ], (u4[u, θ] /. u4sol)];

approxMain[[3]]

(**  Similarly for v2[u,θ] and v3[u,θ] **)

v2sol=Solve[Coefficient[Normal[%], 1/r^4] == 0, 
          v2[u, θ]][[1]]

v3sol=Solve[Coefficient[Normal[%%], 1/r^5] == 0, 
          v3[u, θ]][[1]]

v2rules=FuncRepRules[
          v2[u, θ], (v2[u, θ] /. v2sol)];

v3rules=FuncRepRules[
          v3[u, θ], (v3[u, θ] /. v3sol)];

approxMain[[4]]

(**  Set q1[u] equal to a constant and collect all rules **)

quvRules=Join[FuncRepRules[q1[u], q1c], u4rules, v2rules, v3rules];

allRules=Join[approxRules, quvRules];

approxMain = FacSimp[MainEqs //. allRules]; // Timing

approxMain[[1]]

approxMain[[2]]

approxMain[[3]]

approxMain[[4]]

(**  Solve for Derivative[1, 0][c3][u, θ] and Derivative[1, 0][c4][u, θ] and generate lists of replacement rules **)

c3uSol=Solve[Coefficient[Normal[%], 1/r^4] == 0,Derivative[1, 0][c3][u, θ]][[1]]

c4uSol=Solve[Coefficient[Normal[%%], 1/r^5] == 0,Derivative[1, 0][c4][u, θ]][[1]]

c3uRules=FuncRepRules[Derivative[1, 0][c3][u, θ],(Derivative[1, 0][c3][u, θ]/.c3uSol)];

c4uRules=FuncRepRules[Derivative[1, 0][c4][u, θ],(Derivative[1, 0][c4][u, θ]/.c4uSol)];

(**  Redefine all rules to include above rules **)

allRules=Join[approxRules, quvRules, c3uRules, c4uRules];

FacSimp[approxMain //. allRules]

\(\(\(
\)\(
\)\)\)

(**  Solve remaining equations **)

approxRest = FacSimp[RestEqs //. allRules]; // Timing

approxRest[[1]]

approxRest[[2]]

nNuSol=Solve[Coefficient[Normal[%], 1/r^3] == 0,Derivative[1, 0][nN][u, θ]][[1]]

approxRest[[3]]

MuSol=Solve[Coefficient[Normal[%], 1/r^2] == 0,Derivative[1, 0][M][u, θ]][[1]]

nNuRules=FuncRepRules[Derivative[1, 0][nN][u, θ],(Derivative[1, 0][nN][u, θ]/.nNuSol)];

MuRules=FuncRepRules[Derivative[1, 0][M][u, θ],(Derivative[1, 0][M][u, θ]/.MuSol)];

(**  Final set of rules **)

allRules=Join[allRules, nNuRules, MuRules];

FacSimp[approxRest //. allRules]

(** Main equations satisfied to O[(1/r)]^6, Rest equations to O[(1/r)]^5 **)

\(\(\(
\)\(
\)\)\)

(**  Weyl tensor components **)

approxPSI = 
        FacSimp[Ψ /@ {0, 1, 2, 3, 4} /. 
            approxRules];

approxPSIsol = FacSimp[approxPSI //. allRules];

approxPSIsol[[1]]

approxPSIsol[[2]]

approxPSIsol[[3]]

approxPSIsol[[4]]

approxPSIsol[[5]]

\(
\)

(**  Form of metric functions when q1c=0 **)

FacSimp[(E^(2 βf[u, r, θ])) //. allRules] /. 
      q1c -> 0

FacSimp[(E^(2 γf[u, r, θ])) //. allRules] /. 
      q1c -> 0

FacSimp[U[u, r, θ] //. allRules] /. q1c -> 0

FacSimp[V[u, r, θ] //. allRules] /. q1c -> 0

\(
\)

Off[General::"spell", General::"spell1"];