mirror of
https://github.com/tgorordo/WignerSymbols.jl.git
synced 2026-06-05 15:42:15 -07:00
README and benchmarks
This commit is contained in:
Jutho
parent
1a974193eb
commit
c47faaf768
3 changed files with 160 additions and 8 deletions
10
LICENSE.md
10
LICENSE.md
|
|
@ -1,17 +1,17 @@
|
|||
The WignerSymbols.jl package is licensed under the MIT "Expat" License:
|
||||
|
||||
> Copyright (c) 2017: Jutho.
|
||||
>
|
||||
> Copyright (c) 2017: Jutho Haegeman.
|
||||
>
|
||||
> Permission is hereby granted, free of charge, to any person obtaining a copy
|
||||
> of this software and associated documentation files (the "Software"), to deal
|
||||
> in the Software without restriction, including without limitation the rights
|
||||
> to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
|
||||
> copies of the Software, and to permit persons to whom the Software is
|
||||
> furnished to do so, subject to the following conditions:
|
||||
>
|
||||
>
|
||||
> The above copyright notice and this permission notice shall be included in all
|
||||
> copies or substantial portions of the Software.
|
||||
>
|
||||
>
|
||||
> THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
|
||||
> IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
|
||||
> FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
|
||||
|
|
@ -19,4 +19,4 @@ The WignerSymbols.jl package is licensed under the MIT "Expat" License:
|
|||
> LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
|
||||
> OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
|
||||
> SOFTWARE.
|
||||
>
|
||||
>
|
||||
|
|
|
|||
52
README.md
52
README.md
|
|
@ -1,7 +1,53 @@
|
|||
# WignerSymbols
|
||||
|
||||
[](https://travis-ci.org/jutho/WignerSymbols.jl)
|
||||
[](https://travis-ci.org/jutho/WignerSymbols.jl)
|
||||
[](LICENSE.md)
|
||||
[](https://coveralls.io/github/Jutho/WignerSymbols.jl?branch=master)
|
||||
[](http://codecov.io/github/Jutho/WignerSymbols.jl?branch=master)
|
||||
|
||||
[](https://coveralls.io/github/jutho/WignerSymbols.jl?branch=master)
|
||||
Compute Wigner's 3j and 6j symbols, and related quantities such as Clebsch-Gordan coefficients and Racah's symbols.
|
||||
|
||||
[](http://codecov.io/github/jutho/WignerSymbols.jl?branch=master)
|
||||
## Requirements
|
||||
This requires a recent master edition of Julia (i.e. v0.7.0-DEV), because it depends on some changes in `Base.GMP`. In particular, it uses the mutating functions for reducing allocation overhead while working with `BigInts` (namely [JuliaLang/julia#21654](https://github.com/JuliaLang/julia/pull/21654)). It also depends on `Primes.jl` for generating prime numbers.
|
||||
|
||||
## Installation
|
||||
Until it is registered, install via `Pkg.clone("https://github.com/Jutho/WignerSymbols.jl.git")`.
|
||||
|
||||
## Available functions
|
||||
While the following function signatures are probably self-explanatory, you can query help for them in the Julia REPL to get further details.
|
||||
* `wigner3j(T::Type{<:AbstractFloat} = Float64, j₁, j₂, j₃, m₁, m₂, m₃ = -m₂-m₁) -> ::T`
|
||||
* `wigner6j(T::Type{<:AbstractFloat} = Float64, j₁, j₂, j₃, j₄, j₅, j₆) -> ::T`
|
||||
* `clebschgordan(T::Type{<:AbstractFloat} = Float64, j₁, j₂, j₃, m₁, m₂, m₃ = m₁+m₂) -> ::T`
|
||||
* `racahV(T::Type{<:AbstractFloat} = Float64, j₁, j₂, j₃, m₁, m₂, m₃ = -m₁-m₂) -> ::T`
|
||||
* `racahW(T::Type{<:AbstractFloat} = Float64, j₁, j₂, J, j₃, J₁₂, J₂₃) -> ::T`
|
||||
* `δ(j₁, j₂, j₃) -> ::Bool`
|
||||
* `Δ(T::Type{<:AbstractFloat} = Float64, j₁, j₂, j₃) -> ::T`
|
||||
|
||||
## Implementation
|
||||
Largely based on reading the paper (but not the code):
|
||||
|
||||
[1] [H. T. Johansson and C. Forssén, SIAM Journal on Scientific Compututing 38 (2016) 376-384](https://doi.org/10.1137/15M1021908) ([arXiv:1504.08329](https://arxiv.org/abs/1504.08329))
|
||||
|
||||
with some additional modifications to further improve efficiency for large `j` (angular momenta quantum numbers).
|
||||
|
||||
In particular, 3j and 6j symbols are computed exactly, in the format `√(r) * s` where `r` and `s` are exactly computed as `Rational{BigInt}`,
|
||||
using an intermediate representation based on prime number factorization. As a consequence thereof, all of the above functions can be called
|
||||
requesting `BigFloat` precision for the result. There is currently no convenient syntax for obtaining `r` and `s` directly (see TODO).
|
||||
|
||||
Most intermediate calculations (prime factorizations of numbers and their factorials, conversion between prime powers and `BigInt`s) are cached to improve the efficiency, but this can result in large use of memory when querying Wigner symbols for large values of `j`.
|
||||
|
||||
Also uses ideas from
|
||||
|
||||
[2] [J. Rasch and A. C. H. Yu, SIAM Journal on Scientific Compututing 25 (2003), 1416–1428](https://doi.org/10.1137/S1064827503422932)
|
||||
|
||||
for caching the computed 3j and 6j symbols.
|
||||
|
||||
# Benchmark:
|
||||
|
||||
## Todo:
|
||||
* Wigner 9-j symbols, as explained in [1] and based on
|
||||
|
||||
[3] [L. Wei, New formula for 9-j symbols and their direct calculation, Computers in Physics, 12 (1998), 632–634.](ftp://gravity.physics.uwa.edu.au/pub/Clebsch-Gordan/Papers/Wei9j.pdf)
|
||||
|
||||
* Convenient syntax to get the exact results in the `√(r) * s` format, but in such a way that it can be used by
|
||||
the Julia type system and can be converted afterwards.
|
||||
|
|
|
|||
106
benchmark/benchmarks.jl
Normal file
106
benchmark/benchmarks.jl
Normal file
|
|
@ -0,0 +1,106 @@
|
|||
# some random benchmarks; for j=50000 you need more than 8GB of memory
|
||||
using WignerSymbols
|
||||
|
||||
@time wigner3j(BigFloat, 15, 30, 40, 2, 2, -4)
|
||||
@time wigner3j(BigFloat, 200, 200, 200, -10, 60, -50)
|
||||
@time wigner3j(BigFloat, 50000, 50000, 50000, 1000, -6000, 5000)
|
||||
|
||||
j = 8
|
||||
@time wigner6j(BigFloat, j, j, j, j, j, j)
|
||||
j = 200
|
||||
@time wigner6j(BigFloat, j, j, j, j, j, j)
|
||||
j = 600
|
||||
@time wigner6j(BigFloat, j, j, j, j, j, j)
|
||||
j = 10000
|
||||
@time wigner6j(BigFloat, j, j, j, j, j, j)
|
||||
j = 50000
|
||||
@time wigner6j(BigFloat, j, j, j, j, j, j)
|
||||
|
||||
function compute3jmax(jmax)
|
||||
for j₁ = 0:1//2:jmax
|
||||
for j₂ = 0:1//2:j₁
|
||||
for j₃ = abs(j₁-j₂):j₂
|
||||
if isinteger(j₁)
|
||||
m₁ = 0
|
||||
M₂ = min(j₂, j₃-m₁)
|
||||
for m₂ = M₂:(-1):0
|
||||
wigner3j(j₁,j₂,j₃,0, m₂)
|
||||
end
|
||||
for m₁ = 1:j₁
|
||||
M₂ = min(j₂, j₃-m₁)
|
||||
for m₂ = -M₂:M₂
|
||||
wigner3j(j₁, j₂, j₃, m₁, m₂)
|
||||
end
|
||||
end
|
||||
else
|
||||
for m₁ = 1//2:j₁
|
||||
M₂ = min(j₂, j₃-m₁)
|
||||
for m₂ = -M₂:M₂
|
||||
wigner3j(j₁, j₂, j₃, m₁, m₂)
|
||||
end
|
||||
end
|
||||
end
|
||||
end
|
||||
end
|
||||
end
|
||||
return nothing
|
||||
end
|
||||
|
||||
@time compute3jmax(10)
|
||||
|
||||
function computeinteger3jmax(jmax)
|
||||
for j₁ = 0:jmax
|
||||
for j₂ = 0:j₁
|
||||
for j₃ = abs(j₁-j₂):j₂
|
||||
m₁ = 0
|
||||
M₂ = min(j₂, j₃-m₁)
|
||||
for m₂ = M₂:(-1):0
|
||||
wigner3j(j₁,j₂,j₃,0, m₂)
|
||||
end
|
||||
for m₁ = 1:j₁
|
||||
M₂ = min(j₂, j₃-m₁)
|
||||
for m₂ = -M₂:M₂
|
||||
wigner3j(j₁, j₂, j₃, m₁, m₂)
|
||||
end
|
||||
end
|
||||
end
|
||||
end
|
||||
end
|
||||
return nothing
|
||||
end
|
||||
|
||||
function compute6jmax(jmax)
|
||||
for j₁ = 0:1//2:jmax
|
||||
for j₂ = 0:1//2:j₁
|
||||
for j₃ = abs(j₁-j₂):j₂
|
||||
for j₄ = 0:1//2:jmax
|
||||
for j₅ = abs(j₃-j₄):min(jmax,j₃+j₄)
|
||||
for j₆ = min(abs(j₂-j₄),abs(j₁-j₅)):min(j₂+j₄,j₁+j₅,jmax)
|
||||
wigner6j(j₁, j₂, j₃, j₄, j₅, j₆)
|
||||
end
|
||||
end
|
||||
end
|
||||
end
|
||||
end
|
||||
end
|
||||
return nothing
|
||||
end
|
||||
@time compute6jmax(10)
|
||||
|
||||
|
||||
function computeinteger6jmax(jmax)
|
||||
for j₁ = 0:jmax
|
||||
for j₂ = 0:j₁
|
||||
for j₃ = abs(j₁-j₂):j₂
|
||||
for j₄ = 0:jmax
|
||||
for j₅ = abs(j₃-j₄):min(jmax,j₃+j₄)
|
||||
for j₆ = min(abs(j₂-j₄),abs(j₁-j₅)):min(j₂+j₄,j₁+j₅,jmax)
|
||||
wigner6j(j₁, j₂, j₃, j₄, j₅, j₆)
|
||||
end
|
||||
end
|
||||
end
|
||||
end
|
||||
end
|
||||
end
|
||||
return nothing
|
||||
end
|
||||
Loading…
Add table
Add a link
Reference in a new issue