Comparing discretization techniques
In Periodic problems and plane-wave discretizations we saw how simple 1D problems can be modelled using plane-wave basis sets. This example invites you to work out some details on these aspects yourself using a number of exercises. The solutions are given at the bottom of the page.
For this example we consider the discretization of
\[ H = - \frac12 Δ + V(x) \quad \text{with $V(x) = \cos(x)$}\]
on $[0, 2π]$ with periodic boundary conditions. The $\cos(x)$ takes the role of a lattice-periodic potential. We will be interested in the smallest eigenvalues of this discretized Hamiltonian. Of note the boundary condition matters: The spectrum we will get is different from e.g. the spectrum of $H$ on $\mathbb{R}$.
Finite differences
We approximate functions $ψ$ on $[0, 2\pi]$ by their values at grid points $x_k = 2\pi \frac{k}{N}$, $k=1, \dots, N$. The boundary conditions are imposed by $ψ(x_0) = ψ(x_N), ψ(x_{N+1}) = ψ(x_1)$. We then have
\[ \big(Hψ\big)(x_k) \approx \frac 1 2 \frac{-ψ_{k-1} + 2 ψ_k - ψ_{k+1}}{δx^2} + V(x_k) ψ(x_k)\]
with $δx = \frac{2π}{N}$.
This can be put in matrix form in the following way:
# Finite differences Hamiltonian -1/2 Delta + V on [0, 2pi] with periodic bc.
# Pass it a function V.
using LinearAlgebra
function build_finite_differences_matrix(Vfunction, N::Integer)
δx = 2π/N
# Finite-difference approximation to -½Δ
T = 1/(2δx^2) * Tridiagonal(-ones(N-1), 2ones(N), -ones(N-1))
# The type Tridiagonal is efficient, but to establish the periodic boundary conditions
# we need to add elements not on the three diagonals, so convert to dense matrix
T = Matrix(T)
T[1, N] = T[N, 1] = -1 / (2δx^2)
# Finite-difference approximation to potential: We collect all coordinates ...
x_coords = [k * δx for k=1:N]
V = Diagonal(Vfunction.(x_coords)) # ... and evaluate V on each of the x_coords
T + V
end;Show that the finite-difference approximation of -½Δ is indeed an approximation of the second derivative. Obtain an estimate of the first eigenvalue of $H$. Hint: Take a look at the eigen function from LinearAlgebra.
Plane waves method
In this method, we expand states on the basis
\[ e_G(x) = \frac{1}{\sqrt{2\pi}} e^{iGx} \qquad \text{for $G=-N,\dots,N$}.\]
Show that
\[ \langle e_G, e_{G'}\rangle = ∫_0^{2π} e_G^\ast(x) e_{G'}(x) d x = δ_{G, G'}\]
and (assuming $V(x) = \cos(x)$)
\[ \langle e_G, H e_{G'}\rangle = \frac 1 2 \left(|G|^2 \delta_{G,G'} + \delta_{G, G'+1} + \delta_{G, G'-1}\right).\]
What happens for a more general $V(x)$?
Code this and check the first eigenvalue agrees with the finite-difference case. Compare accuracies at various basis set sizes $N$.
Using DFTK
We now use DFTK to do the same plane-wave discretization in this 1D system. To deal with a 1D case we use a 3D lattice with two lattice vectors set to zero.
using DFTK
a = 2π
lattice = a .* [[1 0 0.]; [0 0 0]; [0 0 0]];Define Hamiltonian: Kinetic + Potential
terms = [Kinetic(),
ExternalFromReal(r -> cos(r[1]))] # r is a vector of size 3
model = Model(lattice; n_electrons=1, terms, spin_polarization=:spinless); # One spinless electronEcut defines the number of plane waves by selecting all those $G$, which satisfy the relationship $½ |G|^2 ≤ \text{Ecut}$.
Ecut = 500
basis = PlaneWaveBasis(model; Ecut, kgrid=(1, 1, 1))PlaneWaveBasis discretization:
architecture : DFTK.CPU()
num. mpi processes : 1
num. julia threads : 1
num. DFTK threads : 1
num. blas threads : 2
num. fft threads : 1
Ecut : 500.0 Ha
fft_size : (135, 1, 1), 135 total points
kgrid : MonkhorstPack([1, 1, 1])
num. red. kpoints : 1
num. irred. kpoints : 1
Estimated memory usage (per MPI process):
single ψ : 3.94KiB
single ρ : 1.05KiB
peak memory (SCF) : 44.2KiB
Discretized Model(custom, 1D):
lattice (in Bohr) : [6.28319 , 0 , 0 ]
[0 , 0 , 0 ]
[0 , 0 , 0 ]
unit cell volume : 6.2832 Bohr
num. electrons : 1
spin polarization : spinless
temperature : 0 Ha
terms : Kinetic()
ExternalFromReal(Main.var"##354".var"#3#4"())We now seek the ground state using the self-consistent field algorithm.
scfres = self_consistent_field(basis; tol=1e-4)
scfres.energiesEnergy breakdown (in Ha):
Kinetic 0.2090845
ExternalFromReal -0.7441493
total -0.535064852288On this simple linear (non-interacting) model, the SCF converges in one step. The ground state energy of is simply the lowest eigenvalue; it should match the smallest eigenvalue of $H$ computed above.
Plotting
We can also get the first eigenvector (in the plane wave basis) and plot it
using Plots
ψ_fourier = scfres.ψ[1][:, 1]; # first k-point, all G components, first eigenvector
# Transform the wave function to real space
ψ = ifft(basis, basis.kpoints[1], ψ_fourier)[:, 1, 1]
# Eigenvectors are only defined up to a phase. We fix it by imposing that psi(0) is real
ψ /= (ψ[1] / abs(ψ[1]))
plot(real(ψ); label="")Again this should match with the result above.
Look at the Fourier coefficients of ψ_fourier and compare with the result above.
The DFTK Hamiltonian
We can ask DFTK for the Hamiltonian
E, ham = energy_hamiltonian(basis, scfres.ψ, scfres.occupation; ρ=scfres.ρ)
H = ham.blocks[1]
typeof(H)DFTK.DftHamiltonianBlockThis is an opaque data structure, which encodes the Hamiltonian. What can we do with it?
using InteractiveUtils
methodswith(typeof(H), supertypes=true)9-element Vector{Method}:- Array(block::HamiltonianBlock) in DFTK at /home/runner/work/DFTK.jl/DFTK.jl/src/terms/Hamiltonian.jl:74
- *(H::HamiltonianBlock, ψ) in DFTK at /home/runner/work/DFTK.jl/DFTK.jl/src/terms/Hamiltonian.jl:63
- (Matrix)(block::HamiltonianBlock) in DFTK at /home/runner/work/DFTK.jl/DFTK.jl/src/terms/Hamiltonian.jl:75
- eltype(block::HamiltonianBlock) in DFTK at /home/runner/work/DFTK.jl/DFTK.jl/src/terms/Hamiltonian.jl:64
- size(block::HamiltonianBlock, i::Integer) in DFTK at /home/runner/work/DFTK.jl/DFTK.jl/src/terms/Hamiltonian.jl:65
- size(block::HamiltonianBlock) in DFTK at /home/runner/work/DFTK.jl/DFTK.jl/src/terms/Hamiltonian.jl:66
- mul!(Hψ::AbstractArray, H::DFTK.DftHamiltonianBlock, ψ::AbstractArray) in DFTK at /home/runner/work/DFTK.jl/DFTK.jl/src/terms/Hamiltonian.jl:137
- PreconditionerNone(::HamiltonianBlock) in DFTK at /home/runner/work/DFTK.jl/DFTK.jl/src/eigen/preconditioners.jl:21
- PreconditionerTPA(ham::HamiltonianBlock; kwargs...) in DFTK at /home/runner/work/DFTK.jl/DFTK.jl/src/eigen/preconditioners.jl:45
This defines a number of methods. For instance, it can be used as a linear operator:
H * DFTK.random_orbitals(basis, basis.kpoints[1], 1)63×1 Matrix{ComplexF64}:
0.005610556324933323 - 0.13713478136307972im
-0.00022716191215306145 - 0.18685765653846337im
0.02742925931596646 - 0.37765063414509314im
0.49287383055752115 - 0.43248706783556357im
-0.1549815782869785 + 0.3456945501759658im
-1.6385020340833576 + 0.2430014591050156im
0.44075367893039186 - 0.6359376404522141im
-3.480818496069009 - 0.42256996943683656im
-0.8925906117043708 - 2.925453568243214im
1.8715513189089634 - 0.35407597291603227im
5.440549149805195 + 4.2166662193827795im
2.369419932415937 - 3.80044014935152im
-6.510457361910034 - 4.102724632331347im
-6.921160769077959 + 9.86594145604858im
7.821067804191677 - 7.137346864415825im
2.6254067148642397 - 20.186504686041825im
10.474699105176393 - 28.21032944588713im
1.6739542891508012 + 12.169090137435264im
13.538079441974116 - 9.797080486439775im
2.087598035735784 + 18.849059316276882im
20.97994941170068 + 2.413975261256023im
7.160678037981892 - 18.950009864279394im
-9.426185400902021 + 10.614623946016586im
10.925528668211625 + 15.370727319596917im
-11.517982279515698 - 0.28082816019553525im
37.09851609603968 - 20.453678637255088im
6.954192984569263 - 52.97055471235189im
27.393233652984453 + 9.192632840116728im
-44.09206324943364 - 17.597111520497844im
36.44856381094992 + 27.334075836115733im
43.98609134601144 + 39.4878992558371im
-74.71261782488551 - 6.139861744420296im
18.372330928118842 + 15.615557891146075im
-54.17640120963555 - 3.0626205757500236im
0.9457938614766469 + 35.085767479459im
7.821914823710978 + 91.41892705271793im
30.692223117456578 + 12.645412657365284im
-1.7323537418836064 + 39.02829737051083im
1.507689399463484 - 12.413824274776843im
32.28169685892509 + 22.923941719115284im
-20.23028852950254 - 32.47187873145238im
-25.049833676850177 + 16.531036753854092im
10.66759155604768 + 19.718607557887527im
27.31887323583357 + 9.077035065763093im
-25.85936593056658 + 6.389334695986996im
-1.7037268617743884 + 3.04217609880579im
6.370775749711867 - 4.900310651081236im
18.473877390125132 - 6.436050685871519im
8.835014342004996 + 11.642370899390574im
-10.840836610458226 + 9.049087301459641im
4.848221305090208 - 5.46522746393955im
0.8235049885924072 - 15.496252421511118im
-1.5139255804535272 + 2.232129387986574im
6.732742479652089 - 4.150638041585512im
-6.553416326883865 - 5.916610504691713im
0.06243264926644415 + 1.4162361975408109im
0.24163340223527544 + 2.595314985382359im
-3.0716559016164804 + 0.8592927769749494im
0.6745703323656937 + 1.9932040556519581im
-0.4817918335883232 - 0.8860396815845123im
0.3363707909976851 - 0.05078963372563038im
0.05505464008491497 - 0.0436739721826989im
-0.02024227354009552 - 0.016433417782910287imWe can also get its full matrix representation:
Array(H)63×63 Matrix{ComplexF64}:
3.21245e-19+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im
0.5+9.98774e-18im 0.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im
8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 2.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im
-2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 4.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im
-1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 8.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im
1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 12.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im
1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 18.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im
-2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 24.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im
-8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 32.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im
4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 40.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im
-1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 50.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im
2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 60.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im
-2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 72.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im
-3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 84.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im
2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 98.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im
7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 112.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im
-2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 128.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im
0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 144.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im
6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 162.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im
-8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 180.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im
3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 200.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -4.10876e-18-1.59967e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im
-3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 220.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 8.50869e-18-3.40232e-18im -4.10876e-18-1.59967e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im
-9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 242.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 0.0+2.88705e-17im 8.50869e-18-3.40232e-18im -4.10876e-18-1.59967e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im
1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 264.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im -1.21529e-17+1.03447e-17im 0.0+2.88705e-17im 8.50869e-18-3.40232e-18im -4.10876e-18-1.59967e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im
-5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 288.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im -1.96473e-17+2.10786e-17im -1.21529e-17+1.03447e-17im 0.0+2.88705e-17im 8.50869e-18-3.40232e-18im -4.10876e-18-1.59967e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im
-2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 312.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -6.4196e-18-3.9968e-18im -1.96473e-17+2.10786e-17im -1.21529e-17+1.03447e-17im 0.0+2.88705e-17im 8.50869e-18-3.40232e-18im -4.10876e-18-1.59967e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im
9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 338.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 7.85256e-20-1.51998e-17im -6.4196e-18-3.9968e-18im -1.96473e-17+2.10786e-17im -1.21529e-17+1.03447e-17im 0.0+2.88705e-17im 8.50869e-18-3.40232e-18im -4.10876e-18-1.59967e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im
6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 364.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 7.00086e-18-3.1146e-18im 7.85256e-20-1.51998e-17im -6.4196e-18-3.9968e-18im -1.96473e-17+2.10786e-17im -1.21529e-17+1.03447e-17im 0.0+2.88705e-17im 8.50869e-18-3.40232e-18im -4.10876e-18-1.59967e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im
-1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 392.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -3.51205e-18+2.34442e-18im 7.00086e-18-3.1146e-18im 7.85256e-20-1.51998e-17im -6.4196e-18-3.9968e-18im -1.96473e-17+2.10786e-17im -1.21529e-17+1.03447e-17im 0.0+2.88705e-17im 8.50869e-18-3.40232e-18im -4.10876e-18-1.59967e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im
3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 420.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im 8.55835e-18+5.02783e-18im -3.51205e-18+2.34442e-18im 7.00086e-18-3.1146e-18im 7.85256e-20-1.51998e-17im -6.4196e-18-3.9968e-18im -1.96473e-17+2.10786e-17im -1.21529e-17+1.03447e-17im 0.0+2.88705e-17im 8.50869e-18-3.40232e-18im -4.10876e-18-1.59967e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im
1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 450.0+0.0im 0.5+1.16569e-17im -2.68068e-18+8.11096e-18im 8.55835e-18+5.02783e-18im -3.51205e-18+2.34442e-18im 7.00086e-18-3.1146e-18im 7.85256e-20-1.51998e-17im -6.4196e-18-3.9968e-18im -1.96473e-17+2.10786e-17im -1.21529e-17+1.03447e-17im 0.0+2.88705e-17im 8.50869e-18-3.40232e-18im -4.10876e-18-1.59967e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im 3.0551e-18-4.69829e-18im
3.0551e-18-4.69829e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 480.5+0.0im -3.18704e-18-2.38796e-17im -2.68068e-18+8.11096e-18im 8.55835e-18+5.02783e-18im -3.51205e-18+2.34442e-18im 7.00086e-18-3.1146e-18im 7.85256e-20-1.51998e-17im -6.4196e-18-3.9968e-18im -1.96473e-17+2.10786e-17im -1.21529e-17+1.03447e-17im 0.0+2.88705e-17im 8.50869e-18-3.40232e-18im -4.10876e-18-1.59967e-17im -4.53964e-18-4.47079e-18im -2.33026e-17+1.55125e-17im 1.48727e-18-1.13135e-17im -9.71381e-19+1.04378e-17im -1.35822e-18+7.74245e-18im -3.63007e-18-1.40327e-17im 2.63164e-17+2.78475e-17im 2.35457e-17+1.83632e-17im 1.13767e-17-3.75559e-17im -1.96104e-17-1.42265e-17im -8.08181e-18+2.00791e-17im 1.59347e-17+3.80377e-18im -5.15297e-19-1.21159e-17im -1.11563e-17+2.11254e-17im -1.06426e-17+4.48441e-19im 9.97107e-18-1.10351e-17im 3.44846e-18+1.33684e-18im -4.68096e-18+5.22695e-18im 4.39275e-18+8.91115e-18im
3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -4.10876e-18+1.59967e-17im 8.50869e-18+3.40232e-18im 1.96473e-17-2.89782e-17im -1.21529e-17-1.03447e-17im -1.82858e-17-3.37562e-17im -6.4196e-18+3.9968e-18im 7.85256e-20+1.51998e-17im 7.00086e-18+3.1146e-18im -3.51205e-18-2.34442e-18im 8.55835e-18-5.02783e-18im -2.68068e-18-8.11096e-18im 8.41322e-19+1.68969e-17im 480.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im
1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -4.10876e-18+1.59967e-17im 8.50869e-18+3.40232e-18im 1.96473e-17-2.89782e-17im -1.21529e-17-1.03447e-17im -1.82858e-17-3.37562e-17im -6.4196e-18+3.9968e-18im 7.85256e-20+1.51998e-17im 7.00086e-18+3.1146e-18im -3.51205e-18-2.34442e-18im 8.55835e-18-5.02783e-18im -2.68068e-18-8.11096e-18im 0.5+9.98774e-18im 450.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im
3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -4.10876e-18+1.59967e-17im 8.50869e-18+3.40232e-18im 1.96473e-17-2.89782e-17im -1.21529e-17-1.03447e-17im -1.82858e-17-3.37562e-17im -6.4196e-18+3.9968e-18im 7.85256e-20+1.51998e-17im 7.00086e-18+3.1146e-18im -3.51205e-18-2.34442e-18im 8.55835e-18-5.02783e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 420.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im
-2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -4.10876e-18+1.59967e-17im 8.50869e-18+3.40232e-18im 1.96473e-17-2.89782e-17im -1.21529e-17-1.03447e-17im -1.82858e-17-3.37562e-17im -6.4196e-18+3.9968e-18im 7.85256e-20+1.51998e-17im 7.00086e-18+3.1146e-18im -3.51205e-18-2.34442e-18im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 392.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im
6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -4.10876e-18+1.59967e-17im 8.50869e-18+3.40232e-18im 1.96473e-17-2.89782e-17im -1.21529e-17-1.03447e-17im -1.82858e-17-3.37562e-17im -6.4196e-18+3.9968e-18im 7.85256e-20+1.51998e-17im 7.00086e-18+3.1146e-18im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 364.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im
1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -4.10876e-18+1.59967e-17im 8.50869e-18+3.40232e-18im 1.96473e-17-2.89782e-17im -1.21529e-17-1.03447e-17im -1.82858e-17-3.37562e-17im -6.4196e-18+3.9968e-18im 7.85256e-20+1.51998e-17im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 338.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im
-2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -4.10876e-18+1.59967e-17im 8.50869e-18+3.40232e-18im 1.96473e-17-2.89782e-17im -1.21529e-17-1.03447e-17im -1.82858e-17-3.37562e-17im -6.4196e-18+3.9968e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 312.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im
-5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -4.10876e-18+1.59967e-17im 8.50869e-18+3.40232e-18im 1.96473e-17-2.89782e-17im -1.21529e-17-1.03447e-17im -1.82858e-17-3.37562e-17im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 288.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im
1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -4.10876e-18+1.59967e-17im 8.50869e-18+3.40232e-18im 1.96473e-17-2.89782e-17im -1.21529e-17-1.03447e-17im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 264.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im
-9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -4.10876e-18+1.59967e-17im 8.50869e-18+3.40232e-18im 1.96473e-17-2.89782e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 242.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im
-3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -4.10876e-18+1.59967e-17im 8.50869e-18+3.40232e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 220.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im
3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -4.10876e-18+1.59967e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 200.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im
-4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -4.53964e-18+4.47079e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 180.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im
6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im -2.33026e-17-1.55125e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 162.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im
-2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 1.48727e-18+1.13135e-17im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 144.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im
-2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -9.71381e-19-1.04378e-17im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 128.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im
7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 5.46006e-18+1.88734e-18im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 112.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im
2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im -3.63007e-18+1.40327e-17im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 98.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im
-3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im 1.35822e-18-1.50482e-17im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 84.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im
-2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 2.35457e-17-1.83632e-17im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 72.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im
2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im 1.13767e-17+3.75559e-17im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 60.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im
1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -1.96104e-17+1.42265e-17im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 50.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im
4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im -8.08181e-18-2.00791e-17im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 40.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im
3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im 1.59347e-17-3.80377e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 32.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im
-2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -5.15297e-19+1.21159e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 24.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im
1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 0.0-1.67108e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 18.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im
1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im -1.06426e-17-4.48441e-19im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 12.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im
-1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im 1.78543e-17-1.72668e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 8.0+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im
-2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.44846e-18-1.33684e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 4.5+0.0im 0.5+1.16569e-17im 8.97271e-18-4.91711e-18im
8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im -4.68096e-18-5.22695e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 2.0+0.0im 0.5+1.16569e-17im
0.5+1.16569e-17im 8.97271e-18-4.91711e-18im -2.04669e-18-7.85756e-18im -1.01526e-17+1.97211e-17im 1.0213e-17+1.52664e-17im 1.23204e-17+9.68665e-21im -2.11585e-18-3.84134e-18im 3.64735e-18+1.55332e-17im 4.20171e-18-1.36319e-18im 1.04751e-17-1.37601e-17im 2.08131e-17+5.36003e-18im -2.55944e-18-3.23188e-17im -3.42569e-17+3.63887e-18im 2.40113e-18+2.87344e-17im 7.23023e-19-4.03029e-18im -2.85334e-18+5.92892e-18im -2.81457e-18+1.0037e-17im 6.26015e-18-6.40017e-20im -4.02022e-18-1.83128e-17im 3.26871e-19-3.07798e-18im -3.28243e-18+3.61775e-18im -9.47944e-18-1.56453e-18im 1.81844e-17+2.83749e-18im -5.03637e-18-1.0721e-18im -2.25513e-17-1.21393e-17im 1.17651e-17+2.33973e-17im 6.59538e-18+2.3907e-18im -2.63164e-17+7.46706e-18im 3.99783e-18-3.16376e-18im 1.17666e-17+6.83303e-18im 3.0551e-18+4.69829e-18im 4.39275e-18-8.91115e-18im 1.17666e-17-6.83303e-18im 3.99783e-18+3.16376e-18im -1.17651e-17-7.53952e-18im 6.59538e-18-2.3907e-18im 9.02289e-18-1.67345e-17im -2.25513e-17+1.21393e-17im -5.03637e-18+1.0721e-18im 1.81844e-17-2.83749e-18im -9.47944e-18+1.56453e-18im -3.28243e-18-3.61775e-18im 3.26871e-19+3.07798e-18im -8.02303e-18+9.48533e-18im 6.26015e-18+6.40017e-20im 0.0-1.14136e-17im -2.85334e-18-5.92892e-18im 7.23023e-19+4.03029e-18im 2.40113e-18-2.87344e-17im -3.42569e-17-3.63887e-18im -2.55944e-18+3.23188e-17im 2.08131e-17-5.36003e-18im -1.03454e-17-1.16715e-17im 4.20171e-18+1.36319e-18im -8.22387e-18-4.9923e-19im -2.11585e-18+3.84134e-18im 1.23204e-17-9.68665e-21im 1.0213e-17-1.52664e-17im -1.01526e-17-1.97211e-17im -2.04669e-18+7.85756e-18im 8.97271e-18+4.91711e-18im 0.5+9.98774e-18im 0.5+0.0imCompare this matrix Array(H) with the one you obtained in Exercise 3, get its eigenvectors and eigenvalues. Try to guess the ordering of $G$-vectors in DFTK.
Increase the size of the problem, and compare the time spent by DFTK's internal diagonalization algorithms to a full diagonalization of Array(H). Hint: The @belapsed and @benchmark macros (from the BenchmarkTools package) are handy for this task. Note that there are some subtleties with global variables (see the BenchmarkTools docs for details). E.g. to use it to benchmark a function like eigen(H) run it as (note the $):
using BenchmarkTools
@benchmark eigen($H)Solutions
Exercise 1
If we consider a function $f : [0, 2π] → \mathbb{R}$, to first order
\[f(x + δx) = f(x) + δx f'(x) + O(δx^2)\]
therefore after rearrangement
\[f'(x) = \frac{f(x + δx) - f(x)}{δx} + O(δx).\]
Similarly
\[f''(x) = \frac{f'(x + δx) - f'(x)}{δx} + O(δx),\]
such that overall
\[f''(x) \simeq \frac{f(x + 2δx) - f(x + δx) - f(x + δx) + f(x)}{δx^2} = \frac{f(x + 2δx) - 2f(x + δx) + f(x)}{δx^2}\]
In finite differences we consider a stick basis of vectors
\[\left\{ e_i = (0, …, 0, \underbrace{δx}_\text{$i$-th position}, 0, …, 0) \middle| i = 1, … N \right\}.\]
Keeping in mind the periodic boundary conditions (i.e. $e_0 = e_N$) projecting the Hamiltonian $H$ onto this basis thus yields the proposed structure.
We start off with $N = 100$ to obtain
Hfd = build_finite_differences_matrix(cos, 100)
L, V = eigen(Hfd)
L[1:5]5-element Vector{Float64}:
-0.5351640695081414
0.3429576842587494
0.8532206682511969
2.053679933442362
2.078512110605267This is already pretty accurate (to about 4 digits) as can be estimated looking at the following convergence plot:
function fconv(N)
L, V = eigen(build_finite_differences_matrix(cos, N))
first(L)
end
Nrange = 10:10:100
plot(Nrange, abs.(fconv.(Nrange) .- fconv(200));
yaxis=:log, legend=false, mark=:x, xlabel="N", ylabel="Absolute error")Exercise 2
We note that
\[\langle e_G, e_{G'}\rangle = ∫_0^{2π} e_G^\ast(x) e_{G'}(x) d x = 1/2π ∫_0^{2π} e^{i(G'-G)x} d x\]
Since $e^{iy}$ is a periodic function with period $2\pi$, $\int_0^{2\pi} e^{i m y} = \delta_{0,m}$. Therefore if $G≠G'$ we have that $\langle e_G, e_{G'}\rangle = 0$, while $G=G'$ implies $\langle e_G, e_{G'}\rangle = 1$. In summary:
\[\langle e_G, e_{G'}\rangle = δ_{G, G'}\]
Next fo $V(x) = \cos(x)$ we obtain
\[\langle e_G, H e_{G'}\rangle = \frac 1 2 ∫_0^{2π} e_G^\ast(x) H e_{G'}(x) d x\]
We start by applying the Hamiltonian to a plane-wave:
\[H e_{G'}(x) = - \frac 1 2 (-|G|^2) \frac 1 {\sqrt{2π}} e^{iG'x) + cos(x) \frac 1 {\sqrt{2π}} e{iG'x}\]
Then, using the result of the first part of the exercise and the fact that $cos(x) = \frac 1 2 \left(e{ix} + e{-ix}\right)$, we get:
\[\begin{align*} ⟨ e_G, H e_{G'}⟩ &= \frac 1 2 G^2 δ_{G, G'} + \frac 1 {4π} \left(∫_0^{2π} e^{ix ⋅ (G'-G+1)} d x + ∫_0^{2π} e^{ix ⋅ (G'-G-1)} d x \right) \\ &= \frac 1 2 \left(|G|^2 \delta_{G,G'} + \delta_{G, G'+1} + \delta_{G, G'-1}\right) \end{align*}\]
In case a more general $V(x)$ was employed, this potential still has to be periodic over $[0, 2\pi]$ to fit our setting. Assuming sufficient regularity in $V$ we can employ a Fourier series:
\[V(x) = \sum_{G=- \infty}^{\infty} \hat{V}_G e_G(x)\]
where
\[\hat{V}_G = \frac{1}{\sqrt{2π}} ∫_0^{2π} V(x) e^{-iGx} dx = ∫_0^{2π} V(x) e_G^\ast dx .\]
Note that one can change of this as a change of basis from the position basis to the plane-wave basis.
Based on this expansion
\[\begin{align*} ⟨ e_G, V e_{G'} ⟩ &= \left\langle e_G, ∑_{G''} \hat{V}_{G''} \, e_{G'+G''} \right\rangle \\ &= \sum_{G''=-\infty}^\infty \hat{V}_{G''} ⟨ e_G, e_{G'+G''} ⟩ \\ &= \sum_{G''=-\infty}^\infty \hat{V}_{G''} \, δ_{G-G', G''} ⟩ \\ &= \hat{V}_{G-G'} \end{align*}\]
and therefore
\[⟨ e_G, H e_{G'} ⟩ = \frac 1 2 |G|^2 \delta_{G,G'} + \hat{V}_{G-G'},\]
i.e. essentially the Fourier transform of $V$ determines the contribution to the matrix elements of the Hamiltonian.
Exercise 3
The Hamiltonian matrix for the plane waves method can be found this way:
# Plane waves Hamiltonian -½Δ + cos on [0, 2pi].
function build_plane_waves_matrix_cos(N::Integer)
# Plane wave approximation to -½Δ
Gsq = [float(i)^2 for i in -N:N]
# Hamiltonian as derived in Exercise 2:
1/2 * Tridiagonal(ones(2N), Gsq, ones(2N))
end;Then we check that the first eigenvalue agrees with the finite-difference case, using $N = 10$:
Hpw_cos = build_plane_waves_matrix_cos(10)
L, V = eigen(Hpw_cos)
L[1:5]5-element Vector{Float64}:
-0.5350648522878154
0.34336012839908336
0.8536343543207997
2.0565044112660984
2.081227363352144We look at the convergence plot to compare the accuracy for various numbers of plane-waves $N$:
function fconv(N)
L = eigvals(build_plane_waves_matrix_cos(N))
first(L)
end
Nrange = 2:10
plot(Nrange, abs.(fconv.(Nrange) .- fconv(200)); yaxis=:log, legend=false,
ylims=(1e-16,Inf), ylabel="Absolute error", xlabel="N", mark=:x)Notice how compared to exercise 1 the considered basis size $N$ is much smaller, indicating that plane-wave methods more quickly lead to accurate solutions than finite-difference methods.
Exercise 4
For efficiency reasons the data in Fourier space is not ordered increasingly with $G$. Therefore to plot the Fourier space representation sensibly, we need to sort by ascending values of the $G$ vectors first. For this we extract the Fourier vector of each plane-wave basis function in the index order:
coords_G_vectors = G_vectors_cart(basis, basis.kpoints[1]) # Get coordinates of first and only k-point
# Only keep first component of each vector (because the others are zero for 1D problems):
coords_Gx = [G[1] for G in coords_G_vectors]
p = plot(coords_Gx, real(ψ_fourier); label="real part", xlims=(-10, 10))
plot!(p, coords_Gx, imag(ψ_fourier); label="imaginary part")The plot is symmetric about the zero (confirming that the orbitals are real) and only takes peaked values, which corresponds to the expected result for a cosine potential.
Exercise 5
To figure out the ordering we consider a small basis and build the Hamiltonian:
basis_small = PlaneWaveBasis(model; Ecut=5, kgrid=(1, 1, 1))
ham_small = Hamiltonian(basis_small)
H_small = Array(ham_small.blocks[1])
H_small[abs.(H_small) .< 1e-12] .= 0 # Drop numerically zero entries
H_small7×7 Matrix{ComplexF64}:
0.0+0.0im 0.5+1.00424e-16im 0.0+0.0im 0.0+0.0im 0.0+0.0im 0.0+0.0im 0.5-1.00424e-16im
0.5-1.00424e-16im 0.5+0.0im 0.5+1.00424e-16im 0.0+0.0im 0.0+0.0im 0.0+0.0im 0.0+0.0im
0.0+0.0im 0.5-1.00424e-16im 2.0+0.0im 0.5+1.00424e-16im 0.0+0.0im 0.0+0.0im 0.0+0.0im
0.0+0.0im 0.0+0.0im 0.5-1.00424e-16im 4.5+0.0im 0.0+0.0im 0.0+0.0im 0.0+0.0im
0.0+0.0im 0.0+0.0im 0.0+0.0im 0.0+0.0im 4.5+0.0im 0.5+1.00424e-16im 0.0+0.0im
0.0+0.0im 0.0+0.0im 0.0+0.0im 0.0+0.0im 0.5-1.00424e-16im 2.0+0.0im 0.5+1.00424e-16im
0.5+1.00424e-16im 0.0+0.0im 0.0+0.0im 0.0+0.0im 0.0+0.0im 0.5-1.00424e-16im 0.5+0.0imThe equivalent version using the build_plane_waves_matrix_cos function is N=3 (both give rice to a 7×7 matrix).
Hother = build_plane_waves_matrix_cos(3)7×7 LinearAlgebra.Tridiagonal{Float64, Vector{Float64}}:
4.5 0.5 ⋅ ⋅ ⋅ ⋅ ⋅
0.5 2.0 0.5 ⋅ ⋅ ⋅ ⋅
⋅ 0.5 0.5 0.5 ⋅ ⋅ ⋅
⋅ ⋅ 0.5 0.0 0.5 ⋅ ⋅
⋅ ⋅ ⋅ 0.5 0.5 0.5 ⋅
⋅ ⋅ ⋅ ⋅ 0.5 2.0 0.5
⋅ ⋅ ⋅ ⋅ ⋅ 0.5 4.5By comparing the entries we find the ordering is 0,1,2,...,-2,-1, which can also be found by inspecting
first.(G_vectors(basis_small, basis_small.kpoints[1]))7-element Vector{Int64}:
0
1
2
3
-3
-2
-1Both matrices have the same eigenvalues:
maximum(abs, eigvals(H_small) - eigvals(Hother))1.7763568394002505e-15and in the eigenvectors we find the same rearrangements in the entries of the eigenvectors of both matrices, matching the DFTK ordering of is 0,1,2,...,-2,-1.
eigvecs(Hother)[:, 1]7-element Vector{Float64}:
-0.008461095314567344
0.08520425570524764
-0.42353481192309905
0.7915641575577494
-0.42353481192310566
0.08520425570524905
-0.008461095314567471eigvecs(H_small)[:, 1]7-element Vector{ComplexF64}:
-0.7915641575577493 + 1.5898485476322146e-16im
0.42353481192310877 - 1.7013307112974292e-16im
-0.08520425570525011 + 5.1339563909648e-17im
0.008461095314567624 - 6.797609502050578e-18im
0.008461095314567069 + 3.3988047510250303e-18im
-0.08520425570524359 - 1.711318796988132e-17im
0.4235348119230965 + 0.0imNotice, that eigenvectors are only defined up to a phase, so the sign may globally be inverted between the two eigenvectors.
Exercise 6
We benchmark the time needed for a full diagonalization (instantiation of the Array plus call of eigen) versus the time needed for running the SCF (i.e. iterative diagonalization using plane waves).
using Printf
for Ecut in 200:200:1600
basis_time = PlaneWaveBasis(model; Ecut, kgrid=(1, 1, 1))
t_eigen = @elapsed eigen(Array(Hamiltonian(basis_time).blocks[1]))
t_scf = @elapsed self_consistent_field(basis_time; tol=1e-6, callback=identity);
@printf "%4i eigen=%8.4fms scf=%8.4fms\n" Ecut 1000t_eigen 1000t_scf
end 200 eigen= 2.4230ms scf= 5.1452ms
400 eigen= 4.4023ms scf= 3.3391ms
600 eigen= 6.4406ms scf= 3.1903ms
800 eigen= 9.9584ms scf= 3.2343ms
1000 eigen= 11.9109ms scf= 3.4995ms
1200 eigen= 14.1403ms scf= 3.2301ms
1400 eigen= 16.4969ms scf= 3.1009ms
1600 eigen= 19.6187ms scf= 3.3643ms