From 5fc391ee7431e684703c4b43241b9d97713f7635 Mon Sep 17 00:00:00 2001 From: Nuwan Yapa Date: Mon, 4 Aug 2025 14:30:36 -0400 Subject: [PATCH] Simpson's rule --- mesons.jl | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/mesons.jl b/mesons.jl index 1481acd..b669710 100644 --- a/mesons.jl +++ b/mesons.jl @@ -25,7 +25,10 @@ greensFunctionP(r, rp) = -1 / (max(r, rp) + r_reg) "Solve the Klein-Gordon equation (or Poisson's equation when m=0) and return an array in MeV for a source array given in fm⁻³ where m is the mass of the meson in MeV/c2." function solveKG(m, source, s::system) - int_measure = ħc .* Δr(s) .* rs(s) .^ 2 + + @assert s.divs % 2 == 0 "Number of mesh divisions must be even for Simpson's rule" + simpsons_weights = (Δr(s)/3) .* [1; repeat([2,4], s.divs ÷ 2)[2:end]; 1] + int_measure = ħc .* simpsons_weights .* rs(s) .^ 2 greensFunction = m == 0 ? greensFunctionP : (r, rp) -> greensFunctionKG(m / ħc, r, rp) greenMat = greensFunction.(rs(s), transpose(rs(s)))