Correct smooth cutoff

Author: BrieucLD

ChainOhmT/quadohmT.jl

@@ -74,13 +74,13 @@ function ohmicspectraldensity_finiteT(x,i,α,s,ωc,β; smooth=false)
         y = 0
     elseif i==2
 	if smooth
-	  y = -2 .*( α*abs.(x).^s ./ ωc^(s-1)) .* (coth.((β/2).*x) .+ 1)./2 .* exp(-abs(x)/wc)
+	  y = -2 .*( α*abs.(x).^s ./ ωc^(s-1)) .* (coth.((β/2).*x) .+ 1)./2 .* exp.(-abs.(x) ./ ωc)
 	else
           y = -2 .*( α*abs.(x).^s ./ ωc^(s-1)) .* (coth.((β/2).*x) .+ 1)./2 
 	end
     elseif i==3
         if smooth
-          y = 2 .*( α*abs.(x).^s ./ ωc^(s-1)) .* (coth.((β/2).*x) .+ 1)./2 .* exp(-abs(x)/wc)
+          y = 2 .*( α*abs.(x).^s ./ ωc^(s-1)) .* (coth.((β/2).*x) .+ 1)./2 .* exp.(-abs.(x) ./ ωc)
         else
           y = 2 .*( α*abs.(x).^s ./ ωc^(s-1)) .* (coth.((β/2).*x) .+ 1)./2 
         end

examples/bath-observables.jl

@@ -34,7 +34,7 @@ tfinal = 60.0           # simulation time
 if β == ∞
     cpars = chaincoeffs_ohmic(N, α, s; ωc=ωc)  # chain parameters, i.e. on-site energies ϵ_i, hopping energies t_i, and system-chain coupling c_0
 else
-    cpars = chaincoeffs_finiteT(N, β; α=α, s=s, J=nothing, ωc=ωc, mc=4, mp=0, AB=nothing, iq=1, idelta=2, procedure=:Lanczos, Mmax=5000, save=false)  # chain parameters, i.e. on-site energies ϵ_i, hopping energies t_i, and system-chain coupling c_0
+    cpars = chaincoeffs_finiteT(N, β; α=α, s=s, J=nothing, ωc=ωc, mc=4, mp=0, AB=nothing, iq=1, idelta=2, procedure=:Lanczos, Mmax=5000, save=false, smooth=false)  # chain parameters, i.e. on-site energies ϵ_i, hopping energies t_i, and system-chain coupling c_0
     #=
     #If cpars is stored in "../ChainOhmT/ohmicT" 
     curdir = @__DIR__

examples/puredephasing.jl

@@ -38,7 +38,7 @@ s = 1 # ohmicity
 if β == ∞
     cpars = chaincoeffs_ohmic(N, α, s; ωc=ωc)  # chain parameters, i.e. on-site energies ϵ_i, hopping energies t_i, and system-chain coupling c_0
 else
-    cpars = chaincoeffs_finiteT(N, β; α=α, s=s, J=nothing, ωc=ωc, mc=4, mp=0, AB=nothing, iq=1, idelta=2, procedure=:Lanczos, Mmax=5000, save=true)  # chain parameters, i.e. on-site energies ϵ_i, hopping energies t_i, and system-chain coupling c_0
+    cpars = chaincoeffs_finiteT(N, β; α=α, s=s, J=nothing, ωc=ωc, mc=4, mp=0, AB=nothing, iq=1, idelta=2, procedure=:Lanczos, Mmax=5000, save=false, smooth=false)  # chain parameters, i.e. on-site energies ϵ_i, hopping energies t_i, and system-chain coupling c_0
     #=#If cpars is stored in "../ChainOhmT/ohmicT" 
     curdir = @__DIR__
     dir_chaincoeff = abspath(joinpath(curdir, "../ChainOhmT/ohmicT"))

examples/sbm_tailored_specdensity.jl

@@ -105,7 +105,7 @@ if Jω_type==:Discrete
 
 elseif Jω_type==:Continuous
 
-    cpars = chaincoeffs_finiteT(N, β, false; α=α, s=s, J=Jω_fct, ωc=ωc, mc=size(AB)[1], mp=0, AB=AB, iq=1, idelta=2, procedure=:Lanczos, Mmax=5000, save=false)  # chain parameters, i.e. on-site energies ϵ_i, hopping energies t_i, and system-chain coupling c_0
+    cpars = chaincoeffs_finiteT(N, β, false; α=α, s=s, J=Jω_fct, ωc=ωc, mc=size(AB)[1], mp=0, AB=AB, iq=1, idelta=2, procedure=:Lanczos, Mmax=5000, save=false, smooth=false)  # chain parameters, i.e. on-site energies ϵ_i, hopping energies t_i, and system-chain coupling c_0
 
 end
 

src/finitetemperature.jl

@@ -27,8 +27,9 @@ Users can provide their own spectral density.
 * AB: component intervals. The defaut intervals are [[-Inf -ωc];[-ωc 0];[0 ωc];[ωc Inf]].
 * Mmax: maximum number of integration points
 * save: if true the coefficients are saved
+* smooth: if true the spectral density is multiplied by a decreasing exponential instead of a Heaviside function.
 """
-function chaincoeffs_finiteT(nummodes, β, ohmic=true; α=1, s=1, J=nothing, ωc=1, mc=4, mp=0, AB=nothing, iq=1, idelta=2, procedure=:Lanczos, Mmax=5000, save=true)
+function chaincoeffs_finiteT(nummodes, β, ohmic=true; α=1, s=1, J=nothing, ωc=1, mc=4, mp=0, AB=nothing, iq=1, idelta=2, procedure=:Lanczos, Mmax=5000, save=true, smooth=false)
 
     N = nummodes #Number of bath modes
 
@@ -45,7 +46,7 @@ function chaincoeffs_finiteT(nummodes, β, ohmic=true; α=1, s=1, J=nothing, ωc
 
     # Express the spectral density according to the intervals
     if ohmic==true
-        wf(x,i) = ohmicspectraldensity_finiteT(x,i,α,s,ωc,β)
+        wf(x,i) = ohmicspectraldensity_finiteT(x,i,α,s,ωc,β; smooth=smooth)
     elseif J==nothing
         throw(ArgumentError("A spectral density should have been provided."))
     else