Skip to content

Commit

Permalink
Add eigenvalue test for rigid rotor and delta potential
Browse files Browse the repository at this point in the history
  • Loading branch information
@ajarifi
ajarifi committed May 24, 2024
1 parent 2844039 commit 914329a
Show file tree
Hide file tree
Showing 4 changed files with 311 additions and 135 deletions.
41 changes: 41 additions & 0 deletions test/DeltaPotential.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -31,3 +31,44 @@ end
println("""```
""")

# <ψₙ|H|ψₙ> = ∫ψₙ*Tψₙdx = Eₙ

println(raw"""
#### Eigen Values
```math
\begin{aligned}
E_n
&= \int_{-\infty}^{\infty} \psi^\ast(x) \hat{H} \psi(x) ~\mathrm{d}x \\
&= \int_{-\infty}^{\infty} \psi^\ast(x) \left[ \hat{V} + \hat{T} \right] \psi(x) ~\mathrm{d}x \\
&= \int_{-\infty}^{\infty} \psi^\ast(x) \left[ -\alpha\delta(x) - \frac{\hbar^2}{2m} \frac{\mathrm{d}^{2}}{\mathrm{d} x^{2}} \right] \psi(x) ~\mathrm{d}x \\
&= \int_{-\infty}^{\infty} \psi^\ast(x) \left[ -\alpha\delta(x) - \frac{\hbar^2}{2m} (\kappa^2 -2\kappa \delta(x))\right]\psi(x) ~\mathrm{d}x \\
\end{aligned}
```
where the $\kappa=m\alpha/\hbar^2$ and the integration with the delta function yeild a term proportional to the wave function at the origin $|\psi(0)|^2$.
```""")

@testset "<ψₙ|H|ψₙ> = ∫ψₙ*Hψₙdx = Eₙ" begin
function ψHψ(DP)
κ = DP.m * DP.α/ DP.ℏ^2
T = -DP.ℏ^2 / (2 * DP.m) *( κ^2 - 2κ * ψ(DP,0)^2 )
V = -DP.α * ψ(DP, 0)^2
return T + V
end

println(" α | m | ℏ | analytical | numerical ")
println("--- | --- | --- | ----------------- | ----------------- ")
for α in [0.1, 0.5, 2.0]
for m in [0.1, 0.5, 2.0]
for ℏ in [0.1, 0.5, 2.0]
DP = DeltaPotential(α=α, m=m, ℏ=ℏ)
analytical = E(DP)
numerical = ψHψ(DP)
acceptance = iszero(analytical) ? isapprox(analytical, numerical, atol=1e-4) : isapprox(analytical, numerical, rtol=1e-4)
@printf("%.1f | %.1f | %.1f | %17.12f | %17.12f %s\n", α, m, ℏ,analytical, numerical, acceptance ? "✔" : "✗")
end
end
end
end

println("""```""")
51 changes: 51 additions & 0 deletions test/RigidRotor.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -123,5 +123,56 @@ Y_{lm}(\theta,\varphi)^* Y_{l'm'}(\theta,\varphi) \sin(\theta)
end
end

println("""```
""")


# <ψₙ|H|ψₙ> = ∫ψₙ*Tψₙdx = Eₙ

println(raw"""
#### Eigen Values
```math
\begin{aligned}
E_n
&= \int_{0}^{2\pi}~\mathrm{d}\varphi \int_{0}^{\pi}~\mathrm{d}\theta ~ \psi^\ast_{lm}(\theta,\varphi) \hat{H} \psi_{lm}(\theta,\varphi) \\
&= \int_{0}^{2\pi}~\mathrm{d}\varphi \int_{0}^{\pi}~\mathrm{d}\theta ~ \psi^\ast_{lm}(\theta,\varphi) \left[ \hat{V} + \hat{T} \right] \psi_{lm}(\theta,\varphi) \\
&= \int_{0}^{2\pi}~\mathrm{d}\varphi \int_{0}^{\pi}~\mathrm{d}\theta ~ \psi^\ast_{lm}(\theta,\varphi) \left[ 0 + \frac{l(l+1)}{2I} \right] \psi_{lm}(\theta,\varphi) \\
&= \int_{0}^{2\pi}~\mathrm{d}\varphi \int_{0}^{\pi}~\mathrm{d}\theta ~ \frac{l(l+1)}{2I} |\psi_{lm}(\theta,\varphi)|^2
\end{aligned}
```
where $l$ is angular momentum quantum number and $I$ is the moment of intertia.
```""")

@testset "<ψₙ|H|ψₙ> = ∫ψₙ*Hψₙdx = Eₙ" begin
function ψHψ(RR;l)
μ = (1/RR.m₁ + 1/RR.m₂)^(-1)
I = μ * RR.R^2
YY = real(quadgk(φ -> quadgk(θ -> conj(Y(RR,θ,φ,l=l,m=0)) * Y(RR,θ,φ,l=l,m=0) * sin(θ), 0, π, maxevals=10)[1], 0, 2π, maxevals=10)[1])
T = l*(l+1)/2/I * YY
return T
end

ℏ=1
println(" m₁ | m₂ | R | l | analytical | numerical ")
println("--- | --- | --- | --- | ----------------- | ----------------- ")
for m₁ in [0.5, 2.0]
for m₂ in [0.5, 2.0]
for R in [0.5, 2.0]
for l in [1, 2, 3]
RR = RigidRotor(m₁=m₁, m₂=m₂, R=R, ℏ=ℏ)
analytical = E(RR,l=l)
numerical = ψHψ(RR,l=l)
acceptance = iszero(analytical) ? isapprox(analytical, numerical, atol=1e-4) : isapprox(analytical, numerical, rtol=1e-4)
@printf("%.1f | %.1f | %.1f | %.1f | %17.12f | %17.12f %s\n", m₁, m₂, R, l,analytical, numerical, acceptance ? "✔" : "✗")
end
end
end
end
end




println("""```
""")
43 changes: 43 additions & 0 deletions test/result/DeltaPotential.log
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,4 +36,47 @@
7.0 | 7.0 | 7.0 | 1.000000000000 | 1.000000000000 ✔
```

#### Eigen Values

```math
\begin{aligned}
E_n
&= \int_{-\infty}^{\infty} \psi^\ast(x) \hat{H} \psi(x) ~\mathrm{d}x \\
&= \int_{-\infty}^{\infty} \psi^\ast(x) \left[ \hat{V} + \hat{T} \right] \psi(x) ~\mathrm{d}x \\
&= \int_{-\infty}^{\infty} \psi^\ast(x) \left[ -\alpha\delta(x) - \frac{\hbar^2}{2m} \frac{\mathrm{d}^{2}}{\mathrm{d} x^{2}} \right] \psi(x) ~\mathrm{d}x \\
&= \int_{-\infty}^{\infty} \psi^\ast(x) \left[ -\alpha\delta(x) - \frac{\hbar^2}{2m} (\kappa^2 -2\kappa \delta(x))\right]\psi(x) ~\mathrm{d}x \\
\end{aligned}
```
where the $\kappa=m\alpha/\hbar^2$ and the integration with the delta function yeild a term proportional to the wave function at the origin $|\psi(0)|^2$.
```
α | m | ℏ | analytical | numerical
--- | --- | --- | ----------------- | -----------------
0.1 | 0.1 | 0.1 | -0.050000000000 | -0.050000000000 ✔
0.1 | 0.1 | 0.5 | -0.002000000000 | -0.002000000000 ✔
0.1 | 0.1 | 2.0 | -0.000125000000 | -0.000125000000 ✔
0.1 | 0.5 | 0.1 | -0.250000000000 | -0.250000000000 ✔
0.1 | 0.5 | 0.5 | -0.010000000000 | -0.010000000000 ✔
0.1 | 0.5 | 2.0 | -0.000625000000 | -0.000625000000 ✔
0.1 | 2.0 | 0.1 | -1.000000000000 | -1.000000000000 ✔
0.1 | 2.0 | 0.5 | -0.040000000000 | -0.040000000000 ✔
0.1 | 2.0 | 2.0 | -0.002500000000 | -0.002500000000 ✔
0.5 | 0.1 | 0.1 | -1.250000000000 | -1.250000000000 ✔
0.5 | 0.1 | 0.5 | -0.050000000000 | -0.050000000000 ✔
0.5 | 0.1 | 2.0 | -0.003125000000 | -0.003125000000 ✔
0.5 | 0.5 | 0.1 | -6.250000000000 | -6.250000000000 ✔
0.5 | 0.5 | 0.5 | -0.250000000000 | -0.250000000000 ✔
0.5 | 0.5 | 2.0 | -0.015625000000 | -0.015625000000 ✔
0.5 | 2.0 | 0.1 | -25.000000000000 | -25.000000000000 ✔
0.5 | 2.0 | 0.5 | -1.000000000000 | -1.000000000000 ✔
0.5 | 2.0 | 2.0 | -0.062500000000 | -0.062500000000 ✔
2.0 | 0.1 | 0.1 | -20.000000000000 | -20.000000000000 ✔
2.0 | 0.1 | 0.5 | -0.800000000000 | -0.800000000000 ✔
2.0 | 0.1 | 2.0 | -0.050000000000 | -0.050000000000 ✔
2.0 | 0.5 | 0.1 | -100.000000000000 | -100.000000000000 ✔
2.0 | 0.5 | 0.5 | -4.000000000000 | -4.000000000000 ✔
2.0 | 0.5 | 2.0 | -0.250000000000 | -0.250000000000 ✔
2.0 | 2.0 | 0.1 | -400.000000000000 | -400.000000000000 ✔
2.0 | 2.0 | 0.5 | -16.000000000000 | -16.000000000000 ✔
2.0 | 2.0 | 2.0 | -1.000000000000 | -1.000000000000 ✔
```

Loading

0 comments on commit 914329a

Please sign in to comment.