Simpson's rule

This commit is contained in:
Nuwan Yapa 2025-08-04 14:30:36 -04:00
parent 0af47ef161
commit 5fc391ee74
1 changed files with 4 additions and 1 deletions

View File

@ -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)))