CosmoLattice output

tanh4_Ex1.h model: V(ϕ) = (1/4)(Λ^4)(Tanh^4)(ϕ/M)

(* Created by Francisco Torrenti - June 2022, updated 2025 *)

(*This notebook reads the output data generated by CosmoLattice (in text format) and creates a series of useful figures. It is adapted to the default tanh4_Ex1.h model implemented in CosmoLattice, but can be easily generalized to other scalar models. *)

Clear["Global`*"] 
(*Clears all variables*)

SetDirectory[NotebookDirectory[]]; 
(*Sets working directory to current directory*)

Scale factor

(*Column structure of file "average_scale_factor.txt"*)

strm=OpenRead["average_scale_factor.txt"];  (*Opens file to read data from*)
(*Skip[strm,String,1];*)  (* If you add headers to the text files, uncomment this to skip the first line *)
NumColumns=4;  (*Number of columns in file*)
scale=ReadList[strm,Table[Number,{NumColumns}]];  (*Reads data and creates table*)

FinalTime=scale[[Length[scale]]][[1]];  (*Final time of simulation*)

a((OverTilde[η]))

sfIn=Interpolation[scale[[All,{1,2}]]];   
(*Creates interpolation from columns 1 and 2*)

Plot[{sfIn[t]},{t,0,FinalTime},
Frame->True,
FrameLabel->{Style[OverTilde[η],FontSize->16,Black],Style["a",FontSize->16,Black]},
FrameTicksStyle->Directive[Black,14],
PlotStyle->{Darker[Blue],Thick},
ImageSize->400
]

a’((OverTilde[η]))

sfdIn=Interpolation[scale[[All,{1,3}]]];  
(*Creates interpolation from columns 1 and 3*)

Plot[{sfdIn[t]},{t,0,FinalTime},
Frame->True,PlotRange->All,
FrameLabel->{Style[OverTilde[η],FontSize->16,Black],Style["a'",FontSize->16,Black]},
FrameTicksStyle->Directive[Black,14],
PlotStyle->{Darker[Blue],Thick},
ImageSize->400
]

a’((OverTilde[η]))/a((OverTilde[η]))

sfddIn=Interpolation[scale[[All,{1,4}]]];  
(*Creates interpolation from columns 1 and 4*)

Plot[{sfddIn[t]},{t,0,FinalTime},
Frame->True,PlotRange->All,
FrameLabel->{Style[OverTilde[η],FontSize->16,Black],Style["a'/a",FontSize->16,Black]},
FrameTicksStyle->Directive[Black,14],
PlotStyle->{Darker[Blue],Thick},
ImageSize->400
]

Volume-averages of field 0 (inflaton (OverTilde[ϕ]))

; (*Column structure of file "average_scalar_0.txt"*)

strm=OpenRead["average_scalar_0.txt"]; (*Opens file to read data from*)
NumColumns=7;  (*Number of columns in file*)

means0=Flatten[ReadList[strm,Table[Number,NumColumns],RecordLists->True],1]; (*Reads data and creates table*)

(*NOTE: Sometimes "nans" get produced in the last two columns, which the ReadList function cannot read. They must be transformed to "0" in the text file to make it work.*)

f0phiIn=Interpolation[means0[[All,{1,2}]]];
(*Creates interpolation from columns 1 and 2*)

(* Volume average of inflation *)
Plot[f0phiIn[t],{t,0,FinalTime},ImageSize->400,
PlotStyle->{{Thick,Darker[Blue]}},Frame->True,PlotRange->All,FrameLabel->{Style["(OverTilde[η])",FontSize->16,Black],Style[OverTilde["<|ϕ|>"],FontSize->16,Black]},FrameTicksStyle->Directive[Black,14]
]

(* Volume average of inflation --> conformal transformation*)

Plot[f0phiIn[t]*sfIn[t],{t,0,FinalTime},ImageSize->400,
PlotStyle->{{Thick,Darker[Blue]}},Frame->True,PlotRange->All,FrameLabel->{Style[(OverTilde[η]),FontSize->16,Black],Style[a OverTilde["<|ϕ|>"],FontSize->16,Black]},FrameTicksStyle->Directive[Black,14]
]

rms0In=Interpolation[Table[{means0[[i]][[1]],means0[[i]][[6]]},{i,1,Length[means0]}]];
(*Creates interpolation from columns 1 and 6*)

(* rms of inflaton (including conformal transformation) *)
LogPlot[sfIn[t]*rms0In[t],{t,0,FinalTime},ImageSize->400,
PlotStyle->{{Thick,Darker[Blue]}},Frame->True,PlotRange->All,FrameLabel->{Style[(OverTilde[η]),FontSize->16,Black],Style["a*rms((OverTilde[ϕ]))",FontSize->16,Black]},FrameTicksStyle->Directive[Black,14]
]

“Energy conservation” (preservation of 1st Friedmann equation)

;

strm=OpenRead["average_energy_conservation.txt"];
NumColumns=4;
EnCons=ReadList[strm,Table[Number,{NumColumns}]];

EnConsIn=Interpolation[Table[{EnCons[[i,1]],EnCons[[i,2]]},{i,1,Length[EnCons]}]];
LogPlot[Abs[EnConsIn[t]],{t,0,FinalTime},PlotRange->{10^(-10),1},ImageSize->400,
PlotStyle->{{Thick,Darker[Blue]}},Frame->True,FrameLabel->{Style[(OverTilde[η]),FontSize->16,Black],Style["(Subscript[Δ, a])",FontSize->16,Black]},FrameTicksStyle->Directive[Black,14]
]

Energy density components

;

strm=OpenRead["average_energies.txt"]; (*Opens file to read data from*)
NumColumns=5; (*Number of columns in file*)
energy=ReadList[strm,Table[Number,{NumColumns}]];

(* Creates an interpolating function for each energy component*)

Kin0=Interpolation[energy[[All,{1,2}]]]; (*Kinetic inflation*)
Grad0=Interpolation[energy[[All,{1,3}]]]; (*Gradient inflaton *)
Pot=Interpolation[energy[[All,{1,4}]]]; (*Inflaton potential *)
ETot=Interpolation[energy[[All,{1,5}]]];(* Sum of components*)

resc=4; (* Coefficient for conformal rescaling *)

(* Plot of energy components with conformal rescaling *)
LogLogPlot[{Kin0[t]*sfIn[t]^resc,Grad0[t]*sfIn[t]^resc,Pot[t]*sfIn[t]^resc,ETot[t]*sfIn[t]^resc},{t,0,FinalTime},

PlotStyle->{{Thick,Darker[Blue]},{Thick,Lighter[Lighter[Blue]]},{Thick,Darker[Green]},{Thick,Black}},

ImageSize->450,PlotRange->{10^(-15),10},Frame->True,FrameLabel->{Style[OverTilde[η],FontSize->16,Black],Style["(a^4)"OverTilde["ρ"],FontSize->16,Black]},FrameTicksStyle->Directive[Black,16]]

(* Same plot as before but in linear scale*)

Plot[{Kin0[t]*sfIn[t]^resc,Grad0[t]*sfIn[t]^resc,Pot[t]*sfIn[t]^resc,ETot[t]*sfIn[t]^resc},{t,0,FinalTime},PlotStyle->{{Thick,Darker[Blue]},{Thick,Lighter[Lighter[Blue]]},{Thick,Darker[Green]},{Thick,Black}},ImageSize->450,PlotRange->{0,0.82},Frame->True,FrameLabel->{Style[OverTilde[η],FontSize->16,Black],Style["(a^4)"OverTilde["ρ"],FontSize->16,Black]},FrameTicksStyle->Directive[Black,16]]

Parameters for spectra outputs

(* Parameters from the simulation *)
  Npoints=32;   (* N: number of lattice points per side *)
  tOutputInfreq=5;  (* time interval between spectra *)
  NumBins=IntegerPart[Sqrt[3]* Npoints/2]; (* number of bins in the spectra *)

(* Parameters for spectra output *)
  s3=5; (* frequency of spectra output. Spectra printed are s3*n+1 for n=1, 2,... *)

(*Number of spectra to be printed*)
NumSpectra=FinalTime/(tOutputInfreq*s3)

colors=Table[ColorData["Rainbow", t],{t,1,0,-1/NumSpectra}]; (*creates a rainbow of colours for spectra figures*)

specTimes=Table[tOutputInfreq*s3*(i-1),{i,1,NumSpectra*s3}]; (*approximate spectra times*)

Spectra Fields

(*Column structure of file "spectra_scalar_X.txt":*)

;

Inflaton

NumColumns=4;
specS0=ReadList["spectra_scalar_0.txt",Table[Number,{NumColumns}]];

Table[
spectraS0[i]=Table[{specS0[[j]][[1]],specS0[[j]][[2]]},{j,NumBins*s3*(i-1) + 1,NumBins*s3*(i-1)+NumBins}],{i,1,NumSpectra}];
Table[
spectraS0d[i]=Table[{specS0[[j]][[1]],specS0[[j]][[3]]},{j,NumBins*s3*(i-1) + 1,NumBins*s3*(i-1)+NumBins}],{i,1,NumSpectra}];
Table[
partnumb0[i]=Table[{specS0[[j]][[1]],specS0[[j]][[4]]},{j,NumBins*s3*(i-1) + 1,NumBins*s3*(i-1)+NumBins}],{i,1,NumSpectra}];

(*Spectrum of inflaton amplitude*)
ListLogLogPlot[Table[spectraS0[i],{i,1,NumSpectra}],PlotStyle->colors,Joined->True,PlotRange->{(10^(-20)),10^(-4)},ImageSize->440,Frame -> True,FrameLabel->{Style["(OverTilde[k])",FontSize->14,Black],Style["(Subscript[(OverTilde[Δ]), (OverTilde[ϕ])])((OverTilde[k]))",FontSize->14,Black]},FrameTicksStyle->Directive[Black,14]]

(*Spectrum of inflaton time-derivative*)
ListLogLogPlot[Table[spectraS0d[i],{i,1,NumSpectra}],PlotStyle->colors,Joined->True,PlotRange->{(10^(-20)),10^(-4)},ImageSize->440,Frame -> True,FrameLabel->{Style["(OverTilde[k])",FontSize->14,Black],Style["(Subscript[(OverTilde[Δ]), (OverTilde[ϕ])'])((OverTilde[k]))",FontSize->14,Black]},FrameTicksStyle->Directive[Black,14]]