TABLE OF CONTENTS


ABINIT/m_builtin_tests [ Modules ]

[ Top ] [ Modules ]

NAME

  m_builtin_tests

FUNCTION

COPYRIGHT

  Copyright (C) 1998-2018 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 .

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_builtin_tests
28 
29  use defs_basis
30  use m_errors
31  use m_abicore
32 
33  use m_io_tools,     only : open_file
34 
35  implicit none
36 
37  private

ABINIT/testfi [ Functions ]

[ Top ] [ 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
  fred(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)

PARENTS

      abinit

CHILDREN

SOURCE

 76 subroutine testfi(builtintest,etotal,filstat,fred,natom,strten,xred)
 77 
 78 
 79 !This section has been created automatically by the script Abilint (TD).
 80 !Do not modify the following lines by hand.
 81 #undef ABI_FUNC
 82 #define ABI_FUNC 'testfi'
 83 !End of the abilint section
 84 
 85  implicit none
 86 
 87 !Arguments ------------------------------------
 88 !scalars
 89  integer,intent(in) :: builtintest,natom
 90  real(dp),intent(in) :: etotal
 91  character(len=fnlen),intent(in) :: filstat
 92 !arrays
 93  real(dp),intent(in) :: fred(3,natom),strten(6),xred(3,natom)
 94 
 95 !Local variables-------------------------------
 96  character(len=fnlen) :: testname(7)='         '
 97  character(len=*), parameter :: format01000="(a,d22.14,a,d12.4)"
 98  character(len=500) :: msg
 99 !scalars
100  integer,parameter :: mtest=7
101  integer :: iatom,ii,problem,tok,temp_unit
102  real(dp) :: etot_mxdev,etot_ref,fred_mxdev,strten_mxdev,xred_mxdev
103 !arrays
104  integer,parameter :: natom_test(mtest)=(/2,1,2,1,1,1,2/)
105  real(dp) :: fred_ref(3,2),strten_ref(6),xred_ref(3,2)
106 
107 ! ***********************************************************************
108 
109  testname(1)='fast'
110  testname(2)='v1'
111  testname(3)='v5'
112  testname(4)='bigdft'
113  testname(5)='etsf_io'
114  testname(6)='libxc'
115  testname(7)='wannier90'
116 
117 !---------------------------------------------------------
118 
119 !Now, open the status file, and either delete it, or produce a report
120  if (open_file(filstat,msg,newunit=temp_unit,form='formatted',status='unknown') /= 0) then
121    MSG_ERROR(msg)
122  end if
123 
124  if(builtintest==0)then
125 
126    close (temp_unit,status='delete')
127 
128  else
129 
130 !  Note : all processors have their own file, so no special
131 !  attention must be paid to the parallel case.
132 
133    write(temp_unit,*)
134    write(temp_unit,*)'Status file, reporting on built-in test ',trim(testname(builtintest))
135    write(temp_unit,*)
136 
137 !  Define reference values, as well as maximum tolerable deviation
138    if(builtintest==1)then
139 
140      etot_ref=-1.05814441948188d+00
141      etot_mxdev=1.0d-9
142      xred_ref(1:3,1)=(/ -0.65048430042634D-01 , 0.0_dp , 0.0_dp /)
143      xred_ref(1:3,2)=(/  0.65048430042634D-01 , 0.0_dp , 0.0_dp /)
144 !    xred(*,*) are reduced coordinates
145      xred_mxdev=1.0d-6
146      fred_ref(1:3,1:2)= 0.0_dp
147 !    fred(*,*) are gradients with respect to reduced coordinates
148      fred_mxdev=5.0d-4
149      strten_ref(1:6)=(/ 0.149D-04  , 0.560D-04 , 0.560D-04 ,&
150 &     0.0_dp , 0.0_dp , 0.0_dp /)
151      strten_mxdev=1.0d-5
152 
153    else if(builtintest==2)then
154 
155 !    This value of etot is accurate to about 1.0d-12
156      etot_ref=-.70811958266295D+02
157 !    Initialisation conditions might give fluctuations
158      etot_mxdev=3.0d-7
159      xred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /)
160      xred_mxdev=1.0d-12
161      fred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /)
162      fred_mxdev=1.0d-12
163 !    This value of strten is accurate to at least 1.0d-8
164      strten_ref(1:3)= 5.09324870E-03
165      strten_ref(4:6)= 0.0_dp
166      strten_mxdev=1.0d-6
167 
168    else if(builtintest==3)then
169 
170      etot_ref=-.892746696311772D+01
171      etot_mxdev=1.0d-8
172      xred_ref(1:3,1)=(/ -0.125_dp , 0.0_dp , 0.0_dp /)
173      xred_ref(1:3,2)=(/  0.125_dp , 0.0_dp , 0.0_dp /)
174      xred_mxdev=1.0d-12
175      fred_ref(1:3,1)=(/ -.140263620278D+00 , 0.0_dp , 0.0_dp /)
176      fred_ref(1:3,2)=(/  .140013483725D+00 , 0.0_dp , 0.0_dp /)
177      fred_mxdev=1.0d-3
178      strten_ref(1:6)=(/  1.3949d-3 ,  1.3643d-3 , 1.3643d-3 ,&
179 &     .0_dp ,  .0_dp  ,  .0_dp    /)
180      strten_mxdev=1.0d-5
181 
182 !    Bigdft
183    else if(builtintest==4)then
184 
185 !    This value of etot is accurate to about 1.0d-12
186      etot_ref=-0.56231990141D+00
187 !    Initialisation conditions might give fluctuations
188      etot_mxdev=3.0d-7
189      xred_ref(1:3,1)=(/ 0.5_dp , 0.5_dp , 0.5_dp /)
190      xred_mxdev=1.0d-12
191      fred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /)
192      fred_mxdev=1.0d-12
193      strten_ref(1)=-0.22537238382D-01
194      strten_ref(2)=-0.22536232141D-01
195      strten_ref(3)=-0.22529038043D-01
196      strten_ref(4)= 0.39104928264D-06
197      strten_ref(5)= 0.62823846408D-06
198      strten_ref(6)= 0.38414125548E-09
199      strten_mxdev=1.0d-7
200 
201 !    ETSF_IO
202    else if(builtintest==5)then
203 
204 !    This value of etot is accurate to about 1.0d-12
205      etot_ref=-0.33307825915795D+02
206 !    Initialisation conditions might give fluctuations
207      etot_mxdev=3.0d-7
208      xred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /)
209      xred_mxdev=1.0d-12
210      fred_ref(1:3,1)=(/ 0.0_dp , -0.000000120877_dp , -0.000000164287_dp /)
211      fred_mxdev=1.0d-7
212 !    This value of strten is accurate to at least 1.0d-8
213      strten_ref(1)= 0.13569940015175D-01
214      strten_ref(2)= 0.33822108610352D-01
215      strten_ref(3)= 0.37991262607028D-01
216      strten_ref(4:6)= 0.0_dp
217      strten_mxdev=1.0d-6
218 
219 !    LibXC
220    else if(builtintest==6)then
221 
222      etot_ref=-0.55196688523897D+01
223      etot_mxdev=5.0d-7
224      xred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /)
225      xred_mxdev=1.0d-12
226      fred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /)
227      fred_mxdev=1.0d-12
228      strten_ref(1:3)= 0.13246699375127D-04
229      strten_ref(4:6)= 0.0_dp
230      strten_mxdev=1.0d-8
231 
232 !    Wannier90
233    else if(builtintest==7)then
234 
235 !    This value of etot is accurate to about 1.0d-12
236      etot_ref=-0.10620085133544D+02
237 !    Initialisation conditions might give fluctuations
238      etot_mxdev=3.0d-7
239      xred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /)
240      xred_ref(1:3,2)=(/ 0.25_dp , 0.25_dp , 0.25_dp /)
241      xred_mxdev=1.0d-12
242      fred_ref(1:3,1)=(/ 0.0_dp , 0.0_dp , 0.0_dp /)
243      fred_ref(1:3,2)=(/ 0.0_dp , 0.0_dp , 0.0_dp /)
244      fred_mxdev=1.0d-12
245 !    This value of strten is accurate to at least 1.0d-8
246      strten_ref(1:3)=0.31413922197317D-03
247      strten_ref(4:6)= 0.0_dp
248      strten_mxdev=1.0d-6
249 
250 !    End of the reference value set up, for different tests
251    end if
252 
253 !  Compare reference and actual values
254 
255    problem=0
256 
257 !  Take care of total energy
258    if(abs(etot_ref-etotal)>etot_mxdev)then
259      problem=1
260      write(temp_unit,'(a)')' Error for total energy : '
261      write(temp_unit,format01000)'        expected ',etot_ref,&
262 &     '  with maximum   deviation',etot_mxdev
263      write(temp_unit,format01000)'        computed ',etotal,&
264 &     '  with effective deviation',abs(etotal-etot_ref)
265    end if
266 
267 !  Take care of nuclei positions
268    tok=1
269    do iatom=1,natom
270      do ii=1,3
271        if(abs(xred_ref(ii,iatom)-xred(ii,iatom))>xred_mxdev)then
272          tok=0
273          write(temp_unit, '(a,i1,a,i1,a)' )&
274 &         ' Error for nuclei position xred(',ii,',',iatom,')'
275          write(temp_unit,format01000)'        expected ',xred_ref(ii,iatom),&
276 &         '  with maximum   deviation',xred_mxdev
277          write(temp_unit,format01000)'        computed ',xred(ii,iatom),&
278 &         '  with effective deviation',&
279 &         abs( xred(ii,iatom)-xred_ref(ii,iatom) )
280        end if
281      end do
282    end do
283    if(tok==0)problem=1
284 
285 !  Take care of forces
286    tok=1
287    do iatom=1,natom
288      do ii=1,3
289        if(abs(fred_ref(ii,iatom)-fred(ii,iatom))&
290 &       >fred_mxdev)then
291          tok=0
292          write(temp_unit, '(a,i1,a,i1,a)' )&
293 &         ' Error for force fred(',ii,',',iatom,')'
294          write(temp_unit,format01000)'        expected ',fred_ref(ii,iatom),&
295 &         '  with maximum   deviation',fred_mxdev
296          write(temp_unit,format01000)'        computed ',fred(ii,iatom),&
297 &         '  with effective deviation',&
298 &         abs( fred(ii,iatom)-fred_ref(ii,iatom) )
299        end if
300      end do
301    end do
302    if(tok==0)problem=1
303 
304 !  Take care of stress
305    tok=1
306    do ii=1,6
307      if(abs(strten_ref(ii)-strten(ii))>strten_mxdev)then
308        tok=0
309        write(temp_unit,'(a,i1,a)')' Error for stress strten(',ii,')'
310        write(temp_unit,format01000)'        expected ',strten_ref(ii),&
311 &       '  with maximum   deviation',strten_mxdev
312        write(temp_unit,format01000)'        computed ',strten(ii),&
313 &       '  with effective deviation',&
314 &       abs( strten(ii)-strten_ref(ii) )
315      end if
316    end do
317    if(tok==0)problem=1
318 
319    if(problem==0)then
320      write(temp_unit,'(a)')' ==> The run finished cleanly.'
321      write(temp_unit,'(a,a)')'     Moreover, comparison of the total energy, and other (few) ',&
322 &     'relevant quantities with reference values has been successful.'
323      write(temp_unit,'(a)')'     This does not mean that no problem is present, however. '
324      write(temp_unit,'(a)')'     Please run the complete set of ABINIT tests to gain a better confidence in your installation.'
325    end if
326 
327    write(temp_unit,*)
328 
329    close(temp_unit)
330 
331  end if !  End of the choice between produce a report, and produce no report
332 
333 end subroutine testfi