Skip to content

Commit

Permalink
Fix real precisions (tblite#204)
Browse files Browse the repository at this point in the history
* Proper type for VC/VS
* Proper type of constants
* Replace dble calls with real(..., kind=wp)
* Fix default_smoothing: it had different values depending on selected real precision
  • Loading branch information
foxtran authored Sep 26, 2024
1 parent 36ff0cd commit e192821
Show file tree
Hide file tree
Showing 5 changed files with 37 additions and 36 deletions.
6 changes: 3 additions & 3 deletions src/tblite/ceh/ceh.f90
Original file line number Diff line number Diff line change
Expand Up @@ -820,7 +820,7 @@ module tblite_ceh_ceh
& 0.1308523996_wp, 0.1635461451_wp, 0.1994570401_wp] ! 101-103

!> Empirical atomic radii for calculation of the coordination number
real(wp), parameter :: ceh_cov_radii(max_elem) = 0.5 * [&
real(wp), parameter :: ceh_cov_radii(max_elem) = 0.5_wp * [&
& 2.4040551903_wp, 1.8947380542_wp, 3.4227634078_wp, 3.5225408137_wp, & ! 1-4
& 3.6150631704_wp, 2.8649682108_wp, 2.4695867541_wp, 2.3533691180_wp, & ! 5-8
& 2.4992147462_wp, 3.3390521781_wp, 4.4665909451_wp, 4.3877250907_wp, & ! 9-12
Expand Down Expand Up @@ -850,7 +850,7 @@ module tblite_ceh_ceh

!> Empirical Pauling EN normalized to EN(F)=1 as start values
!> Used for EN-scaled Coordination number in CEH
real(wp), parameter :: pauling_en_ceh(max_elem) = (1d0/3.98d0) * [&
real(wp), parameter :: pauling_en_ceh(max_elem) = (1e0_wp/3.98e0_wp) * [&
& 1.9435211923_wp, 3.6116085622_wp, 2.4630915335_wp, 2.0658837656_wp, & ! 1-4
& 2.3619778807_wp, 2.9484294262_wp, 3.8753937411_wp, 4.6235054741_wp, & ! 5-8
& 3.9800000000_wp, 3.6615073276_wp, 2.3578254072_wp, 2.4225832022_wp, & ! 9-12
Expand Down Expand Up @@ -1292,7 +1292,7 @@ subroutine get_effective_qat(mol, bas, cn_en, qat)
izp = mol%num(isp)

qat(iat, ispin) = p_ceh_en_to_q(izp)*cn_en(iat) &
& + p_ceh_total_to_q * mol%charge/dble(mol%nat)
& + p_ceh_total_to_q * mol%charge/real(mol%nat, wp)
end do
end do

Expand Down
12 changes: 6 additions & 6 deletions src/tblite/scf/mixer/broyden.f90
Original file line number Diff line number Diff line change
Expand Up @@ -169,12 +169,12 @@ subroutine broyden(n, q, qlast, dq, dqlast, iter, memory, alpha, omega, df, u, a
it1 = mod(itn - 1, memory) + 1

! set parameters
! alpha = 0.25d0
omega0 = 0.01d0
minw = 1.0d0
maxw = 100000.0d0
wfac = 0.01d0
! wfac = 0.05d0
! alpha = 0.25e0_wp
omega0 = 0.01e0_wp
minw = 1.0e0_wp
maxw = 100000.0e0_wp
wfac = 0.01e0_wp
! wfac = 0.05e0_wp

! if case for first iteration: simple damping
if (iter == 1) then
Expand Down
5 changes: 3 additions & 2 deletions src/tblite/solvation/cpcm_dd.f90
Original file line number Diff line number Diff line change
Expand Up @@ -752,7 +752,8 @@ pure subroutine dbasis(self, x, basloc, dbsloc, vplm, vcos, vsin)
real(wp), intent(inout) :: dbsloc(:, :) ! [3, self%nylm]
real(wp), intent(inout) :: vcos(:), vsin(:) ! [self%lmax+1]

integer :: l, m, ind, VC, VS
integer :: l, m, ind
real(wp) :: VC, VS
real(wp) :: cthe, sthe, cphi, sphi, plm, fln, pp1, pm1, pp
real(wp) :: et(3), ep(3)

Expand Down Expand Up @@ -1835,7 +1836,7 @@ pure subroutine rmsvec(n, v, vrms, vmax)
end do

! the much neglected square root
vrms = sqrt(vrms/dble(n))
vrms = sqrt(vrms/real(n, wp))

end subroutine rmsvec

Expand Down
2 changes: 1 addition & 1 deletion src/tblite/solvation/surface.f90
Original file line number Diff line number Diff line change
Expand Up @@ -52,7 +52,7 @@ module tblite_solvation_surface
real(wp), parameter :: tolsesp = 1.e-6_wp

real(wp), parameter :: default_offset = 2.0_wp*aatoau
real(wp), parameter :: default_smoothing = 0.3*aatoau
real(wp), parameter :: default_smoothing = 0.3_wp*aatoau

contains

Expand Down
48 changes: 24 additions & 24 deletions test/unit/test_solvation_surface.f90
Original file line number Diff line number Diff line change
Expand Up @@ -116,12 +116,12 @@ subroutine test_mb01(error)
real(wp), parameter :: probe = 1.4_wp * aatoau
integer, parameter :: nang = 110
real(wp), parameter :: ref(16) = [&
& 1.98249964115938E+2_wp, 9.34967917451750E+1_wp, 7.26746426166825E+1_wp, &
& 3.72308704194178E+1_wp, 1.00057039382294E+2_wp, 8.72703799067120E+1_wp, &
& 1.75563552659476E+1_wp, 5.79324044451170E+1_wp, 9.81702069377092E-3_wp, &
& 1.05256238709304E+2_wp, 6.62363239088925E+1_wp, 1.44944527871945E+2_wp, &
& 3.33346854609299E+1_wp, 5.79746583890200E+1_wp, 6.69252987306994E+0_wp, &
& 4.86484695123678E+1_wp]
& 1.98249964603498E+2_wp, 9.34967918541344E+1_wp, 7.26746425976157E+1_wp, &
& 3.72308705072405E+1_wp, 1.00057039498616E+2_wp, 8.72703799995796E+1_wp, &
& 1.75563553107864E+1_wp, 5.79324044295481E+1_wp, 9.81701754804677E-3_wp, &
& 1.05256238904348E+2_wp, 6.62363240313345E+1_wp, 1.44944528018566E+2_wp, &
& 3.33346853562456E+1_wp, 5.79746582175529E+1_wp, 6.69252984752073E+0_wp, &
& 4.86484694486336E+1_wp]

call get_structure(mol, "MB16-43", "01")

Expand Down Expand Up @@ -152,12 +152,12 @@ subroutine test_mb02(error)
real(wp), parameter :: probe = 1.2_wp * aatoau
integer, parameter :: nang = 230
real(wp), parameter :: ref(16) = [&
& 2.86084867409084E+1_wp, 7.50937554738823E+1_wp, 8.05879870633801E+1_wp, &
& 8.24020440071839E+1_wp, 6.48136730333555E+1_wp, 1.97586790838718E+1_wp, &
& 4.90632286612891E+1_wp, 5.29220735196209E+1_wp, 9.14599027767086E+1_wp, &
& 1.38294851038605E+1_wp, 9.02032750404975E+1_wp, 1.13713659886359E+2_wp, &
& 9.83820273872444E+1_wp, 5.95926059064588E+1_wp, 2.96614651477179E+0_wp, &
& 1.44874751678823E+2_wp]
& 2.86084867868854E+1_wp, 7.50937555534059E+1_wp, 8.05879869880977E+1_wp, &
& 8.24020440962820E+1_wp, 6.48136730299052E+1_wp, 1.97586791688521E+1_wp, &
& 4.90632288004349E+1_wp, 5.29220735596789E+1_wp, 9.14599031786151E+1_wp, &
& 1.38294851260743E+1_wp, 9.02032751808618E+1_wp, 1.13713659875286E+2_wp, &
& 9.83820274680035E+1_wp, 5.95926059359978E+1_wp, 2.96614646358023E+0_wp, &
& 1.44874751490690E+2_wp]

call get_structure(mol, "MB16-43", "02")

Expand Down Expand Up @@ -188,12 +188,12 @@ subroutine test_mb03(error)
real(wp), parameter :: probe = 0.2_wp * aatoau
integer, parameter :: nang = 111
real(wp), parameter :: ref(16) = [&
& 4.93447390062137E+1_wp, 5.42387848008597E+1_wp, 2.58043996919982E+1_wp, &
& 3.26892802882332E+1_wp, 1.27988011902719E+1_wp, 9.45810634328472E+1_wp, &
& 3.43532468932011E+1_wp, 2.76341415226557E+1_wp, 2.74903763271451E+1_wp, &
& 2.85813017380994E+1_wp, 7.99313006099443E+1_wp, 1.26258175368922E+2_wp, &
& 5.38016574327888E+1_wp, 4.16287245551897E+1_wp, 9.95930646163274E+1_wp, &
& 2.36024718010776E+1_wp]
& 4.93447390726497E+1_wp, 5.42387849176901E+1_wp, 2.58043997374119E+1_wp, &
& 3.26892803192176E+1_wp, 1.27988010759842E+1_wp, 9.45810634518707E+1_wp, &
& 3.43532470377123E+1_wp, 2.76341416140764E+1_wp, 2.74903764017798E+1_wp, &
& 2.85813017859723E+1_wp, 7.99313005786035E+1_wp, 1.26258175473983E+2_wp, &
& 5.38016574162998E+1_wp, 4.16287245622076E+1_wp, 9.95930646536509E+1_wp, &
& 2.36024718294637E+1_wp]

call get_structure(mol, "MB16-43", "03")

Expand Down Expand Up @@ -225,12 +225,12 @@ subroutine test_mb04(error)
real(wp), allocatable :: rad(:), surface(:), dsdr(:, :, :)
real(wp), parameter :: probe = 1.4_wp * aatoau
real(wp), parameter :: ref(16) = [&
& 3.17888246831669E+1_wp, 1.10843192647448E+2_wp, 6.88322051937029E+1_wp, &
& 1.14544539434107E+2_wp, 1.70720777045875E+2_wp, 3.13678105884819E+1_wp, &
& 4.58475695169198E+1_wp, 1.93179973322729E+2_wp, 6.00038959489535E+1_wp, &
& 6.11241830519643E+1_wp, 4.51433360077588E+1_wp, 9.79240756407461E+0_wp, &
& 1.11790314128626E+2_wp, 3.26024197541769E+1_wp, 7.04914424301540E+1_wp, &
& 7.70033481249225E+1_wp]
& 3.17888247611432E+1_wp, 1.10843192696983E+2_wp, 6.88322052202638E+1_wp, &
& 1.14544539499935E+2_wp, 1.70720777217273E+2_wp, 3.13678106939536E+1_wp, &
& 4.58475696363698E+1_wp, 1.93179973492600E+2_wp, 6.00038960472752E+1_wp, &
& 6.11241830969292E+1_wp, 4.51433358912678E+1_wp, 9.79240755827738E+0_wp, &
& 1.11790314288253E+2_wp, 3.26024198194955E+1_wp, 7.04914426556603E+1_wp, &
& 7.70033482758105E+1_wp]

call get_structure(mol, "MB16-43", "04")

Expand Down

0 comments on commit e192821

Please sign in to comment.