From c1b477a1b2b796617c4e345a7296a8d429d7a067 Mon Sep 17 00:00:00 2001 From: Ian Jauslin Date: Sun, 26 Feb 2023 18:36:05 -0500 Subject: Update to v0.4 feature: compute the 2-point correlation function in easyeq. feature: compute the Fourier transform of the 2-point correlation function in anyeq and easyeq. feature: compute the local maximum of the 2-point correlation function and its Fourier transform. feature: compute the compressibility for anyeq. feature: allow for linear spacing of rho's. feature: print the scattering length. change: ux and uk now return real numbers. fix: error in the computation of the momentum distribution: wrong definition of delta functions. fix: various minor bugs. optimization: assign explicit types to variables. --- src/potentials.jl | 58 +++++++++++++++++++++++++++++++++++++++++++------------ 1 file changed, 46 insertions(+), 12 deletions(-) (limited to 'src/potentials.jl') diff --git a/src/potentials.jl b/src/potentials.jl index 46cafc0..b480db0 100644 --- a/src/potentials.jl +++ b/src/potentials.jl @@ -1,4 +1,4 @@ -## Copyright 2021 Ian Jauslin +## Copyright 2021-2023 Ian Jauslin ## ## Licensed under the Apache License, Version 2.0 (the "License"); ## you may not use this file except in compliance with the License. @@ -13,10 +13,15 @@ ## limitations under the License. # exponential potential in 3 dimensions -@everywhere function v_exp(k,a) +@everywhere function v_exp( + k::Float64, + a::Float64 +) return 8*pi/(1+k^2)^2*a end -@everywhere function a0_exp(a) +@everywhere function a0_exp( + a::Float64 +) if a>0. return log(a)+2*MathConstants.eulergamma+2*besselk(0,2*sqrt(a))/besseli(0,2*sqrt(a)) elseif a<0. @@ -27,15 +32,24 @@ end end # exp(-x)-a*exp(-b*x) in 3 dimensions -@everywhere function v_expcry(k,a,b) +@everywhere function v_expcry( + k::Float64, + a::Float64, + b::Float64 +) return 8*pi*((1+k^2)^(-2)-a*b*(b^2+k^2)^(-2)) end -@everywhere function a0_expcry(a,b) +@everywhere function a0_expcry( + a::Float64, + b::Float64 +) return 1.21751642717932720441274114683413710125487579284827462 #ish end # x^2*exp(-|x|) in 3 dimensions -@everywhere function v_npt(k) +@everywhere function v_npt( + k::Float64 +) return 96*pi*(1-k^2)/(1+k^2)^4 end @everywhere function a0_npt() @@ -43,7 +57,9 @@ end end # 1/(1+x^4/4) potential in 3 dimensions -@everywhere function v_alg(k) +@everywhere function v_alg( + k::Float64 +) if(k==0) return 4*pi^2 else @@ -53,32 +69,50 @@ end a0_alg=1. #ish # (1+a x^4)/(1+x^2)^4 potential in 3 dimensions -@everywhere function v_algwell(k) +@everywhere function v_algwell( + k::Float64 +) a=4 return pi^2/24*exp(-k)*(a*(k^2-9*k+15)+k^2+3*k+3) end a0_algwell=1. #ish # potential corresponding to the exact solution c/(1+b^2x^2)^2 -@everywhere function v_exact(k,b,c,e) +@everywhere function v_exact( + k::Float64, + b::Float64, + c::Float64, + e::Float64 +) if k!=0 return 48*pi^2*((18+3*sqrt(c)-(4-3*e/b^2)*c-(1-2*e/b^2)*c^1.5)/(4*(3+sqrt(c))^2*sqrt(c))*exp(-sqrt(1-sqrt(c))*k/b)+(-18+3*sqrt(c)+(4-3*e/b^2)*c-(1-2*e/b^2)*c^1.5)/(4*(3-sqrt(c))^2*sqrt(c))*exp(-sqrt(1+sqrt(c))*k/b)+(1-k/b)/2*exp(-k/b)-c*e/b^2*(3*(9-c)*k/b+8*c)/(8*(9-c)^2)*exp(-2*k/b))/k else return 48*pi^2*(-sqrt(1-sqrt(c))/b*(18+3*sqrt(c)-(4-3*e/b^2)*c-(1-2*e/b^2)*c^1.5)/(4*(3+sqrt(c))^2*sqrt(c))-sqrt(1+sqrt(c))/b*(-18+3*sqrt(c)+(4-3*e/b^2)*c-(1-2*e/b^2)*c^1.5)/(4*(3-sqrt(c))^2*sqrt(c))-1/b-c*e/b^2*(27-19*c)/(8*(9-c)^2)) end end -@everywhere function a0_exact(b,c,e) +@everywhere function a0_exact( + b::Float64, + c::Float64, + e::Float64 +) return 1. #ish end # tent potential (convolution of soft sphere with itself): a*pi/12*(2*|x|/b-2)^2*(2*|x|/b+4) for |x|