4
1
mirror of https://github.com/pfloos/quack synced 2024-07-04 18:36:03 +02:00

debug numgrid

This commit is contained in:
Pierre-Francois Loos 2020-03-25 12:56:28 +01:00
parent 1e65d1cf1f
commit ddb940de48
9 changed files with 85 additions and 57 deletions

View File

@ -1,9 +1,30 @@
1 3 1 6
S 3 S 8
1 38.3600000 0.0238090 1 2940.0000000 0.0006800
2 5.7700000 0.1548910 2 441.2000000 0.0052360
3 1.2400000 0.4699870 3 100.5000000 0.0266060
4 28.4300000 0.0999930
5 9.1690000 0.2697020
6 3.1960000 0.4514690
7 1.1590000 0.2950740
8 0.1811000 0.0125870
S 8
1 2940.0000000 -0.0001230
2 441.2000000 -0.0009660
3 100.5000000 -0.0048310
4 28.4300000 -0.0193140
5 9.1690000 -0.0532800
6 3.1960000 -0.1207230
7 1.1590000 -0.1334350
8 0.1811000 0.5307670
S 1 S 1
1 0.2976000 1.0000000 1 0.0589000 1.0000000
P 3
1 3.6190000 0.0291110
2 0.7110000 0.1693650
3 0.1951000 0.5134580
P 1 P 1
1 1.2750000 1.0000000 1 0.0601800 1.0000000
D 1
1 0.2380000 1.0000000

View File

@ -13,7 +13,7 @@
# GGA = 2: # GGA = 2:
# Hybrid = 4: # Hybrid = 4:
# Hartree-Fock = 666 # Hartree-Fock = 666
1 RMFL20 1 RVWN5
# quadrature grid SG-n # quadrature grid SG-n
1 1
# Number of states in ensemble (nEns) # Number of states in ensemble (nEns)
@ -21,4 +21,4 @@
# Ensemble weights: wEns(1),...,wEns(nEns-1) # Ensemble weights: wEns(1),...,wEns(nEns-1)
0.50000 0.00000 0.50000 0.00000
# GOK-DFT: maxSCF thresh DIIS n_diis guess_type ortho_type # GOK-DFT: maxSCF thresh DIIS n_diis guess_type ortho_type
32 0.0000001 T 5 1 1 32 0.00001 T 5 1 1

View File

@ -1,4 +1,4 @@
# nAt nEla nElb nCore nRyd # nAt nEla nElb nCore nRyd
1 1 1 0 0 1 2 2 0 0
# Znuc x y z # Znuc x y z
He 0.0 0.0 0.0 Be 0.0 0.0 0.0

View File

@ -1,3 +1,3 @@
1 1
He 0.0000000000 0.0000000000 0.0000000000 Be 0.0000000000 0.0000000000 0.0000000000

View File

@ -1,9 +1,30 @@
1 3 1 6
S 3 S 8
1 38.3600000 0.0238090 1 2940.0000000 0.0006800
2 5.7700000 0.1548910 2 441.2000000 0.0052360
3 1.2400000 0.4699870 3 100.5000000 0.0266060
4 28.4300000 0.0999930
5 9.1690000 0.2697020
6 3.1960000 0.4514690
7 1.1590000 0.2950740
8 0.1811000 0.0125870
S 8
1 2940.0000000 -0.0001230
2 441.2000000 -0.0009660
3 100.5000000 -0.0048310
4 28.4300000 -0.0193140
5 9.1690000 -0.0532800
6 3.1960000 -0.1207230
7 1.1590000 -0.1334350
8 0.1811000 0.5307670
S 1 S 1
1 0.2976000 1.0000000 1 0.0589000 1.0000000
P 3
1 3.6190000 0.0291110
2 0.7110000 0.1693650
3 0.1951000 0.5134580
P 1 P 1
1 1.2750000 1.0000000 1 0.0601800 1.0000000
D 1
1 0.2380000 1.0000000

View File

@ -53,13 +53,9 @@ subroutine allocate_grid(nNuc,ZNuc,rNuc,nShell,TotAngMomShell,ExpShell,max_ang_m
do iNuc=1,nNuc do iNuc=1,nNuc
context = numgrid_new_atom_grid(radial_precision, & context = numgrid_new_atom_grid(radial_precision,min_num_angular_points,max_num_angular_points, &
min_num_angular_points, & int(ZNuc(iNuc)),max_exponent(iNuc),max_ang_mom(iNuc), &
max_num_angular_points, & min_exponent(iNuc,1:max_ang_mom(iNuc)+1))
int(ZNuc(iNuc)), &
max_exponent(iNuc), &
max_ang_mom(iNuc), &
min_exponent(iNuc,:) )
nGrid = nGrid + numgrid_get_num_grid_points(context) nGrid = nGrid + numgrid_get_num_grid_points(context)

View File

@ -58,28 +58,16 @@ subroutine build_grid(nNuc,ZNuc,rNuc,nShell,TotAngMomShell,ExpShell,max_ang_mom,
do iNuc=1,nNuc do iNuc=1,nNuc
context = numgrid_new_atom_grid(radial_precision, & context = numgrid_new_atom_grid(radial_precision,min_num_angular_points,max_num_angular_points, &
min_num_angular_points, & int(ZNuc(iNuc)),max_exponent(iNuc),max_ang_mom(iNuc), &
max_num_angular_points, & min_exponent(iNuc,1:max_ang_mom(iNuc)+1))
int(ZNuc(iNuc)), &
max_exponent(iNuc), &
max_ang_mom(iNuc), &
min_exponent(iNuc,:) )
center_index = iNuc - 1 center_index = iNuc - 1
num_points = numgrid_get_num_grid_points(context) num_points = numgrid_get_num_grid_points(context)
call numgrid_get_grid(context, & call numgrid_get_grid(context,nNuc,center_index,rNuc(:,1),rNuc(:,2),rNuc(:,3),int(ZNuc(:)), &
nNuc, & root(1,iG+1:num_points),root(2,iG+1:num_points),root(3,iG+1:num_points), &
center_index, & weight(iG+1:num_points))
rNuc(:,1), &
rNuc(:,2), &
rNuc(:,3), &
int(ZNuc(:)), &
root(1,iG+1:num_points), &
root(2,iG+1:num_points), &
root(3,iG+1:num_points), &
weight(iG+1:num_points) )
iG = iG + num_points iG = iG + num_points

View File

@ -133,18 +133,18 @@ program eDFT
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! Construct quadrature grid ! Construct quadrature grid
!------------------------------------------------------------------------ !------------------------------------------------------------------------
call read_grid(SGn,nRad,nAng,nGrid) call read_grid(SGn,radial_precision,nRad,nAng,nGrid)
radial_precision = 1d-7 nGrid = nRad*nAng
call allocate_grid(nNuc,ZNuc,rNuc,nShell,TotAngMomShell,ExpShell,max_ang_mom,min_exponent,max_exponent, & ! call allocate_grid(nNuc,ZNuc,rNuc,nShell,TotAngMomShell,ExpShell,max_ang_mom,min_exponent,max_exponent, &
radial_precision,nRad,nAng,nGrid) ! radial_precision,nRad,nAng,nGrid)
allocate(root(ncart,nGrid),weight(nGrid)) allocate(root(ncart,nGrid),weight(nGrid))
call build_grid(nNuc,ZNuc,rNuc,nShell,TotAngMomShell,ExpShell,max_ang_mom,min_exponent,max_exponent, & call quadrature_grid(nRad,nAng,nGrid,root,weight)
radial_precision,nRad,nAng,nGrid,weight,root)
! call quadrature_grid(nRad,nAng,nGrid,root,weight) ! call build_grid(nNuc,ZNuc,rNuc,nShell,TotAngMomShell,ExpShell,max_ang_mom,min_exponent,max_exponent, &
! radial_precision,nRad,nAng,nGrid,weight,root)
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! Calculate AO values at grid points ! Calculate AO values at grid points

View File

@ -1,4 +1,4 @@
subroutine read_grid(SGn,nRad,nAng,nGrid) subroutine read_grid(SGn,radial_precision,nRad,nAng)
! Read grid type ! Read grid type
@ -10,9 +10,9 @@ subroutine read_grid(SGn,nRad,nAng,nGrid)
! Output variables ! Output variables
double precision,intent(out) :: radial_precision
integer,intent(out) :: nRad integer,intent(out) :: nRad
integer,intent(out) :: nAng integer,intent(out) :: nAng
integer,intent(out) :: nGrid
write(*,*)'----------------------------------------------------------' write(*,*)'----------------------------------------------------------'
write(*,'(A22,I1)')' Quadrature grid: SG-',SGn write(*,'(A22,I1)')' Quadrature grid: SG-',SGn
@ -21,18 +21,22 @@ subroutine read_grid(SGn,nRad,nAng,nGrid)
select case (SGn) select case (SGn)
case(0) case(0)
nRad = 23 radial_precision = 1d-5
nAng = 170 nRad = 23
nAng = 170
case(1) case(1)
radial_precision = 1d-7
nRad = 50 nRad = 50
nAng = 194 nAng = 194
case(2) case(2)
radial_precision = 1d-9
nRad = 75 nRad = 75
nAng = 302 nAng = 302
case(3) case(3)
radial_precision = 1d-11
nRad = 99 nRad = 99
nAng = 590 nAng = 590
@ -42,6 +46,4 @@ subroutine read_grid(SGn,nRad,nAng,nGrid)
end select end select
nGrid = nRad*nAng
end subroutine read_grid end subroutine read_grid