Figure1 for Model

partn[t_,t1_,k_,β_]=Integrate[Cos[k(t-t1)](((t1))^(β-1)),t1];
partn2only[t_,t0_,k_,β_]=(partn[t,t0,k,β]-partn[t,1,k,β]);
Clear[t,kp]
kp[t_]=10(((t))^(-1/4));
n=1;
c=(10^(-4));
magspec[k_,t_]=((c(((k/kp[t]))^4))/(((1+(((k/kp[t]))^(n(4+5/3)))))^(1/n)));
stress[k_,t_]:=NIntegrate[magspec[k1,t]magspec[√((k^2)+(k1^2)-2k k1 μ),t]((k^2)/((k^2)+(k1^2)-2k k1 μ))(1+(μ^2))(1+(((k-k1 μ))^2)/((k^2)+(k1^2)-2k k1 μ))((2π/4π)),{k1,0.001,(10^5)},{μ,-1,1}];
findLogDivisions[{xmin_,xmax_},n_Integer]:=10^FindDivisions[Log10@{xmin,xmax},n];
interval=findLogDivisions[{(10^(-4)),5000},200];
dd=Dimensions[interval];
plotsh1=Table[{interval[[i]],stress[interval[[i]],1]},{i,1,dd[[1]]}];
fitfunc1=Interpolation[plotsh1];
c3[x_]=fitfunc1[x];

t0={2,4,10};
β=-1/4;
kp1=20;

plot1=LogLogPlot[{36c3[k] (partn2only[t0[[1]],t0[[1]],k,β])^2,36c3[k] (partn2only[t0[[2]],t0[[2]],k,β])^2,36c3[k](partn2only[t0[[3]],t0[[3]],k,β])^2},{k,10^-2,10^3},PlotRange->{{0.01,1000},{10^-13,10^-6}},GridLines->{{1,kp1},None},GridLinesStyle->{{Black},None},PlotStyle->{{Red,Opacity[0.5]},{Blue,Opacity[0.5]},{Black,Opacity[0.5]},{Blue,Dashed}},Frame->True,FrameStyle->Directive[Black,Bold,12],FrameLabel->{Style["k",Black,Bold,12],Style["Sp((OverDot[h]))",Black,Bold,12]},PlotLegends->Placed[LineLegend[{Style["t=2",8],Style["t=4",8],Style["t=10",8]},LegendFunction->"Frame",LegendMarkerSize->8],{0.12,0.75}],Epilog->{Text["(a)",{-3,-14.5}]}];
plot2=LogLogPlot[{36c3[k](partn2only[5t0[[1]],5t0[[1]],k,β])^2,36c3[k] (partn2only[5t0[[2]],5t0[[2]],k,β])^2,36c3[k] (partn2only[3t0[[3]],3t0[[3]],k,β])^2},{k,10^-2,10^3},PlotRange->{{0.01,1000},{10^-13,10^-6}},GridLines->{{1,kp1},None},GridLinesStyle->{{Black},None},PlotStyle->{{Red,Opacity[0.5]},{Blue,Opacity[0.5]},{Black,Opacity[0.5]},{Blue,Dashed}},Frame->True,FrameStyle->Directive[Black,Bold,12],FrameLabel->{Style["k",Black,Bold,12],Style["Sp((OverDot[h]))",Black,Bold,12]},PlotLegends->Placed[LineLegend[{Style["t=10",8],Style["t=20",8],Style["t=30",8]},LegendFunction->"Frame",LegendMarkerSize->8],{0.10,0.20}],Epilog->{Text["(b)",{-3,-15}]}];

SetDirectory["/home/ramkishor/Desktop/GWspectra_from_mhd_turbulence (Version 1229)/mathematica/figures"];
Export["model.eps",plot1];
Export["model_fig2.eps",plot2];