From 4007f2eac2d0f977cff074ff25d0ef61a390b8de Mon Sep 17 00:00:00 2001 From: Nuwan Yapa Date: Tue, 28 Jan 2025 18:03:40 -0500 Subject: [PATCH] Fixed some units --- mesons.jl | 49 ++++++++++++++++++++-------------------- nucleons.jl | 8 +++---- test/Pb208_meson_flds.jl | 2 +- 3 files changed, 30 insertions(+), 29 deletions(-) diff --git a/mesons.jl b/mesons.jl index 85c6f0c..9288e3f 100644 --- a/mesons.jl +++ b/mesons.jl @@ -1,34 +1,35 @@ using DifferentialEquations -const ħc = 197.327 # ħc in MeVfm +const ħc = 197.327 # MeVfm # Values defined in C. J. Horowitz and J. Piekarewicz, Phys. Rev. Lett. 86, 5647 (2001) # Values taken from Hartree.f (FSUGarnet) and -const m_ρ = 763.0 -const m_ω = 782.5 -const m_s = 496.939473213388 -const m_γ = 0.000001000 # not defined in paper -const g2_s = 110.349189097820 -const g2_v = 187.694676506801 -const g2_ρ = 192.927428365698 -const g2_g = 0.091701236 # not defined in paper -const κ = 3.260178893447 -const λ = -0.003551486718 # LambdaSS -const ζ = 0.023499504053 # LambdaVV -const Λv = 0.043376933644 # LambdaVR +const m_ρ = 763.0 # MeV/c2 +const m_ω = 782.5 # MeV/c2 +const m_s = 496.939473213388 # MeV/c2 +const m_γ = 0.000001000 # MeV/c2 # not defined in paper +const g2_s = 110.349189097820 # units? +const g2_v = 187.694676506801 # units? +const g2_ρ = 192.927428365698 # units? +const g2_g = 0.091701236 # units? # not defined in paper +const κ = 3.260178893447 # units? +const λ = -0.003551486718 # units? # LambdaSS +const ζ = 0.023499504053 # units? # LambdaVV +const Λv = 0.043376933644 # units? # LambdaVR -const r_reg = 1E-9 # regulator for Green's functions in fm +const r_reg = 1E-9 # fm # regulator for Green's functions -"Green's function for Klein-Gordon equation" +"Green's function for Klein-Gordon equation in natural units" greensFunctionKG(m, r, rp) = -(1 / m) * rp / (r + r_reg) * sinh(m * min(r, rp)) * exp(-m * max(r, rp)) -"Green's function for Poisson's equation" +"Green's function for Poisson's equation in natural units" greensFunctionP(r, rp) = -rp^2 / (max(r, rp) + r_reg) "Solve the Klein-Gordon equation (or Poisson's equation when m=0) for a given a source function (as an array) where - the other parameters are the same from mesons!(...)." + m is the mass of the meson in MeV/c2, + r_max is the r-cutoff in fm." function solveKG(m, source, r_max) - greensFunction = m == 0 ? greensFunctionP : (r, rp) -> greensFunctionKG(m, r, rp) + greensFunction = m == 0 ? greensFunctionP : (r, rp) -> greensFunctionKG(m / ħc, r, rp) mesh = range(0, r_max; length=length(source)) greenMat = greensFunction.(transpose(mesh), mesh) return greenMat * source @@ -50,14 +51,14 @@ function solveMesonWfs(ρ_sp, ρ_vp, ρ_sn, ρ_vn, r_max, divs, iterations=3) (src_Φ0, src_W0, src_B0) = (zeros(1 + divs) for _ in 1:3) for _ in 1:iterations - @. src_Φ0 = g2_s * ((κ/2) * Φ0^2 + (λ/6) * Φ0^3 - (ρ_sp + ρ_sn)) / ħc - Φ0 .= solveKG(m_s / sqrt(ħc), src_Φ0, r_max) + @. src_Φ0 = g2_s * ((κ/2) * Φ0^2 + (λ/6) * Φ0^3 - (ρ_sp + ρ_sn)) + Φ0 .= solveKG(m_s, src_Φ0, r_max) - @. src_W0 = g2_v * ((ζ/6) * W0^3 + 2 * Λv * B0^2 * W0 - (ρ_vp + ρ_vn)) / ħc - W0 .= solveKG(m_ω / sqrt(ħc), src_W0, r_max) + @. src_W0 = g2_v * ((ζ/6) * W0^3 + 2 * Λv * B0^2 * W0 - (ρ_vp + ρ_vn)) + W0 .= solveKG(m_ω, src_W0, r_max) - @. src_B0 = g2_ρ * (2 * Λv * W0^2 * B0 - (ρ_vp - ρ_vn) / 2) / ħc - B0 .= solveKG(m_ρ / sqrt(ħc), src_B0, r_max) + @. src_B0 = g2_ρ * (2 * Λv * W0^2 * B0 - (ρ_vp - ρ_vn) / 2) + B0 .= solveKG(m_ρ, src_B0, r_max) end return (Φ0, W0, B0, A0) diff --git a/nucleons.jl b/nucleons.jl index f69766b..c556318 100644 --- a/nucleons.jl +++ b/nucleons.jl @@ -1,10 +1,10 @@ using DifferentialEquations, Roots -const ħc = 197.327 # ħc in MeVfm -const M_n = 939.5654133 # Neutron mass in MeV/c2 -const M_p = 938.2720813 # Proton mass in MeV/c2 +const ħc = 197.327 # MeVfm +const M_n = 939.5654133 # MeV/c2 +const M_p = 938.2720813 # MeV/c2 -const r_reg = 1E-6 # regulator for the centrifugal term in fm +const r_reg = 1E-6 # fm # regulator for the centrifugal term "The spherical Dirac equation that returns du=[dg, df] in-place where u=[g, f] are the reduced radial components evaluated at r, diff --git a/test/Pb208_meson_flds.jl b/test/Pb208_meson_flds.jl index b94e593..b8b75b4 100644 --- a/test/Pb208_meson_flds.jl +++ b/test/Pb208_meson_flds.jl @@ -13,7 +13,7 @@ xs = test_data[:, 1] r_max = maximum(xs) divs = length(xs) - 1 -(Φ0, W0, B0, A0) = solveMesonWfs(ρ_sp, ρ_vp, ρ_sn, ρ_vn, r_max, divs, 3) +(Φ0, W0, B0, A0) = solveMesonWfs(ρ_sp, ρ_vp, ρ_sn, ρ_vn, r_max, divs, 1) plot(xs, hcat(Φ0, W0, B0, A0), layout=4, label=["Φ0" "W0" "B0" "A0"]) xlabel!("r (fm)")