TABLE OF CONTENTS
ABINIT/m_builtin_tests [ Modules ]
NAME
m_builtin_tests
FUNCTION
COPYRIGHT
Copyright (C) 1998-2024 ABINIT group (DCA,XG,GMR) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_builtin_tests 23 24 use defs_basis 25 use m_errors 26 use m_abicore 27 28 use m_io_tools, only : open_file 29 30 implicit none 31 32 private
ABINIT/testfi [ Functions ]
NAME
testfi
FUNCTION
Routine "Final test" for generation of the test report in the status file: if it appears that the run was a "Built-in Test", then compare the final values of different quantities to the reference values, here hard-coded.
INPUTS
builtintest=number of the builtintest, from the input file. etotal=total energy (sum of 7 contributions) (hartree) filstat=name of the status file gred(3,natom)=symmetrized gradient of etotal with respect to tn natom=number of atoms in cell. strten(6)=components of the stress tensor (hartree/bohr^3) xred(3,natom)=reduced atomic coordinates
OUTPUT
(only writing)
SOURCE
66 subroutine testfi(builtintest,etotal,filstat,gred,natom,strten,xred) 67 68 !Arguments ------------------------------------ 69 !scalars 70 integer,intent(in) :: builtintest,natom 71 real(dp),intent(in) :: etotal 72 character(len=fnlen),intent(in) :: filstat 73 !arrays 74 real(dp),intent(in) :: gred(3,natom),strten(6),xred(3,natom) 75 76 !Local variables------------------------------- 77 character(len=fnlen) :: testname(7)=' ' 78 character(len=*), parameter :: format01000="(a,d22.14,a,d12.4)" 79 character(len=500) :: msg 80 !scalars 81 integer,parameter :: mtest=7 82 integer :: iatom,ii,problem,tok,temp_unit 83 real(dp) :: etot_mxdev,etot_ref,gred_mxdev,strten_mxdev,xred_mxdev 84 !arrays 85 integer,parameter :: natom_test(mtest)=(/2,1,2,1,1,1,2/) 86 real(dp) :: gred_ref(3,2),strten_ref(6),xred_ref(3,2) 87 88 ! *********************************************************************** 89 90 testname(1)='fast' 91 testname(2)='v1' 92 testname(3)='v5' 93 testname(4)='bigdft' 94 testname(5)='etsf_io' 95 testname(6)='libxc' 96 testname(7)='wannier90' 97 98 !--------------------------------------------------------- 99 100 !Now, open the status file, and either delete it, or produce a report 101 if (open_file(filstat,msg,newunit=temp_unit,form='formatted',status='unknown') /= 0) then 102 ABI_ERROR(msg) 103 end if 104 105 if(builtintest==0)then 106 107 close (temp_unit,status='delete') 108 109 else 110 111 ! Note: all processors have their own file, so no special attention must be paid to the parallel case. 112 write(temp_unit,*) 113 write(temp_unit,*)'Status file, reporting on built-in test ',trim(testname(builtintest)) 114 write(temp_unit,*) 115 116 ! Define reference values, as well as maximum tolerable deviation 117 if(builtintest==1)then 118 119 etot_ref=-1.05814441948188d+00 120 etot_mxdev=1.0d-9 121 xred_ref(1:3,1)=(/ -0.65048430042634D-01 , 0.0_dp , 0.0_dp /) 122 xred_ref(1:3,2)=(/ 0.65048430042634D-01 , 0.0_dp , 0.0_dp /) 123 ! xred(*,*) are reduced coordinates 124 xred_mxdev=1.0d-6 125 gred_ref(1:3,1:2)= 0.0_dp 126 ! gred(*,*) are gradients with respect to reduced coordinates 127 gred_mxdev=5.0d-4 128 strten_ref(1:6)=(/ 0.149D-04 , 0.560D-04 , 0.560D-04 ,& 129 & 0.0_dp , 0.0_dp , 0.0_dp /) 130 strten_mxdev=1.0d-5 131 132 else if(builtintest==2)then 133 134 ! This value of etot is accurate to about 1.0d-12 135 etot_ref=-.70811958266295D+02 136 ! Initialisation conditions might give fluctuations 137 etot_mxdev=3.0d-7 138 xred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /) 139 xred_mxdev=1.0d-12 140 gred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /) 141 gred_mxdev=1.0d-12 142 ! This value of strten is accurate to at least 1.0d-8 143 strten_ref(1:3)= 5.09324870E-03 144 strten_ref(4:6)= 0.0_dp 145 strten_mxdev=1.0d-6 146 147 else if(builtintest==3)then 148 149 etot_ref=-.892746696311772D+01 150 etot_mxdev=1.0d-8 151 xred_ref(1:3,1)=(/ -0.125_dp , 0.0_dp , 0.0_dp /) 152 xred_ref(1:3,2)=(/ 0.125_dp , 0.0_dp , 0.0_dp /) 153 xred_mxdev=1.0d-12 154 gred_ref(1:3,1)=(/ -.140263620278D+00 , 0.0_dp , 0.0_dp /) 155 gred_ref(1:3,2)=(/ .140013483725D+00 , 0.0_dp , 0.0_dp /) 156 gred_mxdev=1.0d-3 157 strten_ref(1:6)=(/ 1.3949d-3 , 1.3643d-3 , 1.3643d-3 ,.0_dp , .0_dp , .0_dp /) 158 strten_mxdev=1.0d-5 159 160 ! Bigdft 161 else if(builtintest==4)then 162 163 ! This value of etot is accurate to about 1.0d-12 164 etot_ref=-0.56231990141D+00 165 ! Initialisation conditions might give fluctuations 166 etot_mxdev=3.0d-7 167 xred_ref(1:3,1)=(/ 0.5_dp , 0.5_dp , 0.5_dp /) 168 xred_mxdev=1.0d-12 169 gred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /) 170 gred_mxdev=1.0d-12 171 strten_ref(1)=-0.22537238382D-01 172 strten_ref(2)=-0.22536232141D-01 173 strten_ref(3)=-0.22529038043D-01 174 strten_ref(4)= 0.39104928264D-06 175 strten_ref(5)= 0.62823846408D-06 176 strten_ref(6)= 0.38414125548E-09 177 strten_mxdev=1.0d-7 178 179 ! ETSF_IO 180 else if(builtintest==5)then 181 182 ! This value of etot is accurate to about 1.0d-12 183 etot_ref=-0.33307825915795D+02 184 ! Initialisation conditions might give fluctuations 185 etot_mxdev=3.0d-7 186 xred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /) 187 xred_mxdev=1.0d-12 188 gred_ref(1:3,1)=(/ 0.0_dp , -0.000000120877_dp , -0.000000164287_dp /) 189 gred_mxdev=1.0d-7 190 ! This value of strten is accurate to at least 1.0d-8 191 strten_ref(1)= 0.13569940015175D-01 192 strten_ref(2)= 0.33822108610352D-01 193 strten_ref(3)= 0.37991262607028D-01 194 strten_ref(4:6)= 0.0_dp 195 strten_mxdev=1.0d-6 196 197 ! LibXC 198 else if(builtintest==6)then 199 200 etot_ref=-0.55196688523897D+01 201 etot_mxdev=5.0d-7 202 xred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /) 203 xred_mxdev=1.0d-12 204 gred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /) 205 gred_mxdev=1.0d-12 206 strten_ref(1:3)= 0.13246699375127D-04 207 strten_ref(4:6)= 0.0_dp 208 strten_mxdev=1.0d-8 209 210 ! Wannier90 211 else if(builtintest==7)then 212 213 ! This value of etot is accurate to about 1.0d-12 214 etot_ref=-0.10620085133544D+02 215 ! Initialisation conditions might give fluctuations 216 etot_mxdev=3.0d-7 217 xred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /) 218 xred_ref(1:3,2)=(/ 0.25_dp , 0.25_dp , 0.25_dp /) 219 xred_mxdev=1.0d-12 220 gred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /) 221 gred_ref(1:3,2)=(/ 0.0_dp , 0.0_dp , 0.0_dp /) 222 gred_mxdev=1.0d-12 223 ! This value of strten is accurate to at least 1.0d-8 224 strten_ref(1:3)=0.31413922197317D-03 225 strten_ref(4:6)= 0.0_dp 226 strten_mxdev=1.0d-6 227 228 ! End of the reference value set up, for different tests 229 end if 230 231 ! Compare reference and actual values 232 233 problem=0 234 235 ! Take care of total energy 236 if(abs(etot_ref-etotal)>etot_mxdev)then 237 problem=1 238 write(temp_unit,'(a)')' Error for total energy : ' 239 write(temp_unit,format01000)' expected ',etot_ref,' with maximum deviation',etot_mxdev 240 write(temp_unit,format01000)' computed ',etotal,' with effective deviation',abs(etotal-etot_ref) 241 end if 242 243 ! Take care of nuclei positions 244 tok=1 245 do iatom=1,natom 246 do ii=1,3 247 if(abs(xred_ref(ii,iatom)-xred(ii,iatom))>xred_mxdev)then 248 tok=0 249 write(temp_unit, '(a,i1,a,i1,a)' )' Error for nuclei position xred(',ii,',',iatom,')' 250 write(temp_unit,format01000)' expected ',xred_ref(ii,iatom),' with maximum deviation',xred_mxdev 251 write(temp_unit,format01000)' computed ',xred(ii,iatom),' with effective deviation',& 252 & abs( xred(ii,iatom)-xred_ref(ii,iatom) ) 253 end if 254 end do 255 end do 256 if(tok==0)problem=1 257 258 ! Take care of forces 259 tok=1 260 do iatom=1,natom 261 do ii=1,3 262 if(abs(gred_ref(ii,iatom)-gred(ii,iatom))& 263 & >gred_mxdev)then 264 tok=0 265 write(temp_unit, '(a,i1,a,i1,a)' )' Error for gradients gred(',ii,',',iatom,')' 266 write(temp_unit,format01000)' expected ',gred_ref(ii,iatom),' with maximum deviation',gred_mxdev 267 write(temp_unit,format01000)' computed ',gred(ii,iatom),' with effective deviation',& 268 & abs( gred(ii,iatom)-gred_ref(ii,iatom) ) 269 end if 270 end do 271 end do 272 if(tok==0)problem=1 273 274 ! Take care of stress 275 tok=1 276 do ii=1,6 277 if(abs(strten_ref(ii)-strten(ii))>strten_mxdev)then 278 tok=0 279 write(temp_unit,'(a,i1,a)')' Error for stress strten(',ii,')' 280 write(temp_unit,format01000)' expected ',strten_ref(ii),' with maximum deviation',strten_mxdev 281 write(temp_unit,format01000)' computed ',strten(ii),' with effective deviation',& 282 & abs( strten(ii)-strten_ref(ii) ) 283 end if 284 end do 285 if(tok==0)problem=1 286 287 if(problem==0)then 288 write(temp_unit,'(a)')' ==> The run finished cleanly.' 289 write(temp_unit,'(a,a)')' Moreover, comparison of the total energy, and other (few) ',& 290 & 'relevant quantities with reference values has been successful.' 291 write(temp_unit,'(a)')' This does not mean that no problem is present, however. ' 292 write(temp_unit,'(a)')' Please run the complete set of ABINIT tests to gain a better confidence in your installation.' 293 end if 294 295 write(temp_unit,*) 296 297 close(temp_unit) 298 299 end if ! End of the choice between produce a report, and produce no report 300 301 end subroutine testfi