Numerical analysis - KL model - Heat current - Gapless- for fixed n different temperature

(*In this section we analysed the linear temperature dependence of quantum current in the gapless regime*)

dω=0.00025;ωtab=Table[i,{i,-N[1.8],N[1.8],dω}];
ωntab=Drop[ωtab,{Position[ωtab,0.][[1,1]]}];
zl[m_, r_,ω_,k_,v1_,γ_]:= -m (r^2) (ω^2)+k (v1^2)- I *(γ *ω);
zr[m_, r_,ω_,k_,v2_,γ_]:= -m (r^2) (ω^2)+k (v2^2)- I *(γ *ω);
z[m_,r_,ω_,k_,v1_,v2_]:= -m (r^2) (ω^2)+k ((v1^2) + (v2^2));
m=1.0;r=2.0;a=1.0; k=1.0;γL=0.2;γR=0.2;
f[ω_,T_]:= (1/(Exp[((ℏ ω)/(kb T))]-1));kb=1.0;ℏ=1.0; n=15;

Θ=N[0.0];

v1 = (((2r Sin[Θ]+a)r Cos[Θ])/(Sqrt[(a^2)+ 4 (r^2) ((Cos[Θ])^2)]));
v2=  (((2r Sin[Θ]-a)r Cos[Θ])/(Sqrt[(a^2)+ 4 (r^2) ((Cos[Θ])^2)]));

v1

v2

(*dtl1=0.025;
tl1tab=Table[i,{i,0.0101,0.1001,dtl1}];
dtl2=0.1;
tl2tab=Table[i,{i,0.1001,2.1001,dtl2}];
tltab=Join[tl1tab,tl2tab];
dtr1=0.025;
tr1tab=Table[i,{i,0.01,0.10,dtr1}];
dtr2=0.1;
tr2tab=Table[i,{i,0.10,2.10,dtr2}];
trtab=Join[tr1tab,tr2tab];
Ttab=((tltab+trtab)/2);
tltab-trtab;*)

(*Ttab*)

(*dtl=0.1;tltab=Table[i,{i,2.9,3.0,dtl}];trtab=Table[i,{i,1.9,2.0,dtl}];
Ttab=((tltab+trtab)/2);*)

dtl=0.015;tltab=Table[i,{i,0.051,0.11,dtl}];trtab=Table[i,{i,0.05,0.1,dtl}];
Ttab=((tltab+trtab)/2);

Ttab

GList={};
transcofflist={};
currlist={};
currclasslist={};
currQlist={};
gmatslist={};
mats= Table[0,{n},{n}];
For[l=1,l<=Length[ωntab],l++,
For[i=1,i<=n-1,i++,
mats[[i,i]]=z[m,r,ωntab[[l]],k,v1,v2];mats[[1,1]]=zl[m,r,ωntab[[l]],k,v1,γL];mats[[n,n]]=zr[m,r,ωntab[[l]],k,v2,γR];mats[[i,i+1]]=-k * Abs[v1 v2];mats[[i+1,i]]=-k * Abs[v1 v2]; Z=mats;];GList=Join[GList,{{ωntab[[l]],Inverse[Z]}}];gmatslist=Join[gmatslist,{{ωntab[[l]],((Abs[GList[[l,2,1,n]]])^2)}}];
transcofflist=Join[transcofflist,{{ωntab[[l]],(((Abs[GList[[l,2,1,n]]])^2))*((ωntab[[l]])^2)*γL*γR}}];]

For[p=1,p<=Length[Ttab],p++,For[l=1,l<=Length[ωntab],l++,currlist=Join[currlist,{{Ttab[[p]],((Abs[GList[[l,2,1,n]]])^2)*((ℏ ((ωntab[[l]])^3))/π)*γL*γR*(f[ωntab[[l]],tltab[[p]]]-f[ωntab[[l]],trtab[[p]]])}}];
currclasslist=Join[currclasslist,{{Ttab[[p]],((Abs[GList[[l,2,1,n]]])^2)*(((ωntab[[l]])^2)/π)*γL*γR*(tltab[[p]]-trtab[[p]])}}];
currQlist=Join[currQlist,{{Ttab[[p]],((Abs[GList[[l,2,1,n]]])^2)*(((ℏ^2) ((ωntab[[l]])^4))/((((2*kb*Ttab[[p]]))^2)*π))*γL*γR*kb*(tltab[[p]]-trtab[[p]])*(Csch[((ℏ ωntab[[l]])/(2kb*Ttab[[p]]))])*(Csch[((ℏ ωntab[[l]])/(2kb*Ttab[[p]]))])}}];]]

τ={};
For[j=1,j<=Length[Ttab],j++,τ=Join[τ,{{Ttab[[j]],Sum[currQlist[[i,2]]*dω,{i,(j-1)*Length[ωntab]+1,j*Length[ωntab]}]}}]]

ν={};
For[j=1,j<=Length[Ttab],j++,ν=Join[ν,{{Ttab[[j]],Sum[currlist[[i,2]]*dω,{i,(j-1)*Length[ωntab]+1,j*Length[ωntab]}]}}]]

η={};
For[j=1,j<=Length[Ttab],j++,η=Join[η,{{Ttab[[j]],Sum[currclasslist[[i,2]]*dω,{i,(j-1)*Length[ωntab]+1,j*Length[ωntab]}]}}]]

τ

ν

η

(*This plot shows the linear temperature dependence of the quantum current in the gapless regime*)

ListPlot[τ,PlotRange->{0,All}]

(*Linear fitting*)

Fit[τ,{1,x},x]

p1=Show[{ListPlot[τ,PlotRange->{0,All}],Plot[3.19217319469743*(10^(-6))+0.00013748965511223158 x,{x,0,0.1}]}]

(*ListLinePlot[gmatslist,PlotRange->All]*)

(*ListLinePlot[transcofflist,PlotRange->All]*)

(*ListLinePlot[transcofflist,PlotRange->{{-0.7,0.7},{-0.25,0.25}}]*)

Numerical analysis - KL model - Heat current - Gapless - for fixed temperature and different system size

(*In this section we analysed the system size dependence of quantum current in the gapless regime and showed that the transport is ballistic*)

dω=0.00025;ωtab=Table[i,{i,-N[1.8],N[1.8],dω}];
ωntab=Drop[ωtab,{Position[ωtab,0.][[1,1]]}];
zl[m_, r_,ω_,k_,v1_,γ_]:= -m (r^2) (ω^2)+k (v1^2)- I *(γ *ω);
zr[m_, r_,ω_,k_,v2_,γ_]:= -m (r^2) (ω^2)+k (v2^2)- I *(γ *ω);
z[m_,r_,ω_,k_,v1_,v2_]:= -m (r^2) (ω^2)+k ((v1^2) + (v2^2));
m=1.0;r=2.0;a=1.0; k=1.0;γL=0.2;γR=0.2;
f[ω_,T_]:= (1/(Exp[((ℏ ω)/(kb T))]-1));kb=1.0;ℏ=1.0;

ntab={2,4,8,16,24,32,64};

Θ=N[0.0];

v1 = (((2r Sin[Θ]+a)r Cos[Θ])/(Sqrt[(a^2)+ 4 (r^2) ((Cos[Θ])^2)]));
v2=  (((2r Sin[Θ]-a)r Cos[Θ])/(Sqrt[(a^2)+ 4 (r^2) ((Cos[Θ])^2)]));

v1

v2

(*dtl1=0.025;
tl1tab=Table[i,{i,0.0101,0.1001,dtl1}];
dtl2=0.1;
tl2tab=Table[i,{i,0.1001,2.1001,dtl2}];
tltab=Join[tl1tab,tl2tab];
dtr1=0.025;
tr1tab=Table[i,{i,0.01,0.10,dtr1}];
dtr2=0.1;
tr2tab=Table[i,{i,0.10,2.10,dtr2}];
trtab=Join[tr1tab,tr2tab];
Ttab=((tltab+trtab)/2);
tltab-trtab;*)

(*Ttab*)

(*dtl=0.1;tltab=Table[i,{i,2.9,3.0,dtl}];trtab=Table[i,{i,1.9,2.0,dtl}];
Ttab=((tltab+trtab)/2);*)

dtl=0.015;tltab=Table[i,{i,0.051,0.061,dtl}];trtab=Table[i,{i,0.05,0.06,dtl}];
Ttab=((tltab+trtab)/2);

Ttab=((tltab+trtab)/2);

Ttab

GListtab={};
transcofflisttab={};
currlisttab={};
currclasslisttab={};
currQlisttab={};
gmatslisttab={};
matstab= Table[Table[0,{ntab[[i]]},{ntab[[i]]}],{i,1,Length[ntab]}];
Ztab=Table[Table[0,{ntab[[i]]},{ntab[[i]]}],{i,1,Length[ntab]}];

(*For[l=1,l<=Length[ωntab],l++,For[j=1,j<=Length[ntab],j++,
For[i=1,i<=ntab[[j]]-1,i++,
matstab[[j]][[i,i]]=z[m,r,ωntab[[l]],k,v1,v2];matstab[[j]][[1,1]]=zl[m,r,ωntab[[l]],k,v1,γ];matstab[[j]][[ntab[[j]],ntab[[j]]]]=zr[m,r,ωntab[[l]],k,v2,η];matstab[[j]][[i,i+1]]=-k * Abs[v1 v2];
matstab[[j]][[i+1,i]]=-k * Abs[v1 v2]];Ztab=Join[Ztab, {matstab[[j]]}];];]*)

For[l=1,l<=Length[ωntab],l++,For[j=1,j<=Length[ntab],j++,
For[i=1,i<=ntab[[j]]-1,i++,
matstab[[j]][[i,i]]=z[m,r,ωntab[[l]],k,v1,v2];matstab[[j]][[1,1]]=zl[m,r,ωntab[[l]],k,v1,γL];matstab[[j]][[ntab[[j]],ntab[[j]]]]=zr[m,r,ωntab[[l]],k,v2,γR];matstab[[j]][[i,i+1]]=-k * Abs[v1 v2];
matstab[[j]][[i+1,i]]=-k * Abs[v1 v2]];
Ztab[[j]]=matstab[[j]];
GListtab=Join[GListtab,{{ωntab[[l]],Inverse[Ztab[[j]]]}}];];]

combineByLabel[list_]:={#[[1,1]],#[[All,2]]}&/@GatherBy[list,First];

GListtabnew=combineByLabel[GListtab];

For[p=1,p<=Length[Ttab],p++,For[j=1,j<=Length[ntab],j++,For[l=1,l<=Length[ωntab],l++,currQlisttab=Join[currQlisttab,{{ntab[[j]],{Ttab[[p]],((Abs[GListtabnew[[l,2,j,1,ntab[[j]]]]])^2)*(((ℏ^2) ((ωntab[[l]])^4))/((((2*kb*Ttab[[p]]))^2)*π))*γL*γR*kb*(tltab[[p]]-trtab[[p]])*(Csch[((ℏ ωntab[[l]])/(2kb*Ttab[[p]]))])*(Csch[((ℏ ωntab[[l]])/(2kb*Ttab[[p]]))])}}}];currclasslisttab=Join[currclasslisttab,{{ntab[[j]],{Ttab[[p]],((Abs[GListtabnew[[l,2,j,1,ntab[[j]]]]])^2)*(((ωntab[[l]])^2)/π)*γL*γR*(tltab[[p]]-trtab[[p]])}}}];]]]

currQlisttabnew={#[[1,1]],#[[All,2]]}&/@GatherBy[currQlisttab,First];

currclasslisttabnew= {#[[1,1]],#[[All,2]]}&/@GatherBy[currclasslisttab,First];

(*Definition of quantum current*)

tauForSize[sizeIndex_]:=Table[{Ttab[[j]],Sum[currQlisttabnew[[sizeIndex,2,i,2]]*dω,{i,(j-1)*Length[ωntab]+1,j*Length[ωntab]}]},{j,Length[Ttab]}]

tauAll=Table[{currQlisttabnew[[s,1]],(*system size n*)tauForSize[s]},{s,Length[currQlisttabnew]}];

tauAll

(*Definition of classical current*)

etaForSize[sizeIndex_]:=Table[{Ttab[[j]],Sum[currclasslisttabnew[[sizeIndex,2,i,2]]*dω,{i,(j-1)*Length[ωntab]+1,j*Length[ωntab]}]},{j,Length[Ttab]}]

etaAll=Table[{currclasslisttabnew[[s,1]],(*system size n*)etaForSize[s]},{s,Length[currclasslisttabnew]}];

etaAll

(*This plot shows the ballistic nature of the quantum current in the gapless regime*)

ListPlot[Table[{tauAll[[i]][[1]],tauAll[[i,2,1,2]]},{i,1,Length[ntab]}],PlotRange->{0,All}]

(*This plot shows the ballistic nature of the classical current in the gapless regime*)

ListPlot[Table[{etaAll[[i]][[1]],etaAll[[i,2,1,2]]},{i,1,Length[ntab]}],PlotRange->{0,All}]

Numerical analysis - KL model - Heat current - Gapped- for fixed n different temperature

(*In this section we analysed the linear temperature dependence of quantum current in the gapped regime*)

dω=0.00025;ωtab=Table[i,{i,-N[1.8],N[1.8],dω}];
ωntab=Drop[ωtab,{Position[ωtab,0.][[1,1]]}];
zl[m_, r_,ω_,k_,v1_,γ_]:= -m (r^2) (ω^2)+k (v1^2)- I *(γ *ω);
zr[m_, r_,ω_,k_,v2_,γ_]:= -m (r^2) (ω^2)+k (v2^2)- I *(γ *ω);
z[m_,r_,ω_,k_,v1_,v2_]:= -m (r^2) (ω^2)+k ((v1^2) + (v2^2));
m=1.0;r=2.0;a=1.0; k=1.0;γL=0.2;γR=0.2;
f[ω_,T_]:= (1/(Exp[((ℏ ω)/(kb T))]-1));kb=1.0;ℏ=1.0; n=16;

Θ=N[(π/3)];

v1 = (((2r Sin[Θ]+a)r Cos[Θ])/(Sqrt[(a^2)+ 4 (r^2) ((Cos[Θ])^2)]));
v2=  (((2r Sin[Θ]-a)r Cos[Θ])/(Sqrt[(a^2)+ 4 (r^2) ((Cos[Θ])^2)]));

v1

v2

ω0=(Sqrt[(k/(m (r^2)))])(v1-v2)

(*dtl1=0.025;
tl1tab=Table[i,{i,0.0101,0.1001,dtl1}];
dtl2=0.1;
tl2tab=Table[i,{i,0.1001,2.1001,dtl2}];
tltab=Join[tl1tab,tl2tab];
dtr1=0.025;
tr1tab=Table[i,{i,0.01,0.10,dtr1}];
dtr2=0.1;
tr2tab=Table[i,{i,0.10,2.10,dtr2}];
trtab=Join[tr1tab,tr2tab];
Ttab=((tltab+trtab)/2);
tltab-trtab;*)

(*Ttab*)

(*dtl=0.1;tltab=Table[i,{i,2.9,3.0,dtl}];trtab=Table[i,{i,1.9,2.0,dtl}];
Ttab=((tltab+trtab)/2);*)

dtl=0.015;tltab=Table[i,{i,0.051,0.11,dtl}];trtab=Table[i,{i,0.05,0.1,dtl}];
Ttab=((tltab+trtab)/2);

Ttab

GList={};
transcofflist={};
currlist={};
currclasslist={};
currQlist={};
gmatslist={};
mats= Table[0,{n},{n}];
For[l=1,l<=Length[ωntab],l++,
For[i=1,i<=n-1,i++,
mats[[i,i]]=z[m,r,ωntab[[l]],k,v1,v2];mats[[1,1]]=zl[m,r,ωntab[[l]],k,v1,γL];mats[[n,n]]=zr[m,r,ωntab[[l]],k,v2,γR];mats[[i,i+1]]=-k * Abs[v1 v2];mats[[i+1,i]]=-k * Abs[v1 v2]; Z=mats;];GList=Join[GList,{{ωntab[[l]],Inverse[Z]}}];gmatslist=Join[gmatslist,{{ωntab[[l]],((Abs[GList[[l,2,1,n]]])^2)}}];
transcofflist=Join[transcofflist,{{ωntab[[l]],(((Abs[GList[[l,2,1,n]]])^2))*((ωntab[[l]])^2)*γL*γR}}];]

For[p=1,p<=Length[Ttab],p++,For[l=1,l<=Length[ωntab],l++,currlist=Join[currlist,{{Ttab[[p]],((Abs[GList[[l,2,1,n]]])^2)*((ℏ ((ωntab[[l]])^3))/π)*γL*γR*(f[ωntab[[l]],tltab[[p]]]-f[ωntab[[l]],trtab[[p]]])}}];
currclasslist=Join[currclasslist,{{Ttab[[p]],((Abs[GList[[l,2,1,n]]])^2)*(((ωntab[[l]])^2)/π)*γL*γR*(tltab[[p]]-trtab[[p]])}}];
currQlist=Join[currQlist,{{Ttab[[p]],((Abs[GList[[l,2,1,n]]])^2)*(((ℏ^2) ((ωntab[[l]])^4))/((((2*kb*Ttab[[p]]))^2)*π))*γL*γR*kb*(tltab[[p]]-trtab[[p]])*(Csch[((ℏ ωntab[[l]])/(2kb*Ttab[[p]]))])*(Csch[((ℏ ωntab[[l]])/(2kb*Ttab[[p]]))])}}];]]

τ={};
For[j=1,j<=Length[Ttab],j++,τ=Join[τ,{{Ttab[[j]],Sum[currQlist[[i,2]]*dω,{i,(j-1)*Length[ωntab]+1,j*Length[ωntab]}]}}]]

ν={};
For[j=1,j<=Length[Ttab],j++,ν=Join[ν,{{Ttab[[j]],Sum[currlist[[i,2]]*dω,{i,(j-1)*Length[ωntab]+1,j*Length[ωntab]}]}}]]

η={};
For[j=1,j<=Length[Ttab],j++,η=Join[η,{{Ttab[[j]],Sum[currclasslist[[i,2]]*dω,{i,(j-1)*Length[ωntab]+1,j*Length[ωntab]}]}}]]

τ

ν

η

(*This plot shows the temperature dependence of the quantum current in the gapped regime*)

ListLogPlot[τ,PlotRange->{0,All}]

(*In this plot we compared the temperature dependence of quantum current with the analytic expression*)

p2=Show[LogPlot[0.0000085*((Exp[((-ℏ ω0)/(kb T))])/(Sqrt[T])),{T,0.001,0.1},PlotRange->All],ListLogPlot[τ,PlotRange->{0,All}]]

(*ListLinePlot[gmatslist,PlotRange->All]*)

(*ListLinePlot[transcofflist,PlotRange->All]*)

(*ListLinePlot[transcofflist,PlotRange->{{-0.7,0.7},{-0.25,0.25}}]*)

Numerical analysis - KL model - Heat current - Gapped - for fixed temperature and different system size

(*In this section we analysed the system size dependence of quantum current in the gapped regime and showed that the transport is exponentially suppressed*)

dω=0.00025;ωtab=Table[i,{i,-N[1.8],N[1.8],dω}];
ωntab=Drop[ωtab,{Position[ωtab,0.][[1,1]]}];
zl[m_, r_,ω_,k_,v1_,γ_]:= -m (r^2) (ω^2)+k (v1^2)- I *(γ *ω);
zr[m_, r_,ω_,k_,v2_,γ_]:= -m (r^2) (ω^2)+k (v2^2)- I *(γ *ω);
z[m_,r_,ω_,k_,v1_,v2_]:= -m (r^2) (ω^2)+k ((v1^2) + (v2^2));
m=1.0;r=2.0;a=0.5; k=1.0;γL=0.2;γR=0.2;
f[ω_,T_]:= (1/(Exp[((ℏ ω)/(kb T))]-1));kb=1.0;ℏ=1.0;

ntab={2,4,8,16,24,32,64};

Θ=N[(π/3)];

v1 = (((2r Sin[Θ]+a)r Cos[Θ])/(Sqrt[(a^2)+ 4 (r^2) ((Cos[Θ])^2)]));
v2=  (((2r Sin[Θ]-a)r Cos[Θ])/(Sqrt[(a^2)+ 4 (r^2) ((Cos[Θ])^2)]));

v1

v2

ω0=(Sqrt[(k/(m (r^2)))])(v1-v2)

(*dtl1=0.025;
tl1tab=Table[i,{i,0.0101,0.1001,dtl1}];
dtl2=0.1;
tl2tab=Table[i,{i,0.1001,2.1001,dtl2}];
tltab=Join[tl1tab,tl2tab];
dtr1=0.025;
tr1tab=Table[i,{i,0.01,0.10,dtr1}];
dtr2=0.1;
tr2tab=Table[i,{i,0.10,2.10,dtr2}];
trtab=Join[tr1tab,tr2tab];
Ttab=((tltab+trtab)/2);
tltab-trtab;*)

(*Ttab*)

(*dtl=0.1;tltab=Table[i,{i,2.9,3.0,dtl}];trtab=Table[i,{i,1.9,2.0,dtl}];
Ttab=((tltab+trtab)/2);*)

dtl=0.015;tltab=Table[i,{i,0.011,0.021,dtl}];trtab=Table[i,{i,0.01,0.02,dtl}];
Ttab=((tltab+trtab)/2);

Ttab=((tltab+trtab)/2);

Ttab

GListtab={};
transcofflisttab={};
currlisttab={};
currclasslisttab={};
currQlisttab={};
gmatslisttab={};
matstab= Table[Table[0,{ntab[[i]]},{ntab[[i]]}],{i,1,Length[ntab]}];
Ztab=Table[Table[0,{ntab[[i]]},{ntab[[i]]}],{i,1,Length[ntab]}];

(*For[l=1,l<=Length[ωntab],l++,For[j=1,j<=Length[ntab],j++,
For[i=1,i<=ntab[[j]]-1,i++,
matstab[[j]][[i,i]]=z[m,r,ωntab[[l]],k,v1,v2];matstab[[j]][[1,1]]=zl[m,r,ωntab[[l]],k,v1,γ];matstab[[j]][[ntab[[j]],ntab[[j]]]]=zr[m,r,ωntab[[l]],k,v2,η];matstab[[j]][[i,i+1]]=-k * Abs[v1 v2];
matstab[[j]][[i+1,i]]=-k * Abs[v1 v2]];Ztab=Join[Ztab, {matstab[[j]]}];];]*)

For[l=1,l<=Length[ωntab],l++,For[j=1,j<=Length[ntab],j++,
For[i=1,i<=ntab[[j]]-1,i++,
matstab[[j]][[i,i]]=z[m,r,ωntab[[l]],k,v1,v2];matstab[[j]][[1,1]]=zl[m,r,ωntab[[l]],k,v1,γL];matstab[[j]][[ntab[[j]],ntab[[j]]]]=zr[m,r,ωntab[[l]],k,v2,γR];matstab[[j]][[i,i+1]]=-k * Abs[v1 v2];
matstab[[j]][[i+1,i]]=-k * Abs[v1 v2]];
Ztab[[j]]=matstab[[j]];
GListtab=Join[GListtab,{{ωntab[[l]],Inverse[Ztab[[j]]]}}];];]

combineByLabel[list_]:={#[[1,1]],#[[All,2]]}&/@GatherBy[list,First];

GListtabnew=combineByLabel[GListtab];

For[p=1,p<=Length[Ttab],p++,For[j=1,j<=Length[ntab],j++,For[l=1,l<=Length[ωntab],l++,currQlisttab=Join[currQlisttab,{{ntab[[j]],{Ttab[[p]],((Abs[GListtabnew[[l,2,j,1,ntab[[j]]]]])^2)*(((ℏ^2) ((ωntab[[l]])^4))/((((2*kb*Ttab[[p]]))^2)*π))*γL*γR*kb*(tltab[[p]]-trtab[[p]])*(Csch[((ℏ ωntab[[l]])/(2kb*Ttab[[p]]))])*(Csch[((ℏ ωntab[[l]])/(2kb*Ttab[[p]]))])}}}];currclasslisttab=Join[currclasslisttab,{{ntab[[j]],{Ttab[[p]],((Abs[GListtabnew[[l,2,j,1,ntab[[j]]]]])^2)*(((ωntab[[l]])^2)/π)*γL*γR*(tltab[[p]]-trtab[[p]])}}}];]]]

currQlisttabnew={#[[1,1]],#[[All,2]]}&/@GatherBy[currQlisttab,First];

currclasslisttabnew= {#[[1,1]],#[[All,2]]}&/@GatherBy[currclasslisttab,First];

(*Definition of quantum current*)

tauForSize[sizeIndex_]:=Table[{Ttab[[j]],Sum[currQlisttabnew[[sizeIndex,2,i,2]]*dω,{i,(j-1)*Length[ωntab]+1,j*Length[ωntab]}]},{j,Length[Ttab]}]

tauAll=Table[{currQlisttabnew[[s,1]],(*system size n*)tauForSize[s]},{s,Length[currQlisttabnew]}];

tauAll

(*Definition of classical current*)

etaForSize[sizeIndex_]:=Table[{Ttab[[j]],Sum[currclasslisttabnew[[sizeIndex,2,i,2]]*dω,{i,(j-1)*Length[ωntab]+1,j*Length[ωntab]}]},{j,Length[Ttab]}]

etaAll=Table[{currclasslisttabnew[[s,1]],(*system size n*)etaForSize[s]},{s,Length[currclasslisttabnew]}];

etaAll

(*This plot shows the insulating nature of the quantum current in the gapped regime*)

ListPlot[Table[{tauAll[[i]][[1]],tauAll[[i,2,1,2]]},{i,1,Length[ntab]}],PlotRange->All]

ListLogPlot[Table[{tauAll[[i]][[1]],tauAll[[i,2,1,2]]},{i,1,Length[ntab]}],PlotRange->All]

data=Table[{tauAll[[i]][[1]],tauAll[[i,2,1,2]]},{i,1,Length[ntab]}]

nlm=NonlinearModelFit[data,c* Exp[-b x],{c,b},x];

nlm["BestFitParameters"]

Show[Plot[Log[nlm[x]],{x,0,64},PlotRange->All],ListLogPlot[Table[{tauAll[[i]][[1]],tauAll[[i,2,1,2]]},{i,1,Length[ntab]}]]]

(*This plot shows the ballistic nature of the classical current in the gapped regime*)

ListPlot[Table[{etaAll[[i]][[1]],etaAll[[i,2,1,2]]},{i,1,Length[ntab]}],PlotRange->{0,All}]

(*Plot[0.00015*Exp[((-2(x-1)Abs[v1-v2])/(a(Sqrt[v1 v2])))],{x ,0,25}]*)

(*((-2Abs[v1-v2])/(a(Sqrt[v1 v2])))*)

(*Show[{Plot[Log[0.000085*Exp[((-2(x-1)Abs[v1-v2])/(a(Sqrt[v1 v2])))]],{x ,0,64},PlotRange->All],ListLogPlot[Table[{tauAll[[i]][[1]],tauAll[[i,2,1,2]]},{i,1,Length[ntab]}],PlotRange->{0,All}]}]*)