Plots for the Conformally Invariant Fields plus radiation

<<"C:\\Users\\anderson\\Research\\cosmology\\Order Reduction\\Conformally_invariant_fields\\Plots_for_the_paper\\Solutions_Method_1_corrected.m"

<<"C:\\Users\\anderson\\Research\\cosmology\\Order Reduction\\Conformally_invariant_fields\\Plots_for_the_paper\\Solutions_Method_2.m"

<<"C:\\Users\\anderson\\Research\\cosmology\\Order Reduction\\Conformally_invariant_fields\\Plots_for_the_paper\\Solutions_Method_3_corrected.m"

I am using α and Eqmethod1 for Method 1, α0 for the zeroth iteration of Methods 2 and 3, etc.  I am using Eqmethod2a or b or c for the zeroth and higher iterations of Method 2, and right now I am using Eqa, b, c  for the zeroth and higher iterations of Method 3.  
Further the zeroth iteration of Method 2 is the classical solution.

First I give the plot for the iterations of Method 2.

Plot[{α0[t]/.Eqmethod2a,α1[t-0.76]/.Eqmethod2b,α2[t-0.88]/.Eqmethod2c},{t,-1,1},Frame->True,FrameLabel->{τ,α},FrameTicks->{{{1.4,1.5,1.6,1.7},None},{{-1,-0.5,0,0.5,1},None}},PlotRange->{{-1,1},{1.4,1.75}}]

Next is the plot for the iterations of Method 3.

Plot[{α0[t]/.Eqa,α1[t-0.6]/.Eqb,α2[t-0.665]/.Eqccorrected},{t,-1,1},Frame->True,FrameLabel->{τ,α},FrameTicks->{{{1.4,1.5,1.6,1.7},None},{{-1,-0.5,0,0.5,1},None}},PlotRange->{{-1,1},{1.4,1.7}}]

In the following two plots, the classical solution and the solutions for Method 1, second iteration of Method 2, second iteration of Method 3 are given.

Plot[{α0[t]/.Eqmethod2a,α[t]/.Eqmethod1, α2[t-0.88]/.Eqmethod2c,α2[t-0.665]/.Eqccorrected},{t,-1,1},Frame->True,FrameLabel->{τ,α},FrameTicks->{{{1.4,1.5,1.6,1.7},None},{{-1,-0.5,0,0.5,1},None}},PlotRange->{{-1,1},{1.4,1.75}}]

Plot[{α0[t]/.Eqmethod2a,α[t]/.Eqmethod1, α2[t-0.88]/.Eqmethod2c,α2[t-0.665]/.Eqccorrected},{t,-5,5},Frame->True,FrameLabel->{τ,α},FrameTicks->{{{0,5,10},None},{{-4,-2,0,2,4},None}},PlotRange->{{-5,5},{0,12}}]

The above plots are for the paper.  Next I want to plot the relative differences between the curves for the iterations of Method 2.  First I’ll plot the curves for Method 2 out to larger times.

Plot[{α0[t]/.Eqmethod2a,α1[t-0.76]/.Eqmethod2b,α2[t-0.82]/.Eqmethod2c},{t,5,9}]

Plot[{α0[t]/.Eqmethod2a,α1[t-0.76]/.Eqmethod2b,α2[t-0.82]/.Eqmethod2c},{t,-9,-5}]

So the curves for the 2nd and 3rd iterations cross at some point.

Plot[{((α0[t]/.Eqmethod2a)-(α1[t-0.76]/.Eqmethod2b))/(α0[t]/.Eqmethod2a),((α0[t]/.Eqmethod2a)-(α2[t-0.82]/.Eqmethod2c))/(α0[t]/.Eqmethod2a)},{t,-9,9}]

p1 = Plot[((α1[t-0.76]/.Eqmethod2b)-(α2[t-0.82]/.Eqmethod2c))/(α1[t-0.76]/.Eqmethod2b),{t,-9,9}]

Next I want to plot the relative differences between the curves for the iterations of Method 3.  First I’ll plot the curves for Method 3 out to larger times.

Plot[{α0[t]/.Eqa,α1[t-0.6]/.Eqb,α2[t-0.66]/.Eqccorrected},{t,5,9}]

Plot[{α0[t]/.Eqa,α1[t-0.6]/.Eqb,α2[t-0.66]/.Eqccorrected},{t,-9,-5}]

Again the curves for the 2nd and 3rd iterations cross at some point.

Plot[{((α0[t]/.Eqa)-(α1[t-0.6]/.Eqb))/(α0[t]/.Eqa),((α0[t]/.Eqa)-(α2[t-0.66]/.Eqccorrected))/(α0[t]/.Eqa)},{t,-9,9}]

p2 = Plot[((α1[t-0.6]/.Eqb)-(α2[t-0.66]/.Eqccorrected))/(α1[t-0.6]/.Eqb),{t,-9,9}]

Show[p1,p2]

(* Test to use first not second iteration for M3 *)

Plot[{α0[t]/.Eqmethod2a,α[t]/.Eqmethod1, α2[t-0.88]/.Eqmethod2c,α1[t-0.6]/.Eqb},{t,-1,1},Frame->True,FrameLabel->{τ,α},FrameTicks->{{{1.4,1.5,1.6,1.7},None},{{-1,-0.5,0,0.5,1},None}},PlotRange->{{-1,1},{1.4,1.75}}]

Plot[{α0[t]/.Eqmethod2a,α[t]/.Eqmethod1, α2[t-0.88]/.Eqmethod2c,α2[t-0.66]/.Eqccorrected},{t,-1,1},Frame->True,FrameLabel->{τ,α},FrameTicks->{{{1.4,1.5,1.6,1.7},None},{{-1,-0.5,0,0.5,1},None}},PlotRange->{{-1,1},{1.4,1.75}}]

(* What seems to be happening is that the term proportional to alphaq is smaller than that proportional to betaq and so
one iteration for M3 is very similar to 2 iterations for M2.  If the alpha terms was larger this would presumably not be true.  *)