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"##359".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.028094928642132002 - 0.12273790890839915im
0.027994637898495534 - 0.09946041701137505im
0.1592790311275656 - 0.081400772788719im
0.042211427709555785 + 0.322008849124571im
-0.16760319962168657 - 0.2578554313780062im
-0.6264449397688965 - 1.313033747226755im
-1.0632057198522644 + 0.7557670187632427im
5.014164337415923 + 0.622593962279688im
4.1936146140354245 - 2.012181199575027im
2.0634995006601056 + 1.134664480529986im
-3.2286187776401167 - 3.9955585215852065im
-4.056861544955681 + 2.0807839210867956im
-1.54144553159435 + 3.2359962574130066im
6.91668827608653 + 1.5471686338217006im
3.1607620774665453 - 4.71612046489667im
-18.59481544848212 - 5.79103968630993im
1.4609334896061028 + 16.17419808022199im
-1.408097705231511 + 12.574536868285623im
-10.099418968869916 + 6.229433083693041im
35.39701639702648 + 5.8846078604024985im
-5.3821210638767445 - 34.17894587634089im
37.612151578535716 + 30.58509887700625im
-49.627537477546696 - 27.016630544192658im
5.797553412868598 + 20.369200134708638im
-33.59618564581113 + 4.821625400392442im
25.163051657964242 - 10.643879556965336im
-2.943923136959409 + 12.0609445012144im
14.435280936314646 + 39.243826689525115im
8.763042421525103 + 10.257388665843877im
-16.9529490785266 - 20.779383836149417im
-84.49024065418666 + 6.6346154521691885im
-19.711766328896132 - 21.710560511148994im
69.2013405837688 + 50.61068742522219im
1.2172780562711023 + 30.21467563868524im
-8.576623158897565 + 16.588058873851324im
-48.11870878506111 + 17.20812692436182im
-20.13562337473263 + 5.428444813104489im
-21.680377600496097 + 0.9291053132518675im
36.28452922262856 - 3.4217847527126968im
23.352759709989723 + 49.98284319096013im
-0.4241844184984526 - 31.193294213568866im
-16.554333941220985 + 21.348659491068766im
-2.75465833678347 - 20.81765389649413im
-33.254073012056644 + 31.05981816896098im
-9.631364402633078 + 4.975044093349325im
-14.23552568615056 + 10.559433904526033im
10.912148419557205 - 20.876633011646998im
-8.74700134098445 - 13.457989320735484im
6.521047009716927 - 3.6828228021707816im
12.865303709472148 + 4.481815467116007im
0.507917013825747 + 0.03746643142520515im
-1.5528689342397342 - 1.2437939748270552im
-1.8185792315744471 - 10.80457200020205im
2.543252396299133 + 6.7113659328036785im
-3.6690332852344754 + 1.6394879145626666im
0.842496876615133 + 4.054460203484381im
2.5634432810899423 + 1.6110204528266676im
4.159329283801381 + 1.4059831901981839im
-1.2801727035983963 - 0.19180695149915072im
-0.16791352399109166 + 0.754445963290641im
0.3783974137375403 + 0.3962087850504583im
0.17154129932634743 + 0.28282534335318205im
0.02771557931455349 + 0.08764608729029916imWe can also get its full matrix representation:
Array(H)63×63 Matrix{ComplexF64}:
9.38036e-19+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im
0.5+2.37729e-18im 0.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im
6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 2.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im
-2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 4.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im
9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 8.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im
5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 12.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im
1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 18.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im
-6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 24.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im
6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 32.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im
5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 40.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im
-5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 50.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im
1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 60.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im
-3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 72.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im
-4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 84.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im
6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 98.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im
1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 112.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im
3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 128.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im
4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 144.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im
6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 162.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im
-6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 180.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im
-2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 200.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.96529e-18-1.63017e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im
-4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 220.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 7.46908e-18-4.89191e-18im -3.96529e-18-1.63017e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im
-1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 242.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im 2.3641e-17+8.69902e-18im 7.46908e-18-4.89191e-18im -3.96529e-18-1.63017e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im
1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 264.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im -1.24221e-17+9.18236e-18im 2.3641e-17+8.69902e-18im 7.46908e-18-4.89191e-18im -3.96529e-18-1.63017e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im
-4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 288.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im -5.05513e-18+9.40828e-18im -1.24221e-17+9.18236e-18im 2.3641e-17+8.69902e-18im 7.46908e-18-4.89191e-18im -3.96529e-18-1.63017e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im
-2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 312.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.80457e-18-1.65903e-18im -5.05513e-18+9.40828e-18im -1.24221e-17+9.18236e-18im 2.3641e-17+8.69902e-18im 7.46908e-18-4.89191e-18im -3.96529e-18-1.63017e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im
1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 338.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 2.94354e-19-1.42818e-17im -6.80457e-18-1.65903e-18im -5.05513e-18+9.40828e-18im -1.24221e-17+9.18236e-18im 2.3641e-17+8.69902e-18im 7.46908e-18-4.89191e-18im -3.96529e-18-1.63017e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im
6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 364.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 1.25324e-17+7.3538e-18im 2.94354e-19-1.42818e-17im -6.80457e-18-1.65903e-18im -5.05513e-18+9.40828e-18im -1.24221e-17+9.18236e-18im 2.3641e-17+8.69902e-18im 7.46908e-18-4.89191e-18im -3.96529e-18-1.63017e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im
-7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 392.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 2.40812e-18-1.91415e-17im 1.25324e-17+7.3538e-18im 2.94354e-19-1.42818e-17im -6.80457e-18-1.65903e-18im -5.05513e-18+9.40828e-18im -1.24221e-17+9.18236e-18im 2.3641e-17+8.69902e-18im 7.46908e-18-4.89191e-18im -3.96529e-18-1.63017e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im
0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 420.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im 9.85983e-18+4.45332e-18im 2.40812e-18-1.91415e-17im 1.25324e-17+7.3538e-18im 2.94354e-19-1.42818e-17im -6.80457e-18-1.65903e-18im -5.05513e-18+9.40828e-18im -1.24221e-17+9.18236e-18im 2.3641e-17+8.69902e-18im 7.46908e-18-4.89191e-18im -3.96529e-18-1.63017e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im
1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 450.0+0.0im 0.5-1.91936e-18im 4.3417e-18+1.90195e-17im 9.85983e-18+4.45332e-18im 2.40812e-18-1.91415e-17im 1.25324e-17+7.3538e-18im 2.94354e-19-1.42818e-17im -6.80457e-18-1.65903e-18im -5.05513e-18+9.40828e-18im -1.24221e-17+9.18236e-18im 2.3641e-17+8.69902e-18im 7.46908e-18-4.89191e-18im -3.96529e-18-1.63017e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im -2.32476e-18-1.03486e-17im
-2.32476e-18-1.03486e-17im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 480.5+0.0im -6.21078e-18-8.87719e-18im 4.3417e-18+1.90195e-17im 9.85983e-18+4.45332e-18im 2.40812e-18-1.91415e-17im 1.25324e-17+7.3538e-18im 2.94354e-19-1.42818e-17im -6.80457e-18-1.65903e-18im -5.05513e-18+9.40828e-18im -1.24221e-17+9.18236e-18im 2.3641e-17+8.69902e-18im 7.46908e-18-4.89191e-18im -3.96529e-18-1.63017e-17im -6.5294e-18-5.68244e-18im -2.09005e-17+1.30029e-17im 1.30228e-18-1.05794e-17im 7.50769e-18+1.8538e-17im 7.79391e-18-1.28996e-18im -4.55526e-18-1.45669e-17im -2.07003e-17+3.24918e-17im 2.76834e-17+2.04324e-17im 1.26093e-17-3.72139e-17im -1.9652e-17-2.7462e-17im -5.50336e-18+2.58911e-17im 1.50263e-17+4.32004e-18im 2.08933e-18-9.22629e-18im -2.1533e-17+1.12123e-17im -1.02584e-17-3.88425e-20im 5.20786e-18+2.01367e-18im 1.19815e-17-3.24554e-18im -4.04378e-18+5.31047e-18im -5.18221e-18+5.03797e-20im
-8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.96529e-18+1.63017e-17im 7.46908e-18+4.89191e-18im 2.3641e-17-8.69902e-18im -1.24221e-17-9.18236e-18im -5.05513e-18-9.40828e-18im -6.80457e-18+1.65903e-18im 2.94354e-19+1.42818e-17im 1.25324e-17-7.3538e-18im -4.44267e-19-4.10411e-18im 9.85983e-18-4.45332e-18im 8.49259e-18-1.66418e-17im -6.21078e-18+8.87719e-18im 480.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im
1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.96529e-18+1.63017e-17im 7.46908e-18+4.89191e-18im 2.3641e-17-8.69902e-18im -1.24221e-17-9.18236e-18im -5.05513e-18-9.40828e-18im -6.80457e-18+1.65903e-18im 2.94354e-19+1.42818e-17im 1.25324e-17-7.3538e-18im -4.44267e-19-4.10411e-18im 9.85983e-18-4.45332e-18im 8.49259e-18-1.66418e-17im 0.5+2.37729e-18im 450.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im
-8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.96529e-18+1.63017e-17im 7.46908e-18+4.89191e-18im 2.3641e-17-8.69902e-18im -1.24221e-17-9.18236e-18im -5.05513e-18-9.40828e-18im -6.80457e-18+1.65903e-18im 2.94354e-19+1.42818e-17im 1.25324e-17-7.3538e-18im -4.44267e-19-4.10411e-18im 9.85983e-18-4.45332e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 420.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im
-7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.96529e-18+1.63017e-17im 7.46908e-18+4.89191e-18im 2.3641e-17-8.69902e-18im -1.24221e-17-9.18236e-18im -5.05513e-18-9.40828e-18im -6.80457e-18+1.65903e-18im 2.94354e-19+1.42818e-17im 1.25324e-17-7.3538e-18im -4.44267e-19-4.10411e-18im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 392.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im
6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.96529e-18+1.63017e-17im 7.46908e-18+4.89191e-18im 2.3641e-17-8.69902e-18im -1.24221e-17-9.18236e-18im -5.05513e-18-9.40828e-18im -6.80457e-18+1.65903e-18im 2.94354e-19+1.42818e-17im 1.25324e-17-7.3538e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 364.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im
1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.96529e-18+1.63017e-17im 7.46908e-18+4.89191e-18im 2.3641e-17-8.69902e-18im -1.24221e-17-9.18236e-18im -5.05513e-18-9.40828e-18im -6.80457e-18+1.65903e-18im 2.94354e-19+1.42818e-17im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 338.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im
-2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.96529e-18+1.63017e-17im 7.46908e-18+4.89191e-18im 2.3641e-17-8.69902e-18im -1.24221e-17-9.18236e-18im -5.05513e-18-9.40828e-18im -6.80457e-18+1.65903e-18im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 312.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im
-4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.96529e-18+1.63017e-17im 7.46908e-18+4.89191e-18im 2.3641e-17-8.69902e-18im -1.24221e-17-9.18236e-18im -5.05513e-18-9.40828e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 288.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im
1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.96529e-18+1.63017e-17im 7.46908e-18+4.89191e-18im 2.3641e-17-8.69902e-18im -1.24221e-17-9.18236e-18im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 264.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im
-1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.96529e-18+1.63017e-17im 7.46908e-18+4.89191e-18im 2.3641e-17-8.69902e-18im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 242.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im
-4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.96529e-18+1.63017e-17im 7.46908e-18+4.89191e-18im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 220.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im
-2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.96529e-18+1.63017e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 200.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im
-6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -6.5294e-18+5.68244e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 180.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im
6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im -2.09005e-17-1.30029e-17im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 162.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im
4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.30228e-18+1.05794e-17im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 144.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im
1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 7.50769e-18-1.8538e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 128.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im
1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 3.42423e-18+1.3351e-17im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 112.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im
6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im -4.55526e-18+1.45669e-17im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 98.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im
-4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -7.79391e-18-1.28996e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 84.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im
-3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im 2.76834e-17-2.04324e-17im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 72.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im
1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im 1.26093e-17+3.72139e-17im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 60.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im
-5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.9652e-17+2.7462e-17im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 50.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im
5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im -5.50336e-18-2.58911e-17im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 40.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im
6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im 1.50263e-17-4.32004e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 32.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im
-6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im 2.08933e-18+9.22629e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 24.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im
1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im -2.1533e-17-1.12123e-17im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 18.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im
5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im -1.02584e-17+3.88425e-20im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 12.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im
9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im 5.20786e-18-2.01367e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 8.0+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im
-2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 1.19815e-17+3.24554e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 4.5+0.0im 0.5-1.91936e-18im 6.84621e-18-4.62374e-18im
6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im -4.04378e-18-5.31047e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-18im 2.0+0.0im 0.5-1.91936e-18im
0.5-1.91936e-18im 6.84621e-18-4.62374e-18im -2.10986e-18-7.33721e-18im 9.53855e-19+2.05757e-17im 5.30139e-18+9.96275e-18im 1.15239e-17-1.06252e-18im -6.09606e-18-2.01355e-17im 6.73916e-18-9.6349e-18im 5.08456e-18-5.8805e-19im -5.68698e-18+2.3003e-19im 1.48741e-17+1.10137e-17im -3.44396e-18-3.19427e-17im -4.12338e-17-4.82849e-18im 6.96319e-18+3.06423e-17im 1.41552e-18-4.45966e-18im 1.06844e-17+1.07046e-17im 4.92777e-18+5.57465e-18im 6.87119e-18-1.51362e-19im -6.53743e-18-2.77989e-18im -2.28337e-18-2.65376e-18im -4.23963e-18+3.77063e-18im -1.07755e-17-3.20506e-18im 1.64216e-17-3.33489e-18im -4.8805e-18-2.16383e-18im -2.56782e-17-4.87146e-18im 1.35838e-17+5.69649e-18im 6.55612e-18+3.83412e-18im -7.57355e-18-5.37897e-18im -8.89745e-18+9.6317e-18im 1.13146e-17+6.44402e-18im -8.98722e-18-7.74616e-18im -5.18221e-18-5.03797e-20im 1.13146e-17-6.44402e-18im 0.0+6.78716e-18im -7.57355e-18+5.37897e-18im 6.55612e-18-3.83412e-18im 1.35838e-17-5.69649e-18im -2.56782e-17+4.87146e-18im -4.8805e-18+2.16383e-18im 1.64216e-17+3.33489e-18im -1.07755e-17+3.20506e-18im -4.23963e-18-3.77063e-18im -2.28337e-18+2.65376e-18im -6.53743e-18+2.77989e-18im 6.87119e-18+1.51362e-19im 4.92777e-18-5.57465e-18im 3.61585e-19-1.24616e-17im 1.41552e-18+4.45966e-18im 6.5791e-18-3.57518e-18im -4.12338e-17+4.82849e-18im -3.44396e-18+3.19427e-17im 1.48741e-17-1.10137e-17im -5.68698e-18-2.3003e-19im 5.08456e-18+5.8805e-19im 6.73916e-18+9.6349e-18im -6.09606e-18+2.01355e-17im 1.15239e-17+1.06252e-18im 5.30139e-18-9.96275e-18im 9.53855e-19-2.05757e-17im -2.10986e-18+7.33721e-18im 6.84621e-18+4.62374e-18im 0.5+2.37729e-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.3429576842590336
0.8532206682510264
2.0536799334424938
2.078512110605256This 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.701330711297429e-16im
-0.08520425570525011 + 5.1339563909647985e-17im
0.008461095314567624 - 6.797609502050588e-18im
0.008461095314567069 + 3.398804751025045e-18im
-0.08520425570524359 - 1.711318796988133e-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.4005ms scf= 5.6203ms
400 eigen= 4.4907ms scf= 3.2870ms
600 eigen= 7.3462ms scf= 3.4808ms
800 eigen= 10.7764ms scf= 3.7005ms
1000 eigen= 13.6063ms scf= 3.6196ms
1200 eigen= 17.9297ms scf= 3.6228ms
1400 eigen= 19.9409ms scf= 3.8574ms
1600 eigen= 22.8909ms scf= 3.8303ms