Added nelec as command line argument

This commit is contained in:
Anthony Scemama 2021-03-10 14:52:04 +01:00
commit 04113a93d5
13 changed files with 17971 additions and 130 deletions

View File

@ -1,7 +1,7 @@
IRPF90 = irpf90/bin/irpf90 # --codelet=jastrow_full:1000 #-s nelec:10 -s nnuc:2 -s ncord:5 #-a -d
FC = ifort -xHost -g -mkl=sequential
IRPF90 = irpf90/bin/irpf90 --codelet=factor_een:2 #-s nelec:10 -s nnuc:2 -s ncord:5 #-a -d
FC = ifort -xCORE-AVX512 -g -mkl=sequential
FCFLAGS= -O2 -I .
AR = ar
NINJA = ninja
ARCHIVE = ar crs
RANLIB = ranlib

View File

@ -34,3 +34,22 @@
\sum_{i=1}^{Ne} \sum_{pkl} \sum_a^{Nn}
c_{apkl} R_{ia}^{p-k-l}\, C_{i,a(p-k+l)}^k
$$
* Running
#+begin_src bash :var Ratio=5 Natoms=500
python ./generateData.py -a $Natoms -r $Ratio
#+end_src
Cela genere les trois fichiers:
geometry.txt
elec_coords.txt
jast_coeffs.txt
#+begin_src bash :var Natoms=500
./codelet_factor_een_blas $Natoms
#+end_src

View File

@ -1,24 +1,26 @@
program codelet_factor_een
program codelet_factor_een_blas
implicit none
integer :: i
double precision :: ticks_0, ticks_1, cpu_0, cpu_1
integer, parameter :: irp_imax = 100000
integer, parameter :: irp_imax = 200
call provide_factor_een
PROVIDE factor_een_blas tmp_c
call provide_factor_een_blas
double precision :: irp_rdtsc
call cpu_time(cpu_0)
ticks_0 = irp_rdtsc()
do i=1,irp_imax
call bld_factor_een
call bld_tmp_c
call bld_factor_een_blas
enddo
ticks_1 = irp_rdtsc()
call cpu_time(cpu_1)
print *, 'factor_een'
print *, 'factor_een_blas'
print *, '-----------'
print *, 'Cycles:'
print *, (ticks_1-ticks_0)/dble(irp_imax)
@ -26,4 +28,4 @@ program codelet_factor_een
print *, (cpu_1-cpu_0)/dble(irp_imax)
end

View File

@ -95,24 +95,24 @@ BEGIN_PROVIDER [ double precision, factor_een_deriv_e, (4, nelec) ]
rescale_een_e(i,j,k) * &
rescale_een_n(i,a,m+l)
daccu(1:4) = daccu(1:4) + &
rescale_een_e_deriv_e_t(1:4,i,j,k) * &
rescale_een_e_deriv_e_t(i,1:4,j,k) * &
rescale_een_n(i,a,m)
daccu2(1:4) = daccu2(1:4) + &
rescale_een_e_deriv_e_t(1:4,i,j,k) * &
rescale_een_e_deriv_e_t(i,1:4,j,k) * &
rescale_een_n(i,a,m+l)
enddo
factor_een_deriv_e(1:4,j) = factor_een_deriv_e(1:4,j) + &
(accu * rescale_een_n_deriv_e(1:4,j,a,m+l) + daccu(1:4) * rescale_een_n(j,a,m+l) +&
daccu2(1:4)* rescale_een_n(j,a,m) + accu2*rescale_een_n_deriv_e(1:4,j,a,m)) * cn
(accu * rescale_een_n_deriv_e(j,1:4,a,m+l) + daccu(1:4) * rescale_een_n(j,a,m+l) +&
daccu2(1:4)* rescale_een_n(j,a,m) + accu2*rescale_een_n_deriv_e(j,1:4,a,m)) * cn
factor_een_deriv_e(4,j) = factor_een_deriv_e(4,j) + 2.d0*( &
daccu (1) * rescale_een_n_deriv_e(1,j,a,m+l) + &
daccu (2) * rescale_een_n_deriv_e(2,j,a,m+l) + &
daccu (3) * rescale_een_n_deriv_e(3,j,a,m+l) + &
daccu2(1) * rescale_een_n_deriv_e(1,j,a,m ) + &
daccu2(2) * rescale_een_n_deriv_e(2,j,a,m ) + &
daccu2(3) * rescale_een_n_deriv_e(3,j,a,m ) )*cn
daccu (1) * rescale_een_n_deriv_e(j,1,a,m+l) + &
daccu (2) * rescale_een_n_deriv_e(j,2,a,m+l) + &
daccu (3) * rescale_een_n_deriv_e(j,3,a,m+l) + &
daccu2(1) * rescale_een_n_deriv_e(j,1,a,m ) + &
daccu2(2) * rescale_een_n_deriv_e(j,2,a,m ) + &
daccu2(3) * rescale_een_n_deriv_e(j,3,a,m ) )*cn
enddo
enddo
enddo
@ -152,8 +152,8 @@ BEGIN_PROVIDER [ double precision, factor_een_deriv_e_ref, (4, nelec) ]
rjam_cn = rescale_een_n(j, a, m) * cn
do ii = 1, 4
drjal(ii) = rescale_een_n_deriv_e(ii, j, a, l)
drjam_cn(ii) = rescale_een_n_deriv_e(ii, j, a, m) * cn
drjal(ii) = rescale_een_n_deriv_e(j, ii, a, l)
drjam_cn(ii) = rescale_een_n_deriv_e(j, ii, a, m) * cn
enddo
do i = 1, nelec
@ -162,7 +162,7 @@ BEGIN_PROVIDER [ double precision, factor_een_deriv_e_ref, (4, nelec) ]
rijk = rescale_een_e(i, j, k)
do ii = 1, 4
drijk(ii) = rescale_een_e_deriv_e(ii, j, i, k)
drijk(ii) = rescale_een_e_deriv_e(j, ii, i, k)
enddo
v1 = rijk * rial ! v(x)

View File

@ -1,19 +1,15 @@
BEGIN_PROVIDER [ double precision, factor_een_blas ]
&BEGIN_PROVIDER [ double precision, factor_een_deriv_e_blas, (4, nelec) ]
BEGIN_PROVIDER [ double precision, tmp_c, (nelec,nnuc,0:ncord,0:ncord-1) ]
&BEGIN_PROVIDER [ double precision, dtmp_c, (nelec,4,nnuc,0:ncord,0:ncord-1) ]
implicit none
BEGIN_DOC
! Dimensions 1-3 : dx, dy, dz
! Dimension 4 : d2x + d2y + d2z
! Calculate the intermediate buffers
! tmp_c:
! r_{ij}^k . R_{ja}^l -> tmp_c_{ia}^{kl}
!
! dtmp_c:
! dr_{ij}^k . R_{ja}^l -> dtmp_c_{ia}^{kl}
END_DOC
integer :: i, j, a, p, k, l, lmax, m, n
double precision :: cn(ncord), accu
double precision :: f(nnuc,0:ncord-2,0:ncord-2)
double precision :: tmp_c(nelec,nnuc,0:ncord,0:ncord-1)
double precision :: dtmp_c(4,nelec,nnuc,0:ncord,0:ncord-1)
factor_een_blas = 0.0d0
factor_een_deriv_e_blas(1:4,1:nelec) = 0.0d0
integer :: k
! r_{ij}^k . R_{ja}^l -> tmp_c_{ia}^{kl}
do k=0,ncord-1
@ -26,12 +22,30 @@
! dr_{ij}^k . R_{ja}^l -> dtmp_c_{ia}^{kl}
do k=0,ncord-1
call dgemm('N','N', 4*nelec, nnuc*(ncord+1), nelec, 1.d0, &
rescale_een_e_deriv_e(1,1,1,k), 4*size(rescale_een_e_deriv_e,2),&
rescale_een_e_deriv_e(1,1,1,k), 4*size(rescale_een_e_deriv_e,1),&
rescale_een_n(1,1,0), size(rescale_een_n,1), 0.d0, &
dtmp_c(1,1,1,0,k), 4*size(dtmp_c,2))
dtmp_c(1,1,1,0,k), 4*size(dtmp_c,1))
enddo
END_PROVIDER
BEGIN_PROVIDER [ double precision, factor_een_blas ]
&BEGIN_PROVIDER [ double precision, factor_een_deriv_e_blas, (nelec,4) ]
implicit none
BEGIN_DOC
! Dimensions 1-3 : dx, dy, dz
! Dimension 4 : d2x + d2y + d2z
END_DOC
integer :: i, j, a, p, k, l, lmax, m, n, ii
double precision :: accu, cn
! double precision,dimension(:),allocatable :: cn
factor_een_blas = 0.0d0
factor_een_deriv_e_blas(1:nelec,1:4) = 0.0d0
do n = 1, dim_cord_vect
l = lkpm_of_cindex(1,n)
@ -40,33 +54,37 @@
m = lkpm_of_cindex(4,n)
do a = 1, nnuc
cn(a) = cord_vect_full(n, a)
enddo
cn = cord_vect_full(n, a)
if (cn == 0.d0) cycle
do a = 1, nnuc
accu = 0.d0
do j=1,nelec
accu = accu + rescale_een_n(j,a,m) * tmp_c(j,a,m+l,k)
factor_een_deriv_e_blas(1:4,j) = factor_een_deriv_e_blas(1:4,j) + (&
tmp_c(j,a,m,k) * rescale_een_n_deriv_e(1:4,j,a,m+l) + &
dtmp_c(1:4,j,a,m,k) * rescale_een_n(j,a,m+l) + &
dtmp_c(1:4,j,a,m+l,k) * rescale_een_n(j,a,m) + &
tmp_c(j,a,m+l,k)*rescale_een_n_deriv_e(1:4,j,a,m) &
) * cn(a)
factor_een_deriv_e_blas(4,j) = factor_een_deriv_e_blas(4,j) + (&
dtmp_c(1,j,a,m ,k) * rescale_een_n_deriv_e(1,j,a,m+l) + &
dtmp_c(2,j,a,m ,k) * rescale_een_n_deriv_e(2,j,a,m+l) + &
dtmp_c(3,j,a,m ,k) * rescale_een_n_deriv_e(3,j,a,m+l) + &
dtmp_c(1,j,a,m+l,k) * rescale_een_n_deriv_e(1,j,a,m ) + &
dtmp_c(2,j,a,m+l,k) * rescale_een_n_deriv_e(2,j,a,m ) + &
dtmp_c(3,j,a,m+l,k) * rescale_een_n_deriv_e(3,j,a,m ) &
)*cn(a)*2.d0
enddo
factor_een_blas = factor_een_blas + accu * cn(a)
factor_een_blas = factor_een_blas + accu * cn
do ii=1,4
do j=1,nelec
factor_een_deriv_e_blas(j,ii) = factor_een_deriv_e_blas(j,ii) + (&
tmp_c(j,a,m,k) * rescale_een_n_deriv_e(j,ii,a,m+l) + &
dtmp_c(j,ii,a,m,k) * rescale_een_n(j,a,m+l) + &
dtmp_c(j,ii,a,m+l,k) * rescale_een_n(j,a,m) + &
tmp_c(j,a,m+l,k)*rescale_een_n_deriv_e(j,ii,a,m) &
) * cn
enddo
enddo
cn = cn+cn
do j=1,nelec
factor_een_deriv_e_blas(j,4) = factor_een_deriv_e_blas(j,4) + (&
dtmp_c(j,1,a,m ,k) * rescale_een_n_deriv_e(j,1,a,m+l) + &
dtmp_c(j,2,a,m ,k) * rescale_een_n_deriv_e(j,2,a,m+l) + &
dtmp_c(j,3,a,m ,k) * rescale_een_n_deriv_e(j,3,a,m+l) + &
dtmp_c(j,1,a,m+l,k) * rescale_een_n_deriv_e(j,1,a,m ) + &
dtmp_c(j,2,a,m+l,k) * rescale_een_n_deriv_e(j,2,a,m ) + &
dtmp_c(j,3,a,m+l,k) * rescale_een_n_deriv_e(j,3,a,m ) &
)*cn
enddo
enddo
enddo

File diff suppressed because it is too large Load Diff

View File

@ -3,7 +3,16 @@ BEGIN_PROVIDER [ integer, nelec ]
BEGIN_DOC
! Number of electrons
END_DOC
nelec = 10
character*(32) :: buffer
integer, external :: iargc
if (iargc() == 0) then
nelec = 10
else
call getarg(1,buffer)
read(buffer,*)nelec
endif
END_PROVIDER
BEGIN_PROVIDER [ integer, nelec_up ]
@ -11,7 +20,7 @@ BEGIN_PROVIDER [ integer, nelec_up ]
BEGIN_DOC
! Number of alpha and beta electrons
END_DOC
nelec_up = 5
nelec_up = nelec/2
END_PROVIDER

60
generateData.py Executable file
View File

@ -0,0 +1,60 @@
#!/usr/bin/env python
import numpy as np
import sys, getopt
from jastlib import getCoefList, get_sphere_distribution, generateElsAroundPoints
def main(argv):
Natom = 2
Ratio = 5
inputfile = ''
outputfile = ''
try:
opts, args = getopt.getopt(argv,"ha:r:",["atom=","ratio="])
except getopt.GetoptError:
print('test.py -a <natoms> -r <ratio>')
sys.exit(2)
for opt, arg in opts:
if opt == '-h':
print('test.py -a <inputfile> -r <outputfile>')
sys.exit()
elif opt in ("-a", "--atom"):
Natom = int(arg)
elif opt in ("-r", "--ratio"):
Ratio = int(arg)
Nelec = Ratio
Nord = 5
L1 = 20.0
n = Natom # number of points
dmin = 0.1 # min dist between points
Ls = np.array([L1,L1,L1]) # lengths of the box
shift = -10.0
kappa = 2.0
filename_atom = "geometry.txt"
filename_coeffsa = "jast_coeffs.txt"
(coeffsa, coeffsb, coeffsc) = getCoefList(Nord,n)
coeffsall = np.concatenate((coeffsa,coeffsb,coeffsc))
atomList,_,_ = get_sphere_distribution(n, dmin, Ls, maxiter=1e4, allow_wall=True)
L1 = 15.0
n = Nelec # number of points
dmin = 0.1 # min dist between points
Ls = np.array([L1,L1,L1]) # lengths of the box
shift = -10.0
kappa = 2.0
filename_elec = "elec_coord.txt"
rlist = generateElsAroundPoints(n,atomList,dmin)
# Save file
np.savetxt(filename_elec,rlist)
np.savetxt(filename_atom,atomList)
np.savetxt(filename_coeffsa,coeffsall)
if __name__ == "__main__":
main(sys.argv[1:])

View File

@ -1,2 +1,500 @@
0.000000 0.000000 0.000000
0.000000 0.000000 2.059801
9.169699233477617284e+00 1.577356709974558768e+01 1.026790877215631070e+01
1.709096115690725526e+01 1.887196960617206543e+01 9.799886585728481592e+00
2.759679637155061371e+00 9.340816860053413606e-01 1.427913730962294281e+01
8.466339359065724324e+00 9.231835516571540445e-01 1.640641208255591366e+01
7.239083748223058556e+00 1.725722749662392630e+01 1.787551375153963562e+01
1.474534281245462353e+01 1.833951609869788513e+01 6.010741396284879912e-02
1.513197838597113787e+01 1.124136766210814820e-01 8.422615822432087285e+00
1.565537492181737633e+01 1.748573069812977554e+00 6.284695229046888265e+00
1.723692999612131871e+01 7.388450412650704457e+00 4.172286160665874988e+00
1.715795128560187877e+01 1.450030322639299385e+01 1.976514987287575309e+01
1.102229137605018749e+01 5.401791402001889786e+00 9.274795594529734899e+00
1.615006086012495956e+00 1.633069527897562523e+01 1.910553292984471341e+01
6.385605171265195779e+00 1.316589720661558260e+01 5.436216733982728755e+00
1.292626807035237491e+01 5.886232823823796423e+00 1.715754087186913068e+01
1.208625197377544858e+01 9.786394515437207176e+00 1.717854136000790533e+01
1.494633023918510517e+01 1.603464738388862543e+01 1.406542597486318869e+01
1.667900812370894670e+01 1.251829711646267995e+01 4.270300401191871487e+00
1.586198715643285517e+01 2.875672063917811272e-01 8.891552400836387093e+00
9.419404525581425602e+00 1.389052206291107971e+01 1.108900939317432588e+01
1.963965760373481828e+01 1.146902995603590725e+01 1.282907948294069911e+01
1.198047511787546071e+01 5.562236588668721282e+00 1.915349450607307347e+01
1.488104353517301348e+01 5.804103516304268240e+00 6.937216682401476930e+00
1.081824688317853855e+01 1.816125330714973174e+01 3.402481046984189295e+00
1.719374956854118608e+01 4.884797944310266260e+00 1.920695796409663814e+01
1.793960241583939208e+01 1.653024768810333889e+00 3.298821317147142551e+00
1.897808829921623364e+00 8.140060824774359105e+00 1.961610203110009820e+01
1.522322197862856541e+01 1.141399000556607923e+01 9.491306080383932198e+00
1.885179661581180710e+01 1.492261282193172889e+01 6.484741121232293182e+00
1.116624933697523758e+01 6.358757303413402617e+00 6.916852334485867893e+00
1.641047798611942454e+01 8.078339459447054338e+00 3.974595146357355890e+00
1.839623071904149398e+01 1.747807382087439976e+01 1.224414893914888935e+01
1.013635602912024680e+01 4.994384310269368576e+00 1.919873788632172307e+01
3.532729893816910494e+00 1.177599028381978385e+01 1.364746894301321234e+01
6.364475283610449452e+00 2.756442052672813947e+00 1.993988542501183847e+00
1.986281358109961559e+01 1.591076223093103437e+01 8.678810962712169896e+00
1.175127911665890146e+01 7.772089039972498448e+00 7.187138704784985066e+00
1.193125527485955395e+01 1.431390781406945933e+01 1.265541473195264288e+01
1.133185649869511025e+01 7.714119831904899804e+00 1.988089589887153608e+01
3.970252576323316518e+00 1.945022550869681766e+01 1.250876195144565095e+01
8.316736617812782839e+00 1.175130684010653859e+01 1.051450652667101338e+01
2.363693618139217634e+00 2.468837030925852272e+00 1.005198216788683041e-02
6.296064183630758038e+00 1.379169961939031452e+01 1.087964194835752529e+01
1.404224396433654753e+01 1.676940411433356370e+00 1.264252268440966454e+01
7.748241981417103297e-01 1.214872603120482353e+01 8.330412247841345597e-01
1.482580840274279943e+01 1.223347895464069168e+01 3.415821930923017558e+00
1.583665650690541682e+01 1.353650167723072428e+01 2.697684980051406889e+00
7.228036394543392973e+00 1.856489259785730184e+01 5.412992838630334980e+00
1.250523870392217773e+01 1.930800659024364307e+01 1.323894588602757771e+01
1.476425588956215051e+01 1.973578521040820988e+01 1.845229986910515763e+01
1.453346429546398255e+00 8.164114196818228919e+00 1.096993647314665132e+01
1.645780314591159055e+01 3.609174059640158916e+00 1.214909757599990314e+01
1.503594622796937230e+01 8.361488631403306115e+00 4.885904156983071900e+00
4.607709491539058178e+00 1.711603235809798917e+01 1.077737795893560246e+01
5.444611672020791104e+00 4.650442422860228575e+00 1.909067263900083233e+01
5.313792477754875065e-01 8.018457454298415499e+00 3.294528391053614946e+00
2.008248149095013257e+00 5.871637319393618881e-01 1.274877050725333660e+01
4.236439749525723997e+00 1.622219492454691547e+01 1.673052983670413330e+00
4.863853221184999853e+00 3.211812899760422280e+00 1.196906268597565060e+01
1.423284986943139252e+01 5.125867755205759657e-01 1.957274773922198818e+01
8.788927972275073941e+00 1.690304861735086561e+01 1.986478748330291566e+01
1.873007149410335614e+01 1.202269311771925153e+01 1.238139871086807275e+01
1.332404831179675497e+01 2.430077255068527897e+00 5.774860805419890220e+00
2.860148112180942448e+00 5.699127286772864842e-01 3.785091713927268842e+00
3.169014233130353908e+00 1.187067878487777506e+01 5.904637425708301635e+00
6.485116764957301605e+00 1.449966765623678455e+01 1.093038093646896591e+01
1.638799794515688291e+01 5.828648438494672845e+00 6.770382163340869397e+00
6.488517205429760182e+00 9.709343068893831585e+00 1.716889983416628596e+01
7.335968100444041795e+00 8.477396516109221736e+00 1.244480223455498091e+01
1.287827118726432296e+01 1.462681359150602489e+01 4.404176750645225624e+00
1.177381305500459518e+01 6.986630747738691305e+00 8.980641095899795090e+00
7.780325482042506735e+00 1.464235462946512101e+01 1.190153382636764334e+01
7.363964503444961451e+00 9.378804346496965039e+00 1.963274323415326705e+01
5.776655941050856669e+00 5.339409256429639150e+00 6.151279230717858759e+00
3.116758779185213601e+00 1.742245812485627354e+01 1.569327886491995905e+01
8.567387735155199024e+00 1.208934780212089777e+00 5.196919625657359099e+00
1.415037591844026466e+01 2.847983412165777661e+00 1.233328815020159297e+01
1.636468986461008512e+01 1.451079484590765034e+01 1.207738244144391793e+01
7.417186920223405977e+00 1.310703527774940014e+01 6.177255611674361546e+00
7.707056020550233200e+00 1.551121691246043000e+01 1.191495277078642445e+00
6.066579737290198615e+00 1.398026616942766864e+01 1.670256441921655899e+00
9.404773402560927309e+00 8.163372380824096552e-02 9.563611400839750587e+00
3.572392081411224218e-01 1.349107317825030705e+01 1.442921281054204741e+01
1.046956357964869078e+01 1.720237551282311728e+01 1.976912018536546611e+01
3.587378947280008834e+00 1.982017106695849051e+01 1.800465657877865766e+01
1.233057899646373912e+01 5.017501465684251372e+00 6.245401926461646269e+00
1.223960883456811644e+00 4.752974531379106082e+00 6.328475477314010611e-01
9.772396890948849446e+00 5.571050679520041626e+00 5.226614381600134251e+00
6.975527861247479144e+00 1.870471225537538373e+01 9.192423017315331180e+00
9.240112207612600770e+00 1.896763728569224128e+01 1.101673199128806857e+01
1.341972145993210219e+01 1.870280684111147096e+01 1.795990955568402470e+01
1.035326854711128064e+01 4.878469341937064385e+00 3.751485854957772315e+00
1.870795449017081324e+01 1.831914660403077821e+00 8.373263276422404644e+00
8.087921672582453425e+00 1.950401904800463271e+01 1.454718045256741199e+01
1.027044331124545806e+01 3.508664047516378837e+00 1.445507762924474804e+01
1.554986522192784903e+01 1.020126281964536297e+01 1.835357595579026579e+01
7.168074592017841695e+00 1.100460329646218582e+01 1.428017210088657407e+01
6.485423174327742402e+00 1.322683084053914371e+01 1.992940432667338335e+00
6.553663472634834619e+00 4.229294242312933605e+00 1.569147564716142362e+00
1.481466695650968113e+00 7.831868024012984542e+00 1.127861742440441617e+01
1.335657649179298012e+01 1.647583063010893767e+01 1.138630873736699201e+01
1.379151973752265015e+01 1.871206752648178195e+01 1.427418709907454719e+01
1.832833852719470968e+01 9.370541708219858990e+00 1.306342934326463201e+01
4.831030605796211574e+00 1.953813384789268781e+01 1.357433411806962908e+01
1.424415448690996833e+01 1.767406449802044932e+01 6.249021198885172268e+00
3.905369047332567511e+00 3.269990121459864785e+00 7.224162243389821825e+00
1.358237647699499284e+01 7.315184490172310205e+00 2.425126693107051423e+00
1.212372619074920443e+01 7.230526657229456866e+00 6.173809719590686029e+00
1.189462711607931489e+00 1.436258062726703599e+01 1.388908238946449991e+01
5.791282582461516171e+00 1.154811124642241005e+01 4.853987867958833746e+00
1.580896164277743488e+01 1.336431055791168987e+01 1.464733318205402313e-01
1.461920345301095558e+01 2.868597875720253487e+00 6.843543709812403009e+00
6.466614211741039675e-01 2.591729448327484420e-01 9.027239297261210993e+00
9.797766150642786442e+00 1.073713633697784431e+01 6.555964178001707054e+00
7.720982694389178391e+00 7.626541579051428599e+00 3.378125038476549324e+00
5.629253920787151699e+00 2.939295942554405183e+00 7.733903921451464214e+00
6.583043287686914269e+00 1.922035337287248957e+00 3.584219282118321637e+00
1.489921294361408677e+01 1.407384696930397183e+01 1.706132384981712136e+01
1.884796602727871928e+01 1.412843909897030237e+01 4.981866736853817201e+00
5.257804257464737674e+00 9.596920031787581351e-01 4.521353899301237433e+00
8.035825356431727684e+00 5.683074702560486635e+00 7.592413901629657680e+00
8.674200728650191650e+00 1.989106949092978738e+01 4.798674908083784274e+00
1.196880775847672940e+00 1.881420583740770169e+01 5.170873776613637673e+00
1.273187955638934987e+01 1.194740895765314548e+01 1.332625470554491542e+00
1.135098960471650287e+01 2.938928551817328039e+00 7.838721421271126033e+00
1.716530594140219890e+01 1.648801429086308801e+01 1.242794559996758075e+01
1.059808339130593069e+01 4.094109684779065894e+00 1.252581464872687889e+01
1.630451438891242777e+00 1.581507340905474557e+01 1.241343681908335661e-01
1.153722834561401100e+01 2.198244146456629355e+00 6.327048268880295367e+00
1.273252708789018151e+01 1.974450378774391268e+01 6.961717815856731661e+00
5.526234739443327548e+00 2.096001487241014871e+00 8.073834430983314547e+00
8.083971380937386542e+00 2.941693099249573784e+00 9.457602991575523532e+00
6.242035997829011862e+00 8.171402498789335667e+00 1.457729081929404202e+01
1.889743919435310460e+01 1.215973328705928047e+01 7.291307883355515607e+00
1.066803974906196650e+01 6.740766286962682763e+00 1.468358927758984578e+01
2.502460941375919123e+00 5.335195575479789731e-01 1.692966502639559678e+01
1.048201056139653886e+01 6.589474791250122365e+00 2.907510843947571644e+00
1.284394055156609227e+01 8.013476693909149517e-01 7.799321414157689425e+00
1.275682590461077126e+01 5.346479388172761915e+00 5.472946928985457760e+00
1.697267337457670422e+01 1.088446665380526568e+01 9.908254242298493608e+00
1.428632679128270411e+01 1.162533321997503144e+00 1.629234164043686928e+00
1.434638452207311587e+01 1.864403200474162858e+01 9.149899528133031978e+00
9.359057449688751973e+00 1.510912487929017090e+01 7.663348029212091461e-01
1.813803852725210941e+00 8.753440112058431311e+00 7.845304719674784266e+00
2.392105698224156196e+00 1.005849153808980390e+01 1.303589885117864000e+00
7.240038161511195725e+00 6.211595784793457753e+00 7.528096338530554377e-01
1.575097193267283302e+01 1.425308429393408893e+01 3.492565385916086562e+00
6.852006210223235527e+00 6.432083373615052224e+00 3.118185210866770696e+00
1.507958164934602507e+01 2.830615204387856298e+00 1.903946411331689603e+01
4.445188836935606957e+00 1.275376891718884131e+01 1.687684535924550389e+01
1.067327282018373502e+01 1.705875157642136042e+00 1.195729025803453638e+01
1.700310157905415842e-01 1.682608960965131928e+01 1.853309823648882215e+01
2.914889981242891892e+00 8.472259320255947301e+00 1.991803554830299561e+01
1.470187229496327674e+01 1.125965725688386776e+01 5.604931196696478324e+00
1.384619611683915785e+01 6.193769753542673229e+00 9.318303329658837342e+00
1.687633223865343268e+01 9.998018666633090135e+00 6.508778650551381162e+00
1.894965617366820965e+01 8.764780954071175856e+00 7.076551841313538205e+00
1.252764425472890508e+01 4.148018359600489369e+00 4.545009340633208161e+00
9.660902813795068766e+00 1.721618200071324267e+00 1.040153036026642752e+01
5.655425052058092206e+00 7.712607150642014631e+00 1.664238063140947332e+01
1.562937078944793967e+01 1.060325627156740680e+01 5.860328603923159463e-01
1.817102462266575102e+01 1.250817643177721195e+01 1.968098402390315727e+01
1.377079690279232338e+01 1.539053487312786217e+01 8.550198777069393685e+00
3.250551032209525459e+00 1.375558510358362518e+01 1.806495464619420588e+01
1.207307934606351729e+01 5.141087676208591084e+00 1.312080964285153151e+01
1.437148000527093927e+01 1.837315409754416606e+01 1.612365761140965859e+01
1.502778688073879820e+01 1.173909706723641300e+01 4.438244572471181648e+00
1.055758228475412608e+01 8.009826550848726967e+00 6.030390763814763133e+00
5.157143775218013815e+00 1.052749940283834817e+01 1.599094933836748833e+00
9.787922256488773343e+00 1.677935984755819376e+01 1.337584035126054438e+01
2.521657104383914216e+00 1.254899130385144801e+01 3.864098256185941604e+00
1.208656097083387060e+00 1.391123129787090207e+01 1.591501088437121858e+01
1.055480606855721426e+01 1.323199991241144957e+01 7.818182720033275857e+00
1.544403842232301827e+01 1.738902728097084349e+01 4.583533941105946141e+00
1.757202739302711336e+01 2.980108701569501584e+00 1.103273534184709526e+01
1.643728387023264759e+01 9.769287967799678007e-01 8.060568127292810914e+00
1.278840209801471239e+01 1.693639903025687232e+01 1.382404603068302329e+01
1.169669811206154009e+01 2.695099306954551466e-02 1.991042129157159479e+01
1.541744203695378523e+01 1.215086142188359197e+01 9.796895312052478211e-01
1.543379942320465936e+01 1.909963614070471749e+01 5.803827847132508211e+00
1.463271756041105576e+01 3.962347469628362262e+00 1.562629601219827258e+01
1.311019235754942169e+00 1.806350703350956621e+01 1.691897225244660774e+01
7.150391605065522072e+00 1.099888167710175679e+01 9.641049574578634207e+00
3.145362264789162143e+00 4.758730510699178851e+00 6.295283862467170977e+00
1.017060012717325712e+00 1.929594661322981253e+01 2.751458013493301991e-01
3.181873213230266284e+00 9.611839323944655789e+00 1.578636379598616735e+01
1.596304347011204428e+01 6.861583375044146393e+00 8.411436081657491570e+00
1.907578939273505014e+01 1.847632882492380446e+00 1.406320446220798281e+01
1.119728213760762614e+01 1.359780421190231969e+01 1.673256084377378627e+01
6.535523375795470713e+00 1.938031781879219206e+00 1.468663570434848076e+01
4.326896634117655793e+00 1.360523117160250095e+01 1.823387660933931187e+01
1.041626180439516602e+01 1.332330615971187271e+01 9.258117227669540483e+00
4.564590549733757996e+00 1.908114977586326333e+01 7.312327081900251713e+00
1.779454862026031492e+01 9.817898212787985912e+00 6.904928105639305258e+00
7.437721417815799541e+00 1.410557107763569640e+01 1.864438672650637940e+01
5.656565917369198004e+00 9.857522133432643940e+00 1.133238109297281859e+01
1.121355730341712231e+01 5.977450085190683104e-01 3.205617411419101614e+00
3.647214150558852896e+00 1.086846157058459816e+01 9.922834663035908775e+00
7.793463204232464747e+00 1.568222676824082562e+00 1.073223310005559661e+01
4.336842923939356886e+00 1.846081931805071363e+01 6.292163080181540913e+00
9.088426316570345165e+00 1.228361123709773750e+01 1.325651904055651897e+01
1.727857858344762221e+01 1.602937772276546724e+01 4.209985410461292155e+00
1.151942410175748499e+01 1.176001307649075400e+01 1.674440900361202633e+00
1.546399444772657006e+00 1.716802515638108417e+01 3.580105343112962579e+00
5.898703450689948724e+00 1.685608625726840160e+01 1.604633626058270934e+01
7.175481518346358278e+00 7.580601818862851005e+00 9.179994488728883795e+00
3.604382120423810232e+00 8.563465904527696182e+00 2.631529676708994625e+00
4.250283669469377301e+00 1.253027564972482644e+01 1.882070910096640404e+01
9.740525878953356198e+00 1.727570903050724738e+01 4.492743788056213106e+00
3.693192317030589145e+00 6.023030532782868818e+00 6.375073520935861104e+00
1.105803627605270378e+01 6.074927402221177886e+00 1.383583337766531152e+01
1.616569070477469694e+01 1.066719534745342912e+01 8.216518553597676799e+00
2.516905025276066077e+00 1.862618004395644888e+01 1.613017395783256092e+01
3.213196906933981634e+00 3.748722888848616819e+00 4.460985010115381399e-01
1.966014335741365926e+01 1.393809446621066073e+01 3.643561262848484805e+00
1.417102957996105594e+01 3.775794039091937560e+00 1.721168929006050519e+01
9.887642490799372297e+00 4.287989921782076053e-01 1.948145960504241359e+01
1.952226811650194094e+01 1.541244440198597410e+01 9.277429026528391631e+00
1.612333232399475946e+01 1.101106998443748886e+01 3.527639941285400926e+00
1.613864510998473989e+01 1.983189145666934827e+00 1.416171063278793874e+01
3.613618941471716806e+00 1.250754671063615753e+01 1.192691576606816461e+01
1.486017314112947574e+01 8.644803421197210014e+00 1.327183703369135870e+01
1.252619458099862904e+01 1.199893750175939910e+00 5.573239821230129287e+00
1.724319733856466286e+01 1.471600899572874965e+01 3.906241106931320761e+00
1.382366611819664826e+01 1.762058486150377234e+01 1.220997562373410972e+01
1.171178652411146537e+01 6.381320997427815556e+00 4.087364107663155721e+00
1.888821619247920225e+01 1.781782568894866259e+01 1.619415155824892238e+01
1.700812652296330185e+01 3.679640756773159271e+00 1.856671058339995284e+01
8.850029219321672969e+00 1.809486181237064173e+01 6.264059085325723863e-01
1.459103302255055468e+01 1.642738962536993341e+01 6.789633676739534884e+00
1.208304459237176331e+01 1.049002525105023054e+01 1.065752865159921470e+01
1.871463013746178206e+01 1.943360703768160391e+01 7.144703470916016208e+00
4.548650799665536759e+00 1.832255934794095875e+00 3.218582839415500274e+00
4.374273787707599226e+00 1.152343984508856600e+01 3.390100931205208834e+00
6.860549105292356842e+00 5.397034062245793962e+00 1.863558197710208830e+01
1.788705434578027464e+01 1.593883220930704070e+01 1.944189561195890903e+01
8.041121364293671192e+00 1.123577737912157204e+01 1.752625820573212589e+00
2.190924701186898194e+00 5.359200802521124629e-01 1.415401012565485495e+01
1.425140191352204866e+01 1.945315963481521493e+01 7.173271848286288943e-01
7.574508967392055148e+00 1.466943229300573748e+01 9.560737208799583442e+00
8.683210787959783516e+00 1.636071948760724126e+01 2.973351028185446943e+00
9.838778387330165742e+00 1.266870873029016131e+01 1.656151853945769048e+01
6.201487756087347591e+00 6.842761873310987397e+00 5.408515942382068786e+00
1.957943680058414770e+01 9.435128419604396299e+00 8.504901610318773564e+00
4.922790991485563872e+00 1.182431712748640074e+01 6.007083857379427627e+00
1.375511158803277212e+01 1.818450581307167013e+01 1.294061095284070539e+01
1.505834046826122474e+01 1.606426904659734234e+01 1.295832000371466286e+00
1.139900749432539229e+01 1.080550136192624144e+01 2.148159016403747845e+00
5.140949712612369993e+00 7.848217044362382211e-01 6.259986265313672860e+00
1.591732060288471828e+01 1.135438277855769940e+01 7.717596285640384579e+00
1.159804752691622909e+01 9.970702357772918134e+00 3.740838786608680078e+00
7.819214671258782445e+00 1.474692532406405121e+01 1.599763849799937798e+00
1.183588723443627089e+01 5.230688823295944090e+00 5.610299551811772645e+00
1.403170107465347272e+01 6.259679124919326654e+00 1.266197025549406874e+01
9.996432011853704225e+00 1.898762363866631020e+01 2.688494534691205917e+00
1.192060300806565287e+01 7.491809468189725152e+00 1.055981505824936484e+01
6.573122764497137283e+00 1.203749261839902118e+01 1.975730110879115031e+01
1.344417738990171962e+01 5.124552024014925600e+00 1.167506286585825848e+01
1.859305526500111938e+01 1.979887438072610095e+01 1.792313159572520576e+00
1.918181828964280555e+01 1.426409918571056679e+01 1.760341972177295844e+01
9.575270009881180044e+00 1.918761001448747194e+01 1.928623123099356107e+01
4.873725067145000978e+00 9.199366532460725665e+00 2.531704541945858367e+00
1.564035728411033510e+01 1.252302258707919158e+01 1.505181646230310832e+01
1.621279947189046666e+01 1.624441165448580549e+01 1.528540931575910911e+01
8.379356826536280778e+00 1.168773515224531678e+01 1.787758853742731091e+01
4.416571257052471111e+00 9.174278024635979634e+00 1.213313047383849153e+01
7.084761271774819846e-02 1.167908362701079561e+01 8.472104743632808521e+00
1.497432981447423472e+01 1.355864920688349962e+01 3.210673209057224309e+00
3.686163728562852349e+00 1.492188218167430414e+01 3.721359880901264905e+00
1.412776119573074851e+01 4.257434626616181106e+00 1.383721648358000600e+01
6.682963364623352831e+00 1.043920103302009572e+01 1.469523318137093248e+01
1.001565812571090852e+01 3.129634825270593002e+00 1.864825595533753955e+01
1.462851177669342206e+01 3.276207936145394406e+00 1.152881124395525880e+01
1.875632791702631508e+01 1.477048738222013569e+01 1.033559927088666797e+01
3.778486259716378193e+00 4.066352116376210191e+00 1.563465890382935264e+01
1.906494520449296104e+01 1.943870630782797448e+01 1.504521734889994455e+01
1.636378701606963659e+01 1.071675497675586897e+01 2.583844457303485775e-01
2.165653843046921878e+00 1.626574602399584535e+01 8.054189853689381451e+00
1.382119584208396290e+01 8.246367444278554615e+00 1.717841162611012962e+01
1.263974053490238170e+01 1.832827651411862391e+01 2.269823822614678299e+00
7.740323277769951105e+00 1.839637983935513432e+01 2.776148143075061192e+00
2.086282483451702419e+00 1.651797579436778562e+01 1.703858963968762197e+01
6.075145313221064214e-01 1.556135225865963356e+01 1.965173021111040441e+01
5.670268815950906927e+00 1.048903246785755172e+01 1.622772691770148157e+01
1.849266189042910469e+01 1.188751960699073962e+01 1.597109991753374558e+01
3.252231616064027442e+00 1.284032228348912241e+01 7.286524767690507609e+00
1.291105429487849676e+01 6.657892301134986646e+00 2.702236745246722194e-03
1.644175614139331287e+01 8.194118621433430505e+00 4.619907135767813422e+00
1.406422077532010029e+01 1.246805348694896765e+01 3.522920640951345828e+00
1.204883716664414806e+01 1.432721976382515372e+01 3.382605026347320631e+00
1.145086253670984000e+01 1.598479098208942872e+01 1.801224166181521724e+01
1.295741354267908640e+01 1.654270050442148143e+00 1.711636862466958675e+01
7.894732517944349937e+00 1.233371541027474816e+00 1.582319957908883978e+01
4.491767420267085420e+00 1.533607740451308032e+01 2.024289899479694288e+00
5.243001409788283773e-01 1.437133840905082849e+01 6.295622896450538519e+00
1.201605095631214049e+01 8.007914059888483038e+00 1.268201147991653066e+01
1.480438236357707105e+01 9.917860369484790439e+00 1.290935790405736761e+01
1.609105371737595647e+01 1.307582670188235952e+01 8.433258002441252899e+00
9.739545676964016963e+00 6.536081100275565881e+00 7.776810642826754716e+00
1.699295128448666503e+01 1.101459659569650640e+01 2.833516304861345425e+00
1.686645928129948757e+01 7.408879691236367471e+00 3.393384256342297789e+00
5.395859856051473002e+00 2.533105487131386280e+00 8.958995425975171401e+00
6.673056257009816861e+00 7.707233704030890920e+00 1.493996872606269477e+01
9.878975405044121416e-02 5.217253952423677177e+00 1.254033916023494122e+01
1.461128691102269350e+01 1.145788418144584497e+01 9.676687949951978673e-01
5.519197535996100967e-01 5.586111681257596828e+00 6.274444259195758988e+00
3.145191851978612618e+00 1.380797994761278247e+01 1.965478809996008547e+00
1.550090335191627133e+01 4.141003213289384810e-01 1.322327712828933421e+01
4.296475938928385752e+00 1.337600293195997025e+01 1.713530232802409259e+01
1.098049761285559356e+01 1.558834429003534261e+01 4.206763093004495246e-01
1.247038295107616790e+01 1.591462482593096439e+01 2.389730152102316207e+00
9.644263775972312658e+00 1.858831173933583258e+01 1.584190961914308815e+01
1.987601136368350652e+01 3.563648295642329877e-01 1.218486587117110531e+01
1.749899381701510137e+01 3.782779742567288217e+00 1.608574888702913341e+01
4.002362609191980169e+00 2.534559626380721298e+00 1.305270826154579566e+01
1.604909862808269239e+01 3.052993926770311006e+00 2.122357353580592854e+00
2.200985056849560362e+00 1.633372870391161769e+01 1.625632745905396348e+01
1.726133807466638359e+01 8.326886490832736243e+00 1.872394545778708164e+01
6.789331031250672055e+00 1.274398555534647670e+01 8.806489572560941781e+00
5.597456073012196498e+00 1.882279806375989040e+01 1.844155392401813742e+01
1.284826126024520843e+00 1.303992105498211274e+01 4.082339996386066261e+00
6.294339505372774646e+00 2.182714931405311809e+00 1.677847820224420161e+01
2.451324015378193444e+00 7.979195079683267799e+00 1.594253083674845151e+01
9.908172339326572597e+00 5.565212882370551561e+00 1.961760245430402705e+01
1.786011996445774841e+01 6.028290166295851016e+00 1.559253349100063701e+01
1.789484663831146349e+01 1.727321569440691107e+01 1.519559567849651138e+01
4.191447956158955712e-01 8.018907057028917151e+00 1.807513658920371213e+01
1.642091737348711078e+01 1.053080083393326305e+01 1.703342872249784534e+00
5.922847010751340235e+00 1.918857360095121578e+01 1.653362773606214375e+01
9.536818433797416006e+00 1.333892336395517830e+01 7.512399838910754468e+00
1.680341679903424890e+01 7.529930575159847095e+00 1.995409275699890728e+01
4.149754194465000268e+00 7.194338574685308352e+00 5.945781266663825448e+00
1.185102601175690928e+01 1.877389159368834726e+01 1.000037132070595902e+01
1.749233864963638396e+01 5.485012275413119198e+00 5.085538213868856516e-01
4.833907077491808479e+00 1.155565956254654303e+01 3.911124806017960420e+00
4.159825480777978868e-01 1.184411847099284998e+01 6.661269787834713796e+00
1.621267998846927583e+01 3.816351978493306429e+00 1.827854224177770348e+01
4.757870608430780912e-01 1.206571547830803670e+00 6.703627844375183997e+00
1.533390341419469571e+01 1.240473277888194836e+01 1.227383933249375225e+01
1.904768384223537225e+01 1.569521779573192077e+01 3.111574311669593751e+00
7.299771494006224160e-01 1.232175875934019871e+01 5.805957984230789570e+00
1.398811330220331151e+01 1.404192146256841767e+01 6.506161523542928649e+00
1.526238902524621643e+01 1.290589864564890732e+01 1.604116950580530698e+01
2.288454067766856337e+00 1.818224525289178928e+01 4.150189207217820808e+00
6.667093317345251791e+00 1.484316707165325155e+01 9.989457536206741040e+00
1.630826253146604454e+00 1.716781061523413854e+01 9.447789472906769959e+00
6.195985191129147474e+00 1.322995079792692330e+01 1.459017586795142307e+01
3.477449830177639090e+00 3.478218691316947719e+00 1.426305218767998717e+01
1.930223244723388021e+01 2.035225319390843168e+00 7.668663419862808972e-01
1.559991400958519492e+01 1.156841870072824641e+01 2.970472413292839509e+00
1.908496799876049010e+01 8.946840101451552840e+00 8.068000786817538739e+00
1.272188666045470917e+01 1.101470247198306218e+01 4.623091744216756283e+00
4.096968285998925374e+00 1.698743840488505086e+01 3.503630134729796008e+00
8.624185051487650355e+00 1.130897493659643871e-01 1.066059498060944399e+01
7.371316569484371861e+00 1.073058695585202926e+01 1.498934628051900653e+01
5.382317744155873385e+00 9.164719027817831432e+00 1.995190573388237709e+01
9.188511126192429046e+00 5.823244027332073358e+00 1.079127523405837685e+01
1.178503550620674645e+01 6.131118327784750299e+00 9.056675554485074997e+00
4.615843839105789215e+00 1.813512218088646488e+00 8.790551174139071833e-02
4.458993825633552177e+00 1.034114235853223640e+01 1.957059132159483283e+01
5.711694920743766168e+00 4.628644087549844244e+00 1.398890524212875519e+01
1.026184143257599501e+01 9.283642659267927755e+00 7.054226875109705475e+00
1.171702620882371804e+01 2.276401260494658185e+00 1.810096945087098419e+01
7.313860623077066414e+00 1.136536331465612548e+01 4.333435179470209064e+00
1.161897344231803864e+01 1.519603298759204435e+01 3.382249439070197372e+00
1.156896961606582153e+01 8.965684216027224096e+00 1.827736372208834581e+01
1.167370077791304439e+01 1.051009905261099675e+01 8.807622667959403628e+00
7.606672067973761386e+00 1.265725343657043211e+01 3.979380870847071261e+00
8.892335704184702649e+00 5.104533485456581587e+00 2.453042442841448789e+00
1.672416470286611556e+01 7.264109419389366273e+00 2.155153508137672702e+00
1.300313228937827859e+01 5.226936495928229398e+00 7.224382256964512194e+00
1.526044201138132372e+01 7.523756127501057378e+00 1.253176100436383944e+00
1.503140840626626584e+01 2.281189329522508302e+00 1.566487391834201759e+01
1.725284145996638685e+01 1.729211335181137787e+00 1.102979596681249852e+01
6.760573654131977861e+00 1.293522977412822073e+01 7.663858940594245439e+00
1.621239281076305971e+01 1.713461380275644075e+01 2.214231878610852799e+00
5.697872101122003485e-01 7.850066496542680738e+00 1.270615299210492033e+01
1.385054130936484285e+01 4.748277023904399208e+00 5.879718410124894490e+00
1.947694008113235498e+01 9.848355111751477153e+00 1.037328615483129113e+01
1.096948629944599851e+01 2.456217269844021001e+00 2.394115664983085257e+00
7.238305823828159902e+00 1.929903683394915692e+00 1.077155935473052040e+01
1.166342212104290077e+01 1.407302366702650787e+01 7.389359595792202029e+00
5.092331392100135190e+00 1.952424587496226493e+01 1.923154861807247329e+01
9.366425180191411570e+00 1.615164355335112845e+01 1.911607876834722219e+01
3.315452033839434431e+00 1.473878392392635561e+01 1.812393203845930145e+01
1.900480110888890906e+01 7.074939479502684314e+00 1.223207699986319952e+01
1.112062688587192039e+00 1.350382159552646399e+01 1.788244529975350972e+01
1.208543849838354234e+01 1.917538945159808250e+01 2.065915709806946943e+00
9.560953012010797281e+00 7.998329148897871832e+00 7.450779254324775280e+00
1.083078929122250500e+01 6.773165572482309216e-02 4.050027038129977441e+00
1.520035442423352201e+01 2.823423909758457917e+00 1.660600955629038467e+01
7.457920666287114031e+00 1.589040576353029088e+01 8.032331946108683951e+00
6.218894466525117792e+00 1.902130382216541449e+01 1.903258778343918678e+01
1.185352537044184373e+01 2.131383638373813838e+00 1.705043604149870617e+01
7.907877423623723701e+00 2.355629672603098967e-01 1.768134544461914004e+01
9.180573509978016133e+00 1.056672002888144668e+01 1.615673900391059448e+01
9.758028775418940981e+00 1.897910123530181536e+01 3.232856588131445275e+00
9.529556848976369565e+00 5.533775831517859345e+00 1.005084109974087703e+01
1.888463110367292686e+01 1.628556394454996692e+01 1.803847414451672293e+01
6.238548605600606756e+00 1.732812294171118594e+01 1.390202976478331909e+00
1.703968926632573044e+00 5.406395310068439208e+00 7.791335319289476757e+00
5.318212580207630324e+00 1.959798809581453582e+01 1.623296528008798489e+01
1.030974595776082481e+01 3.345118938173370360e+00 1.876368434371347149e+00
6.726524294819553162e+00 1.710588877924181261e+01 3.587736141512962718e+00
4.672488187259027370e-01 1.330291362956780787e+01 1.617364891779558178e+01
1.823219818623249111e+01 6.472742230971137545e+00 2.236150071190512456e+00
1.338693627635828776e+01 1.345878971051173156e+00 1.234222008199759202e+01
1.947357313764631215e+01 4.948194784649516009e+00 3.232011569633046477e+00
7.125923551155466384e-01 1.162057903825268568e+01 2.197833471575816322e+00
5.765106905171362150e+00 2.344048545419590290e+00 1.009968143530120166e+01
1.702982426306183328e+01 6.072617927489281087e+00 1.344762850718726988e+01
1.953240139022602762e+01 6.355949293339637762e+00 5.018385607396793624e+00
1.245453658156568366e+01 3.228496948834442826e+00 9.246036724649957961e+00
1.874147434626193842e+01 1.288676320323209890e+00 1.084079430093304275e+01
3.727357917213356409e+00 7.890228993383980871e+00 1.887220717306394491e+01
7.321766473129764918e+00 1.447808428304094619e+01 3.497951812525244009e+00
7.865472306426224414e-01 1.175719643886713506e+01 1.517215169072994208e+01
1.813191620677649940e+00 1.420133879752781070e+01 2.130895718843435205e+00
4.951315637270963244e+00 1.270505767711414435e+01 1.503488503815237820e+01
1.968294423243641589e+01 8.070047226443377575e-01 1.569001423339116563e+01
1.566772562717018147e+01 6.727297417688427039e-01 1.644914485797703207e+01
1.872485440230902043e+01 5.507240933872878941e+00 1.370100621573215349e+01
1.161658169643436622e+01 8.393490346501719657e-01 5.285574035829931461e+00
2.735011847701449561e+00 9.285366500191569727e+00 3.563701714451903424e+00
1.595953564417274961e+01 4.982221599869225415e+00 9.866112924259621053e+00
8.236398937164604916e+00 1.225732164001447977e+01 1.073870054404710839e+00
1.673230079118365410e+01 1.707138375264127816e+00 1.186127956517092308e+01
1.133157190343477083e+01 6.748273940311153574e+00 1.735995103299493536e+01
1.123515339456784190e+01 3.358381384637665334e-01 4.451683544357505795e+00
4.637006946943795782e+00 1.053287841873853736e+01 1.890483257401232109e+01
8.316004687068261347e+00 1.015832009654642221e+01 9.997387582388792993e+00
1.155453639981338654e+00 1.119721790588304700e+01 1.015230889964557903e-01
1.781615549334679827e+01 6.965508912105047301e+00 9.275527641105492549e+00
1.080805166589502875e+01 1.058441706397450766e+01 8.304530153792464731e+00
1.501464242245880598e+01 1.007311195971194273e+01 1.119912106401039775e+00
1.419991605053262163e+01 1.775556257620701572e+01 1.286914118149863029e+01
4.276250814546123102e+00 1.627273238694782975e+01 1.338681989661718585e+00
1.854777207937439343e+01 1.925368107712428412e+01 1.792532699376242533e+01
1.223975203151594915e+00 1.509404945054635405e+00 9.011561993634691348e+00
2.194943875494046459e+00 1.877640576542419737e+01 9.705328951731640785e+00
5.216227060122527348e+00 1.486080380567930703e+01 1.688088707737033545e+01
8.212309994247849687e+00 1.287695777335280667e+01 6.582461320676003425e+00
1.251332426363976680e+01 1.395650004243432818e+01 1.712772751477069022e+01
1.276409123650453559e+01 9.267092658243401004e+00 1.117797986214748462e+01
6.113755862360137350e+00 4.078389036120655398e-01 4.247092748998301914e+00
1.242978019441407156e+01 1.326256672807486048e+01 1.169215723874168944e+01
9.042791400820338055e+00 4.356422507665575594e+00 6.449680275291116871e+00
1.099070526776987222e+01 7.993714068441293286e+00 6.805343837475432700e+00
1.589809313358526488e+01 7.419391712002973982e+00 2.740259265591922233e+00
1.403169879421245447e+01 8.972298589150595305e+00 8.031561069202838965e+00
1.344854912285986792e+01 1.337958972288058135e+01 1.745916261721233198e+01
1.929539884202950262e+01 5.880174791470014206e+00 6.137258875502791255e+00
1.481705838303003375e+01 1.859961296602317127e+00 1.070658133606104556e+01
1.747793859169851771e+01 1.802051836868805523e+01 3.046717052164580330e+00
8.166905019318310366e+00 7.107870788028291642e+00 3.094851515184415813e+00
2.911619256801625255e+00 1.748206940353133731e+01 1.667585480238292206e+01
9.841559860317685704e+00 4.407466849272965170e+00 4.015772381534159052e+00
8.534641826538900133e+00 7.001124362981947513e-03 1.083859604224306494e+01
1.179950797532610274e+01 7.733610875448507116e+00 3.715706806233098902e+00
1.031814661896667751e+01 1.504274861610019087e+01 8.641941454933366629e+00
8.487179669857354725e+00 2.779087131125754784e+00 1.771467618721269588e+01
1.817694302030358244e+01 1.368297972028561738e+01 1.459998433404883755e+01
1.523108736137761809e+01 1.914783005037991259e+01 1.271830247155052618e+01
1.676565671838121574e+01 1.137791753488018465e+01 1.149470907912997220e-01
4.713786035411972719e+00 1.502631037162001082e+00 1.781772595544743876e+00
1.496393319315559722e+01 1.907873074229365074e+01 9.267905967302020542e+00
3.420605072294604643e+00 1.946138302136004228e+01 1.554521849678835999e+01
5.059396683262584737e+00 1.179111420257443754e+01 1.694676442505262770e+01
1.034505581870227076e+01 1.047187126578368321e+01 1.513149670850754447e+01
5.766942056752768053e+00 9.343210501397724244e+00 7.221102977105225307e+00
1.859697591806138917e+01 2.139702594727030949e+00 1.926127136089304059e+01
1.758948960639834524e+01 8.716762176901982073e-01 1.994986541769686283e+01
9.791295245901567412e+00 7.476110121679734988e+00 1.719085699736939787e+01
5.237537596468319734e+00 1.114089935463471903e+01 2.413380171945278541e-01
2.947436828339069503e+00 4.101110664128533756e+00 7.002043945902258315e+00
1.432988014056335402e+01 1.921806287020138271e+01 3.934162851604201538e+00
1.154680355918341128e+01 3.270462716187296781e+00 1.951812629060489712e+01
4.846414428820264853e+00 1.087854833453446446e+01 7.050825938023857375e+00
1.212728915878844838e+01 9.183471762621229217e+00 3.013213345939611543e+00
8.256669085996799495e+00 2.639001216334704303e+00 1.964759874273240570e+01
5.578484768657336446e+00 5.250179413170503295e+00 5.202500442144026849e+00
3.628998141330312954e+00 5.424377769518844872e+00 1.685060418027228124e+01
1.675031186141319139e+01 2.738918567328580966e+00 4.941011954075147372e-02
7.800310628073869879e+00 1.009891445538099575e+01 1.608966164039095403e+01
1.734493420202401737e+01 1.961265577635142066e+01 1.363516602159723412e+01
6.630357772344597223e+00 3.263443697742078875e+00 4.760493490091127100e+00
8.559244331671484574e+00 9.128667856159415450e+00 1.602140160223702026e+01
2.679660835208141911e+00 5.250720355636633307e+00 1.341015448463152016e+00
6.038176216444695044e+00 6.157414420908484232e+00 9.374931235597072643e-02
1.794816126171273751e+01 1.999483621425812885e+01 2.852150505341961573e+00
4.257699770888214275e+00 2.783806729552100734e+00 1.521570975371453471e+01
3.969559505967985569e+00 1.893103996289927693e+01 1.641274544068463115e+01
6.616243527471028507e+00 1.860375951597969291e+01 1.485333710912073002e+01
1.906086127887575543e+01 4.294385959116819862e+00 1.636802165767673856e+01
5.054294927043523344e-01 5.301083774118295899e+00 7.469612132794376080e+00
1.012444166562786663e+01 1.499032934216872892e+00 1.588788869604790754e+01
1.160365945173237812e+01 1.965862931449910889e+01 5.536512241907749043e+00
6.199859829526282340e+00 9.282932393847531216e+00 8.980042228686686556e+00
1.303740662239231085e+01 9.235762115290071961e+00 1.420503195579945910e+01
1.104860128164335720e+01 5.683132552589700737e+00 1.026562735946136051e+01
1.516211207459281773e+01 1.315441397925543754e+01 1.735247591384620947e+01

File diff suppressed because it is too large Load Diff

271
jastlib.py Normal file
View File

@ -0,0 +1,271 @@
import numpy as np
import matplotlib.pyplot as plt
def get_sphere_distribution(n, dmin, Ls, maxiter=1e4, allow_wall=True):
"""Get random points in a box with given dimensions and minimum separation.
Parameters:
- n: number of points
- dmin: minimum distance
- Ls: dimensions of box, shape (3,) array
- maxiter: maximum number of iterations.
- allow_wall: whether to allow points on wall;
(if False: points need to keep distance dmin/2 from the walls.)
Return:
- ps: array (n, 3) of point positions,
with 0 <= ps[:, i] < Ls[i]
- n_iter: number of iterations
- dratio: average nearest-neighbor distance, divided by dmin.
Note: with a fill density (sphere volume divided by box volume) above about
0.53, it takes very long. (Random close-packed spheres have a fill density
of 0.64).
Author: Han-Kwang Nienhuys (2020)
Copying: BSD, GPL, LGPL, CC-BY, CC-BY-SA
See Stackoverflow: https://stackoverflow.com/a/62895898/6228891
"""
Ls = np.array(Ls).reshape(3)
if not allow_wall:
Ls -= dmin
# filling factor; 0.64 is for random close-packed spheres
# This is an estimate because close packing is complicated near the walls.
# It doesn't work well for small L/dmin ratios.
sphere_vol = np.pi/6*dmin**3
box_vol = np.prod(Ls + 0.5*dmin)
fill_dens = n*sphere_vol/box_vol
if fill_dens > 0.64:
msg = f'Too many to fit in the volume, density {fill_dens:.3g}>0.64'
raise ValueError(msg)
# initial try
ps = np.random.uniform(size=(n, 3)) * Ls
# distance-squared matrix (diagonal is self-distance, don't count)
dsq = ((ps - ps.reshape(n, 1, 3))**2).sum(axis=2)
dsq[np.arange(n), np.arange(n)] = np.infty
for iter_no in range(int(maxiter)):
# find points that have too close neighbors
close_counts = np.sum(dsq < dmin**2, axis=1) # shape (n,)
n_close = np.count_nonzero(close_counts)
if n_close == 0:
break
# Move the one with the largest number of too-close neighbors
imv = np.argmax(close_counts)
# new positions
newp = np.random.uniform(size=3)*Ls
ps[imv]= newp
# update distance matrix
new_dsq_row = ((ps - newp.reshape(1, 3))**2).sum(axis=-1)
dsq[imv, :] = dsq[:, imv] = new_dsq_row
dsq[imv, imv] = np.inf
else:
raise RuntimeError(f'Failed after {iter_no+1} iterations.')
if not allow_wall:
ps += dmin/2
dratio = (np.sqrt(dsq.min(axis=1))/dmin).mean()
return ps, iter_no+1, dratio
def generateElsAroundPoints(n,LS,dmin):
"""
Parameters:
- n: number of points
- LS: list of position of all atoms
- dmin: minimum intra block distance
- shift: inter block distance
Return:
- r: array (n, 3) of point positions,
"""
xs = None
for Ls in LS:
# Get list of random points around Ls
distrib,a,b = get_sphere_distribution(n,dmin,Ls)
if xs is None:
xs = distrib[:,0]
ys = distrib[:,1]
zs = distrib[:,2]
else:
xs = np.concatenate((xs,distrib[:,0]))
ys = np.concatenate((ys,distrib[:,1]))
zs = np.concatenate((zs,distrib[:,2]))
return((np.array((xs,ys,zs))).T)
def getCoefList(Nord,Natom):
assert(Nord < 11)
dict = {
0 : lambda x,y:x-y-2,
1 : lambda x,y:x-y,
2 : lambda x,y:x-y,
3 : lambda x,y:x-y,
4 : lambda x,y:x-y,
5 : lambda x,y:x-y,
6 : lambda x,y:x-y,
7 : lambda x,y:x-y,
8 : lambda x,y:x-y,
9 : lambda x,y:x-y,
10 : lambda x,y:x-y,
11 : lambda x,y:x-y,
}
count = 0
for p in range(2,Nord+1):
for k in range(p-1,-1,-1):
lmax = dict[k](p,k)
for l in range(lmax,-1,-1):
if (p-k-l) & 1 is 0:
count += 1
coeflista = np.random.rand(Nord+1,Natom)
coeflistb = np.random.rand(Nord+1)
coeflistc = np.random.rand(count,Natom)
return (coeflista.reshape((Nord+1)*Natom),coeflistb,coeflistc.reshape(count*Natom))
#return (coeflista,coeflistb,coeflistc)
def get_sphere_distribution(n, dmin, Ls, maxiter=1e4, allow_wall=True):
"""Get random points in a box with given dimensions and minimum separation.
Parameters:
- n: number of points
- dmin: minimum distance
- Ls: dimensions of box, shape (3,) array
- maxiter: maximum number of iterations.
- allow_wall: whether to allow points on wall;
(if False: points need to keep distance dmin/2 from the walls.)
Return:
- ps: array (n, 3) of point positions,
with 0 <= ps[:, i] < Ls[i]
- n_iter: number of iterations
- dratio: average nearest-neighbor distance, divided by dmin.
Note: with a fill density (sphere volume divided by box volume) above about
0.53, it takes very long. (Random close-packed spheres have a fill density
of 0.64).
Author: Han-Kwang Nienhuys (2020)
Copying: BSD, GPL, LGPL, CC-BY, CC-BY-SA
See Stackoverflow: https://stackoverflow.com/a/62895898/6228891
"""
Ls = np.array(Ls).reshape(3)
if not allow_wall:
Ls -= dmin
# filling factor; 0.64 is for random close-packed spheres
# This is an estimate because close packing is complicated near the walls.
# It doesn't work well for small L/dmin ratios.
sphere_vol = np.pi/6*dmin**3
box_vol = np.prod(Ls + 0.5*dmin)
fill_dens = n*sphere_vol/box_vol
if fill_dens > 0.64:
msg = f'Too many to fit in the volume, density {fill_dens:.3g}>0.64'
raise ValueError(msg)
# initial try
ps = np.random.uniform(size=(n, 3)) * Ls
# distance-squared matrix (diagonal is self-distance, don't count)
dsq = ((ps - ps.reshape(n, 1, 3))**2).sum(axis=2)
dsq[np.arange(n), np.arange(n)] = np.infty
for iter_no in range(int(maxiter)):
# find points that have too close neighbors
close_counts = np.sum(dsq < dmin**2, axis=1) # shape (n,)
n_close = np.count_nonzero(close_counts)
if n_close == 0:
break
# Move the one with the largest number of too-close neighbors
imv = np.argmax(close_counts)
# new positions
newp = np.random.uniform(size=3)*Ls
ps[imv]= newp
# update distance matrix
new_dsq_row = ((ps - newp.reshape(1, 3))**2).sum(axis=-1)
dsq[imv, :] = dsq[:, imv] = new_dsq_row
dsq[imv, imv] = np.inf
else:
raise RuntimeError(f'Failed after {iter_no+1} iterations.')
if not allow_wall:
ps += dmin/2
dratio = (np.sqrt(dsq.min(axis=1))/dmin).mean()
return ps, iter_no+1, dratio
def scalingee(r,kappa=1.0):
return (numpy.ones_like(r) - numpy.exp(-kappa*r))/kappa
def scalingen(r,kappa=1.0):
return numpy.exp(-kappa*r)
if False:
Nord = 5
L1 = 2.0
n = 2 # number of points
dmin = 0.1 # min dist between points
Ls = np.array([L1,L1,L1]) # lengths of the box
shift = -10.0
kappa = 2.0
filename_atom = str(n) + "_geometry.txt"
filename_elec = str(n)
filename_coeffs = str(n) + "_jast_coeffs.txt"
(coeffsa, coeffsb, coeffsc) = getCoefList(Nord,n)
coeffsall = np.concatenate((coeffsa,coeffsb,coeffsc))
print(coeffsa.shape,coeffsb.shape,coeffsc.shape)
atomList,_,_ = get_sphere_distribution(n, dmin, Ls, maxiter=1e4, allow_wall=True)
#print(atomList)
L1 = 5.0
n = 5 # number of points
dmin = 0.1 # min dist between points
Ls = np.array([L1,L1,L1]) # lengths of the box
shift = -10.0
kappa = 2.0
filename_elec = filename_elec + "_" + str(n) + "_elec_coord.txt"
#rlist = generateBlockRandomPointsAtShftApart(n,L1,dmin,shift)
rlist = generateElsAroundPoints(n,atomList,dmin)
# Save file
np.savetxt(filename_elec,rlist)
np.savetxt(filename_atom,atomList)
np.savetxt(filename_coeffs,coeffsall)
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
xs = rlist.T[0]
ys = rlist.T[1]
zs = rlist.T[2]
ax.scatter(xs, ys, zs, marker='o')
plt.show()
rijScaled = np.array([[(lambda xval, yval: np.linalg.norm(xval-yval))(xval,yval) for yval in rlist] for xval in rlist])
plt.imshow(rijScaled)
plt.colorbar()
plt.show()

View File

@ -3,7 +3,9 @@ BEGIN_PROVIDER [ integer, nnuc ]
BEGIN_DOC
! Number of nuclei
END_DOC
nnuc = 2
!nnuc = 2
! read(*,*)nnuc
nnuc = nelec/5
END_PROVIDER
@ -14,7 +16,8 @@ BEGIN_PROVIDER [ integer, typenuc ]
! Type of the nuclei
END_DOC
typenuc = 1
typenuc_arr = (/1, 1/)
!typenuc_arr = (/1, 1/)
typenuc_arr = 1
END_PROVIDER

View File

@ -186,7 +186,7 @@ BEGIN_PROVIDER [double precision, rescale_een_n, (nelec, nnuc, 0:ncord)]
END_PROVIDER
BEGIN_PROVIDER [double precision, rescale_een_n_deriv_e, (4, nelec, nnuc, 0:ncord)]
BEGIN_PROVIDER [double precision, rescale_een_n_deriv_e, (nelec, 4, nnuc, 0:ncord)]
implicit none
BEGIN_DOC
! Derivative of the scaled distance J_{een} wrt R_{ia}
@ -198,23 +198,23 @@ BEGIN_PROVIDER [double precision, rescale_een_n_deriv_e, (4, nelec, nnuc, 0:ncor
kappa_l = - dble(l) * kappa
do a = 1, nnuc
do i = 1, nelec
! r'(x) \lor r''(x)
do ii = 1, 4
rescale_een_n_deriv_e(ii, i, a, l) = &
! r'(x) \lor r''(x)
rescale_een_n_deriv_e(i, ii, a, l) = &
kappa_l * elnuc_dist_deriv_e(ii, i, a)
!print *, "pp", ii, i, a, elnuc_dist_deriv_e(ii, i, a)
enddo
! \left(r''(x)+r'(x)^2\right)
rescale_een_n_deriv_e(4, i, a, l) = rescale_een_n_deriv_e(4, i, a, l) + &
rescale_een_n_deriv_e(1, i, a, l) * rescale_een_n_deriv_e(1, i, a, l) + &
rescale_een_n_deriv_e(2, i, a, l) * rescale_een_n_deriv_e(2, i, a, l) + &
rescale_een_n_deriv_e(3, i, a, l) * rescale_een_n_deriv_e(3, i, a, l)
rescale_een_n_deriv_e(i, 4, a, l) = rescale_een_n_deriv_e(i, 4, a, l) + &
rescale_een_n_deriv_e(i, 1, a, l) * rescale_een_n_deriv_e(i, 1, a, l) + &
rescale_een_n_deriv_e(i, 2, a, l) * rescale_een_n_deriv_e(i, 2, a, l) + &
rescale_een_n_deriv_e(i, 3, a, l) * rescale_een_n_deriv_e(i, 3, a, l)
! \times e^{r(x)}
do ii = 1, 4
rescale_een_n_deriv_e(ii, i, a, l) = &
rescale_een_n_deriv_e(ii, i, a, l) * rescale_een_n(i, a, l)
rescale_een_n_deriv_e(i, ii, a, l) = &
rescale_een_n_deriv_e(i, ii, a, l) * rescale_een_n(i, a, l)
enddo
enddo
enddo
@ -243,7 +243,7 @@ BEGIN_PROVIDER [double precision, elnuc_dist_deriv_e, (4, nelec, nnuc)]
end do
END_PROVIDER
BEGIN_PROVIDER [double precision, rescale_een_e_deriv_e, (4, nelec, nelec, 0:ncord)]
BEGIN_PROVIDER [double precision, rescale_een_e_deriv_e, (nelec, 4, nelec, 0:ncord)]
BEGIN_DOC
! Derivative of the scaled distance J_{een} wrt R_{ia}
END_DOC
@ -259,27 +259,27 @@ BEGIN_PROVIDER [double precision, rescale_een_e_deriv_e, (4, nelec, nelec, 0:nco
do i = 1, nelec
! r'(x) \lor r''(x)
do ii = 1, 4
rescale_een_e_deriv_e(ii, i, j, l) = &
rescale_een_e_deriv_e(i, ii, j, l) = &
kappa_l * elec_dist_deriv_e(ii, i, j)
enddo
! \left(r''(x)+r'(x)^2\right)
rescale_een_e_deriv_e(4, i, j, l) = rescale_een_e_deriv_e(4, i, j, l) + &
rescale_een_e_deriv_e(1, i, j, l) * rescale_een_e_deriv_e(1, i, j, l) + &
rescale_een_e_deriv_e(2, i, j, l) * rescale_een_e_deriv_e(2, i, j, l) + &
rescale_een_e_deriv_e(3, i, j, l) * rescale_een_e_deriv_e(3, i, j, l)
rescale_een_e_deriv_e(i, 4, j, l) = rescale_een_e_deriv_e(i, 4, j, l) + &
rescale_een_e_deriv_e(i, 1, j, l) * rescale_een_e_deriv_e(i, 1, j, l) + &
rescale_een_e_deriv_e(i, 2, j, l) * rescale_een_e_deriv_e(i, 2, j, l) + &
rescale_een_e_deriv_e(i, 3, j, l) * rescale_een_e_deriv_e(i, 3, j, l)
! \times e^{r(x)}
do ii = 1, 4
rescale_een_e_deriv_e(ii, i, j, l) = &
rescale_een_e_deriv_e(ii, i, j, l) * rescale_een_e(i, j, l)
rescale_een_e_deriv_e(i, ii, j, l) = &
rescale_een_e_deriv_e(i, ii, j, l) * rescale_een_e(i, j, l)
enddo
enddo
enddo
enddo
END_PROVIDER
BEGIN_PROVIDER [double precision, rescale_een_e_deriv_e_t, (4, nelec, nelec, 0:ncord)]
BEGIN_PROVIDER [double precision, rescale_een_e_deriv_e_t, (nelec, 4, nelec, 0:ncord)]
implicit none
BEGIN_DOC
! Transposed rescale_een_e_deriv_e
@ -288,7 +288,7 @@ BEGIN_PROVIDER [double precision, rescale_een_e_deriv_e_t, (4, nelec, nelec, 0:n
do l=0,ncord
do j=1,nelec
do i=1,nelec
rescale_een_e_deriv_e_t(1:4,j,i,l) = rescale_een_e_deriv_e(1:4,i,j,l)
rescale_een_e_deriv_e_t(j,1:4,i,l) = rescale_een_e_deriv_e(i,1:4,j,l)
enddo
enddo
enddo