Elegant method to attach two partitions

This commit is contained in:
Nuwan Yapa 2025-04-07 19:24:42 -04:00
parent 9f8cea12a3
commit a143edb496
1 changed files with 12 additions and 9 deletions

View File

@ -1,4 +1,4 @@
using DifferentialEquations, Interpolations
using LinearAlgebra, DifferentialEquations, Interpolations
include("bisection.jl")
include("common.jl")
include("system.jl")
@ -73,17 +73,20 @@ function solveNucleonWf(κ, p::Bool, E, s::system; shooting=true, normalize=true
wf = hcat(sol.u...)
if shooting # join two segments
rescale_factor_g = wf[1, end] / wf_right[1, 1]
rescale_factor_f = wf[2, end] / wf_right[2, 1]
@assert isfinite(rescale_factor_g) && isfinite(rescale_factor_f) "Cannot rescale the right partition"
isapprox(rescale_factor_g, rescale_factor_f; rtol=0.03) || @warn "Discontinuity between two partitions"
wf_right_rescaled = wf_right .* ((rescale_factor_g + rescale_factor_f) / 2)
wf = hcat(wf[:, 1:(end - 1)], wf_right_rescaled)
u1 = wf[:, end]
u2 = wf_right[:, 1]
if norm(u2) < 1e-10
@warn "Right partition too small to rescale, setting to zero"
wf_right .= 0.0
else
proj = only(u1' * u2) / norm(u2)^2
wf_right .*= proj
end
wf = hcat(wf[:, 1:(end - 1)], wf_right)
end
if normalize
norm = sum(wf .* wf) * Δr(s) # integration by Reimann sum
wf = wf ./ sqrt(norm)
wf ./= norm(wf) * Δr(s) # integration by Reimann sum
end
return wf