TABLE OF CONTENTS


ABINIT/driver [ Functions ]

[ Top ] [ Functions ]

NAME

 driver

FUNCTION

 Driver for ground state, response function, screening
 and sigma calculations. The present routine drives the following operations.
 An outer loop allows computation related to different data sets.
 For each data set, either a GS calculation, a RF calculation,
 a SUS calculation, a SCR calculation or a SIGMA calculation is made.
 In both cases, the input variables are transferred in the proper variables,
 selected big arrays are allocated, then the gstate, respfn, ...  subroutines are called.

COPYRIGHT

 Copyright (C) 1999-2018 ABINIT group (XG,MKV,MM,MT,FJ)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

INPUTS

 codvsn= code version
 cpui=initial CPU time
 filnam(5)=character strings giving file names
 filstat=character strings giving name of status file
 mpi_enregs=information about MPI parallelization
 ndtset=number of datasets
 ndtset_alloc=number of datasets, corrected for allocation of at least one data set.
 npsp=number of pseudopotentials
 pspheads(npsp)=<type pspheader_type>all the important informatio from the
   pseudopotential file header, as well as the psp file name

OUTPUT

SIDE EFFECTS

  Input/Output
 results_out(0:ndtset_alloc)=<type results_out_type>contains the results
   needed for outvars, including evolving variables
   Default values are set up in the calling routine
 dtsets(0:ndtset_alloc)=<type datasets_type>
   input: all input variables initialized from the input file.
   output: the effective set of variables used in the different datasets.
           Some variables, indeed, might have been redefined in one of the children.
           of this routine.

NOTES

 The array filnam is used for the name of input and output files,
 and roots for generic input, output or temporary files.
 Pseudopotential file names are set in pspini and pspatm,
 using another name. The name filstat will be needed beyond gstate to check
 the appearance of the "exit" flag, to make a hasty exit, as well as
 in order to output the status of the computation.

PARENTS

      abinit

CHILDREN

      abi_io_redirect,abi_linalg_finalize,abi_linalg_init,bethe_salpeter
      chkdilatmx,destroy_results_out,destroy_results_respfn,dtfil_init
      dtfil_init_img,dtset_copy,dtset_free,echo_xc_name,eph,exit_check
      f_malloc_set_status,fftw3_cleanup,fftw3_init_threads,find_getdtset
      gather_results_out,gstateimg,gwls_sternheimer,init_results_respfn
      libpaw_write_comm_set,libxc_functionals_end,libxc_functionals_init
      mkrdim,mpi_environment_set,nonlinear,nullify_wvl_data,pawang_free
      pawrad_free,pawtab_free,pawtab_nullify,psps_free,psps_init_from_dtset
      psps_init_global,respfn,screening,sigma,status,timab,wfk_analyze,wrtout
      xc_vdw_done,xc_vdw_init,xc_vdw_libxc_init,xc_vdw_memcheck,xc_vdw_read
      xc_vdw_show,xc_vdw_trigger,xc_vdw_write,xcart2xred,xg_finalize
      xgscalapack_config,xmpi_bcast,xred2xcart

SOURCE

111 subroutine driver(codvsn,cpui,dtsets,filnam,filstat,&
112 &                 mpi_enregs,ndtset,ndtset_alloc,npsp,pspheads,results_out)
113 
114  use defs_basis
115  use defs_datatypes
116  use defs_abitypes
117  use defs_wvltypes
118  use m_errors
119  use m_results_out
120  use m_results_respfn
121  use m_xmpi
122  use m_xomp
123  use m_abi_linalg
124  use m_abicore
125  use m_exit
126  use libxc_functionals
127 #if defined DEV_YP_VDWXC
128  use m_xc_vdw
129 #endif
130  use m_xgScalapack
131 
132  use m_time,         only : timab
133  use m_xg,           only : xg_finalize
134  use m_libpaw_tools, only : libpaw_write_comm_set
135  use m_geometry,     only : mkrdim, xcart2xred, xred2xcart, chkdilatmx
136  use m_pawang,       only : pawang_type, pawang_free
137  use m_pawrad,       only : pawrad_type, pawrad_free
138  use m_pawtab,       only : pawtab_type, pawtab_nullify, pawtab_free
139  use m_fftw3,        only : fftw3_init_threads, fftw3_cleanup
140  use m_psps,         only : psps_init_global, psps_init_from_dtset, psps_free
141  use m_dtset,        only : dtset_copy, dtset_free, find_getdtset
142  use m_mpinfo,       only : mpi_distrib_is_ok
143  use m_dtfil,        only : dtfil_init, dtfil_init_img, status
144  use m_respfn_driver,    only : respfn
145  use m_screening_driver, only : screening
146  use m_sigma_driver,     only : sigma
147  use m_bethe_salpeter,   only : bethe_salpeter
148  use m_eph_driver,       only : eph
149  use m_wfk_analyze,      only : wfk_analyze
150  use m_gstateimg,        only : gstateimg
151  use m_gwls_sternheimer, only : gwls_sternheimer
152  use m_nonlinear,        only : nonlinear
153  use m_drivexc,          only : echo_xc_name
154 
155 #if defined HAVE_BIGDFT
156  use BigDFT_API,   only: xc_init, xc_end, XC_MIXED, XC_ABINIT,&
157 &                        mpi_environment_set,bigdft_mpi, f_malloc_set_status
158 #endif
159 
160 !This section has been created automatically by the script Abilint (TD).
161 !Do not modify the following lines by hand.
162 #undef ABI_FUNC
163 #define ABI_FUNC 'driver'
164 !End of the abilint section
165 
166  implicit none
167 
168  !Arguments ------------------------------------
169  !scalars
170  integer,intent(in) :: ndtset,ndtset_alloc,npsp
171  real(dp),intent(in) :: cpui
172  character(len=6),intent(in) :: codvsn
173  character(len=fnlen),intent(in) :: filstat
174  type(MPI_type),intent(inout) :: mpi_enregs(0:ndtset_alloc)
175  !arrays
176  character(len=fnlen),intent(in) :: filnam(5)
177  type(dataset_type),intent(inout) :: dtsets(0:ndtset_alloc)
178  type(pspheader_type),intent(in) :: pspheads(npsp)
179  type(results_out_type),target,intent(inout) :: results_out(0:ndtset_alloc)
180 
181  !Local variables-------------------------------
182  !scalars
183 #if defined DEV_YP_VDWXC
184  character(len=fnlen) :: vdw_filnam
185 #endif
186  integer,parameter :: level=2,mdtset=9999
187  integer,save :: paw_size_old=-1
188  integer :: idtset,ierr,iexit,iget_cell,iget_occ,iget_vel,iget_xcart,iget_xred
189  integer :: ii,iimage,iimage_get,jdtset,jdtset_status,jj,kk
190  integer :: mtypalch,mu,mxnimage,nimage,openexit,paw_size,prtvol
191  real(dp) :: etotal
192  character(len=500) :: message
193  character(len=500) :: dilatmx_errmsg
194  logical :: converged,results_gathered,test_img,use_results_all
195  type(dataset_type) :: dtset
196  type(datafiles_type) :: dtfil
197  type(pawang_type) :: pawang
198  type(pseudopotential_type) :: psps
199  type(results_respfn_type) :: results_respfn
200  type(wvl_data) :: wvl
201 #if defined DEV_YP_VDWXC
202  type(xc_vdw_type) :: vdw_params
203 #endif
204  !arrays
205  integer :: mkmems(3)
206  integer,allocatable :: jdtset_(:),npwtot(:)
207  real(dp) :: acell(3),rprim(3,3),rprimd(3,3),tsec(2)
208  real(dp),allocatable :: acell_img(:,:),amu_img(:,:),rprim_img(:,:,:)
209  real(dp),allocatable :: fcart_img(:,:,:),fred_img(:,:,:)
210  real(dp),allocatable :: etotal_img(:),mixalch_img(:,:,:),strten_img(:,:),miximage(:,:)
211  real(dp),allocatable :: occ(:),xcart(:,:),xred(:,:),xredget(:,:)
212  real(dp),allocatable :: occ_img(:,:),vel_cell_img(:,:,:),vel_img(:,:,:),xred_img(:,:,:)
213  type(pawrad_type),allocatable :: pawrad(:)
214  type(pawtab_type),allocatable :: pawtab(:)
215  type(results_out_type),pointer :: results_out_all(:)
216 
217 !******************************************************************
218 
219  DBG_ENTER("COLL")
220 
221  call timab(640,1,tsec)
222  call timab(641,3,tsec)
223  call status(0,filstat,iexit,level,'enter         ')
224 
225 !Structured debugging if prtvol==-level
226  prtvol=dtsets(1)%prtvol
227  if(prtvol==-level)then
228    write(message,'(80a,a,a)')  ('=',ii=1,80),ch10,' driver : enter , debug mode '
229    call wrtout(std_out,message,'COLL')
230  end if
231 
232  if(ndtset>mdtset)then
233    write(message,'(a,i2,a,i5,a)')'  The maximal allowed ndtset is ',mdtset,' while the input value is ',ndtset,'.'
234    MSG_BUG(message)
235  end if
236 
237  mtypalch=dtsets(1)%ntypalch
238  do ii=1,ndtset_alloc
239    mtypalch=max(dtsets(ii)%ntypalch,mtypalch)
240  end do
241  call psps_init_global(mtypalch, npsp, psps, pspheads)
242 
243  ABI_ALLOCATE(jdtset_,(0:ndtset))
244  if(ndtset/=0)then
245    jdtset_(:)=dtsets(0:ndtset)%jdtset
246  else
247    jdtset_(0)=0
248  end if
249 
250  do idtset=1,ndtset_alloc
251    call init_results_respfn(dtsets,ndtset_alloc,results_respfn)
252  end do
253 
254  call timab(641,2,tsec)
255 
256 !*********************************************************************
257 !Big loop on datasets
258 
259 !Do loop on idtset (allocate statements are present)
260  do idtset=1,ndtset_alloc
261 
262    if(mpi_enregs(idtset)%me<0) cycle
263 
264    call timab(642,1,tsec)
265    call abi_io_redirect(new_io_comm=mpi_enregs(idtset)%comm_world)
266    call libpaw_write_comm_set(mpi_enregs(idtset)%comm_world)
267 
268    jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=1
269 
270    if(ndtset>=2)then
271      jdtset_status=jdtset
272    else
273      jdtset_status=0
274    end if
275 
276    call status(jdtset_status,filstat,iexit,level,'loop jdtset   ')
277 
278 !  Copy input variables into a local dtset.
279    call dtset_copy(dtset, dtsets(idtset))
280 
281 !  Set other values
282    dtset%jdtset = jdtset
283    dtset%ndtset = ndtset
284 
285 !  Print DATASET number
286    write(message,'(83a)') ch10,('=',mu=1,80),ch10,'== DATASET'
287    if (jdtset>=100) then
288      write(message,'(2a,i4,65a)') trim(message),' ',jdtset,' ',('=',mu=1,64)
289    else
290      write(message,'(2a,i2,67a)') trim(message),' ',jdtset,' ',('=',mu=1,66)
291    end if
292    write(message,'(3a,i5)') trim(message),ch10,'-   nproc =',mpi_enregs(idtset)%nproc
293    if (.not.mpi_distrib_is_ok(mpi_enregs(idtset),dtset%mband,dtset%nkpt,&
294 &   dtset%mkmem,dtset%nsppol)) then
295      write(message,'(2a)') trim(message),&
296 &     '   -> not optimal: autoparal keyword recommended in input file'
297    end if
298    write(message,'(3a)') trim(message),ch10,' '
299    call wrtout(ab_out,message,'COLL')
300    call wrtout(std_out,message,'PERS')     ! PERS is choosen to make debugging easier
301 
302    if ( dtset%np_slk == 0 ) then
303      call xgScalapack_config(SLK_DISABLED,dtset%slk_rankpp)
304    else if ( dtset%np_slk == 1000000 ) then
305      call xgScalapack_config(SLK_AUTO,dtset%slk_rankpp)
306    else if ( dtset%np_slk > 1 ) then
307      call xgScalapack_config(dtset%np_slk,dtset%slk_rankpp)
308    else
309      call xgScalapack_config(SLK_AUTO,dtset%slk_rankpp)
310    end if
311 
312 !  Copy input values
313    mkmems(1) = dtset%mkmem
314    mkmems(2) = dtset%mkqmem
315    mkmems(3) = dtset%mk1mem
316 
317 !  Initialize MPI data for the parallelism over images
318    nimage=mpi_enregs(idtset)%my_nimage
319 
320 !  Retrieve evolving arrays (most of them depend on images)
321    ABI_ALLOCATE(acell_img,(3,nimage))
322    ABI_ALLOCATE(amu_img,(dtset%ntypat,nimage))
323    ABI_ALLOCATE(mixalch_img,(dtset%npspalch,dtset%ntypalch,nimage))
324    ABI_ALLOCATE(occ_img,(dtset%mband*dtset%nkpt*dtset%nsppol,nimage))
325    ABI_ALLOCATE(rprim_img,(3,3,nimage))
326    ABI_ALLOCATE(vel_img,(3,dtset%natom,nimage))
327    ABI_ALLOCATE(vel_cell_img,(3,3,nimage))
328    ABI_ALLOCATE(xred_img,(3,dtset%natom,nimage))
329    do iimage=1,nimage
330      ii=mpi_enregs(idtset)%my_imgtab(iimage)
331      acell_img   (:  ,iimage) = dtset%acell_orig(:,ii)
332      amu_img   (:  ,iimage)   = dtset%amu_orig(1:dtset%ntypat,ii)
333      mixalch_img   (:,:,iimage) =dtset%mixalch_orig(1:dtset%npspalch,1:dtset%ntypalch,ii)
334      rprim_img   (:,:,iimage) = dtset%rprim_orig(:,:,ii)
335      vel_img     (:,:,iimage) = dtset%vel_orig(:,1:dtset%natom,ii)
336      vel_cell_img(:,:,iimage) = dtset%vel_cell_orig(:,:,ii)
337      xred_img (:,:,iimage) = dtset%xred_orig(:,1:dtset%natom,ii)
338      occ_img  (:  ,iimage) = dtset%occ_orig(1:dtset%mband*dtset%nkpt*dtset%nsppol,ii)
339    end do
340 
341 !  ****************************************************************************
342 !  Treat the file names (get variables)
343 
344 !  In the case of multiple images, the file names will be overwritten later (for each image)
345    call dtfil_init(dtfil,dtset,filnam,filstat,idtset,jdtset_,mpi_enregs(idtset),ndtset)
346    if (dtset%optdriver==RUNL_GSTATE.and.dtset%nimage>1) then
347      call dtfil_init_img(dtfil,dtset,dtsets,idtset,jdtset_,ndtset,ndtset_alloc)
348    end if
349 
350 !  ****************************************************************************
351 !  Treat other get variables
352 
353 !  If multi dataset mode, and already the second dataset,
354 !  treatment of other get variables.
355    if( ndtset>1 .and. idtset>1 )then
356 
357 !    Check if parallelization over images is activated
358      mxnimage=maxval(dtsets(1:ndtset_alloc)%nimage)
359      ABI_ALLOCATE(miximage,(mxnimage,mxnimage))
360      test_img=(mxnimage/=1.and.maxval(dtsets(:)%npimage)>1)
361      use_results_all=.false.
362      if (test_img.and.mpi_enregs(idtset)%me_cell==0) then
363        use_results_all=.true.
364        ABI_DATATYPE_ALLOCATE(results_out_all,(0:ndtset_alloc))
365      else
366        results_out_all => results_out
367      end if
368      results_gathered=.false.
369 
370      call find_getdtset(dtsets,dtset%getocc,'getocc',idtset,iget_occ,miximage,mxnimage,ndtset_alloc)
371      if(iget_occ/=0)then
372 !      Gather contributions to results_out from images, if needed
373        if (test_img.and.mpi_enregs(idtset)%me_cell==0.and.(.not.results_gathered)) then
374          call gather_results_out(dtsets,mpi_enregs,results_out,results_out_all,use_results_all, &
375 &         allgather=.true.,only_one_per_img=.true.)
376          results_gathered=.true.
377        end if
378        if ((.not.test_img).or.mpi_enregs(idtset)%me_cell==0) then
379          do iimage=1,nimage
380            ii=mpi_enregs(idtset)%my_imgtab(iimage)
381            occ_img(:,iimage)=zero
382            do iimage_get=1,dtsets(iget_occ)%nimage
383              do jj=1,dtset%mband*dtset%nkpt*dtset%nsppol
384                occ_img(jj,iimage)=occ_img(jj,iimage)+miximage(ii,iimage_get) &
385 &               *results_out_all(iget_occ)%occ(jj,iimage_get)
386              end do
387            end do
388          end do
389        end if
390      end if
391 
392 !    Getcell has to be treated BEFORE getxcart since acell and rprim will be used
393      call find_getdtset(dtsets,dtset%getcell,'getcell',idtset,iget_cell,miximage,mxnimage,ndtset_alloc)
394      if(iget_cell/=0)then
395 !      Gather contributions to results_out from images, if needed
396        if (test_img.and.mpi_enregs(idtset)%me_cell==0.and.(.not.results_gathered)) then
397          call gather_results_out(dtsets,mpi_enregs,results_out,results_out_all,use_results_all, &
398 &         allgather=.true.,only_one_per_img=.true.)
399          results_gathered=.true.
400        end if
401        if ((.not.test_img).or.mpi_enregs(idtset)%me_cell==0) then
402          do iimage=1,nimage
403            ii=mpi_enregs(idtset)%my_imgtab(iimage)
404            acell_img(:,iimage)=zero
405            rprim_img(:,:,iimage)=zero
406            do iimage_get=1,dtsets(iget_cell)%nimage
407              do jj=1,3
408                acell_img(jj,iimage)=acell_img(jj,iimage)+&
409 &               miximage(ii,iimage_get)*results_out_all(iget_cell)%acell(jj,iimage_get)
410                do kk=1,3
411                  rprim_img(kk,jj,iimage)=rprim_img(kk,jj,iimage)+&
412 &                 miximage(ii,iimage_get)*results_out_all(iget_cell)%rprim(kk,jj,iimage_get)
413                end do
414              end do
415              dtset%rprimd_orig(:,:,iimage)=dtsets(iget_cell)%rprimd_orig(:,:,iimage)
416 !            Check that the new acell and rprim are consistent with the input dilatmx
417              call mkrdim(acell_img(:,iimage),rprim_img(:,:,iimage),rprimd)
418              call chkdilatmx(dtset%chkdilatmx,dtset%dilatmx,rprimd,&
419 &             dtset%rprimd_orig(1:3,1:3,iimage), dilatmx_errmsg)
420              _IBM6("dilatxm_errmsg: "//TRIM(dilatmx_errmsg))
421              if (LEN_TRIM(dilatmx_errmsg) /= 0) then
422                acell_img(1,iimage) = sqrt(sum(rprimd(:,1)**2))
423                acell_img(2,iimage) = sqrt(sum(rprimd(:,2)**2))
424                acell_img(3,iimage) = sqrt(sum(rprimd(:,3)**2))
425                rprim_img(:,1,iimage) = rprimd(:,1) / acell_img(1,iimage)
426                rprim_img(:,2,iimage) = rprimd(:,2) / acell_img(1,iimage)
427                rprim_img(:,3,iimage) = rprimd(:,3) / acell_img(1,iimage)
428                MSG_WARNING(dilatmx_errmsg)
429              end if
430            end do
431          end do
432        end if
433      end if
434 
435      call find_getdtset(dtsets,dtset%getxred,'getxred',idtset,iget_xred,miximage,mxnimage,ndtset_alloc)
436      if(iget_xred/=0)then
437 !      Gather contributions to results_out from images, if needed
438        if (test_img.and.mpi_enregs(idtset)%me_cell==0.and.(.not.results_gathered)) then
439          call gather_results_out(dtsets,mpi_enregs,results_out,results_out_all,use_results_all, &
440 &         allgather=.true.,only_one_per_img=.true.)
441          results_gathered=.true.
442        end if
443        if ((.not.test_img).or.mpi_enregs(idtset)%me_cell==0) then
444          do iimage=1,nimage
445            ii=mpi_enregs(idtset)%my_imgtab(iimage)
446            xred_img(:,:,iimage)=zero
447            do iimage_get=1,dtsets(iget_xred)%nimage
448              do jj=1,dtset%natom
449                do kk=1,3
450                  xred_img(kk,jj,iimage)=xred_img(kk,jj,iimage)+&
451 &                 miximage(ii,iimage_get)*results_out_all(iget_xred)%xred(kk,jj,iimage_get)
452                end do
453              end do
454            end do
455          end do
456        end if
457      end if
458 
459      call find_getdtset(dtsets,dtset%getxcart,'getxcart',idtset,iget_xcart,miximage,mxnimage,ndtset_alloc)
460      if(iget_xcart/=0)then
461 !      Gather contributions to results_out from images, if needed
462        if (test_img.and.mpi_enregs(idtset)%me_cell==0.and.(.not.results_gathered)) then
463          call gather_results_out(dtsets,mpi_enregs,results_out,results_out_all,use_results_all, &
464 &         allgather=.true.,only_one_per_img=.true.)
465          results_gathered=.true.
466        end if
467        if ((.not.test_img).or.mpi_enregs(idtset)%me_cell==0) then
468          ABI_ALLOCATE(xcart,(3,dtset%natom))
469          ABI_ALLOCATE(xredget,(3,dtset%natom))
470          do iimage=1,nimage
471            ii=mpi_enregs(idtset)%my_imgtab(iimage)
472            xred_img(:,:,iimage)=zero
473            do iimage_get=1,dtsets(iget_xcart)%nimage
474 !            Compute xcart of the previous dataset
475              call  mkrdim(results_out_all(iget_xcart)%acell(:,iimage_get),&
476 &             results_out_all(iget_xcart)%rprim(:,:,iimage_get),rprimd)
477              do jj=1,dtset%natom
478                do kk=1,3
479                  xredget (kk,jj)=results_out_all(iget_xcart)%xred(kk,jj,iimage_get)
480                end do
481              end do
482              call xred2xcart(dtset%natom,rprimd,xcart,xredget)
483 !            xcart from previous dataset is computed. Now, produce xred for the new dataset,
484 !            with the new acell and rprim ...
485              call mkrdim(acell_img(:,iimage),rprim_img(:,:,iimage),rprimd)
486              call xcart2xred(dtset%natom,rprimd,xcart,xredget(:,:))
487              do jj=1,dtset%natom
488                do kk=1,3
489                  xred_img(kk,jj,iimage)=xred_img(kk,jj,iimage)+miximage(ii,iimage_get)*xredget(kk,jj)
490                end do
491              end do
492            end do
493          end do
494          ABI_DEALLOCATE(xcart)
495          ABI_DEALLOCATE(xredget)
496        end if
497      end if
498 
499      call find_getdtset(dtsets,dtset%getvel,'getvel',idtset,iget_vel,miximage,mxnimage,ndtset_alloc)
500      if(iget_vel/=0)then
501 !      Gather contributions to results_out from images, if needed
502        if (test_img.and.mpi_enregs(idtset)%me_cell==0.and.(.not.results_gathered)) then
503          call gather_results_out(dtsets,mpi_enregs,results_out,results_out_all,use_results_all, &
504 &         allgather=.true.,only_one_per_img=.true.)
505          results_gathered=.true.
506        end if
507        if ((.not.test_img).or.mpi_enregs(idtset)%me_cell==0) then
508          do iimage=1,nimage
509            ii=mpi_enregs(idtset)%my_imgtab(iimage)
510            vel_img(:,:,iimage)=zero
511            vel_cell_img(:,:,iimage)=zero
512            do iimage_get=1,dtsets(iget_vel)%nimage
513              do jj=1,dtset%natom
514                do kk=1,3
515                  vel_img(kk,jj,iimage)=vel_img(kk,jj,iimage)+&
516 &                 miximage(ii,iimage_get)*results_out_all(iget_vel)%vel(kk,jj,iimage_get)
517                  vel_cell_img(kk,jj,iimage)=vel_cell_img(kk,jj,iimage)+&
518 &                 miximage(ii,iimage_get)*results_out_all(iget_vel)%vel_cell(kk,jj,iimage_get)
519                end do
520              end do
521            end do
522          end do
523        end if
524      end if
525 
526 !    In the case of parallelization over images, has to distribute data
527      if (test_img) then
528        if (iget_occ/=0) then
529          call xmpi_bcast(occ_img,0,mpi_enregs(idtset)%comm_cell,ierr)
530        end if
531        if (iget_cell/=0) then
532          call xmpi_bcast(acell_img,0,mpi_enregs(idtset)%comm_cell,ierr)
533          call xmpi_bcast(rprim_img,0,mpi_enregs(idtset)%comm_cell,ierr)
534        end if
535        if (iget_vel/=0) then
536          call xmpi_bcast(vel_img,0,mpi_enregs(idtset)%comm_cell,ierr)
537          call xmpi_bcast(vel_cell_img,0,mpi_enregs(idtset)%comm_cell,ierr)
538        end if
539        if (iget_xred/=0.or.iget_xcart/=0) then
540          call xmpi_bcast(xred_img,0,mpi_enregs(idtset)%comm_cell,ierr)
541        end if
542      end if
543 
544 !    Clean memory
545      ABI_DEALLOCATE(miximage)
546      if (test_img.and.mpi_enregs(idtset)%me_cell==0) then
547        if (results_gathered) then
548          call destroy_results_out(results_out_all)
549        end if
550        ABI_DATATYPE_DEALLOCATE(results_out_all)
551      else
552        nullify(results_out_all)
553      end if
554 
555    end if
556 
557 !  ****************************************************************************
558 !  Treat the pseudopotentials : initialize the psps/PAW variable
559 
560    call status(jdtset_status,filstat,iexit,level,'init psps     ')
561    call psps_init_from_dtset(dtset, idtset, psps, pspheads)
562 
563 !  The correct dimension of pawrad/tab is ntypat. In case of alchemical psps
564 !  pawrad/tab(ipsp) is invoked with ipsp<=npsp. So, in order to avoid any problem,
565 !  declare pawrad/tab at paw_size=max(ntypat,npsp).
566    paw_size=0;if (psps%usepaw==1) paw_size=max(dtset%ntypat,dtset%npsp)
567    if (paw_size/=paw_size_old) then
568      if (paw_size_old/=-1) then
569        call pawrad_free(pawrad)
570        call pawtab_free(pawtab)
571        ABI_DATATYPE_DEALLOCATE(pawrad)
572        ABI_DATATYPE_DEALLOCATE(pawtab)
573      end if
574      ABI_DATATYPE_ALLOCATE(pawrad,(paw_size))
575      ABI_DATATYPE_ALLOCATE(pawtab,(paw_size))
576      call pawtab_nullify(pawtab)
577      paw_size_old=paw_size
578    end if
579 
580 !  ****************************************************************************
581 !  WVL allocations.
582 
583 !  Nullify wvl_data. It is important to do so irregardless of the value of usewvl:
584    call nullify_wvl_data(wvl)
585 
586 !  Set up mpi information from the dataset
587    if (dtset%usewvl == 1) then
588 #if defined HAVE_BIGDFT
589      call f_malloc_set_status(iproc=mpi_enregs(idtset)%me_wvl)
590      call mpi_environment_set(bigdft_mpi,mpi_enregs(idtset)%me_wvl,&
591 &     mpi_enregs(idtset)%nproc_wvl,mpi_enregs(idtset)%comm_wvl,&
592 &     mpi_enregs(idtset)%nproc_wvl)
593 #endif
594    end if
595 
596 !  ****************************************************************************
597 !  At this stage, all the data needed for the treatment of one dataset
598 !  have been transferred from multi-dataset arrays.
599 
600    iexit=0
601 
602 !  Smaller integer arrays
603    etotal=zero
604    ABI_ALLOCATE(npwtot,(dtset%nkpt))
605    npwtot = 0
606    ABI_ALLOCATE(xred,(3,dtset%natom))
607    ABI_ALLOCATE(occ,(dtset%mband*dtset%nkpt*dtset%nsppol))
608    if(dtset%optdriver/=RUNL_GSTATE)then
609      occ(:)=occ_img(:,1)
610      acell(:)=acell_img(:,1)
611      rprim(:,:)=rprim_img(:,:,1)
612      xred(:,:)=xred_img(:,:,1)
613    end if
614 
615 !  ****************************************************************************
616 !  Exchange-correlation
617 
618    call echo_xc_name (dtset%ixc)
619 
620    if (dtset%ixc<0) then
621      call libxc_functionals_init(dtset%ixc,dtset%nspden)
622 
623 #if defined DEV_YP_VDWXC
624      if ( (dtset%vdw_xc > 0) .and. (dtset%vdw_xc < 3) ) then
625        vdw_params%functional = dtset%vdw_xc
626        vdw_params%acutmin = dtset%vdw_df_acutmin
627        vdw_params%aratio = dtset%vdw_df_aratio
628        vdw_params%damax = dtset%vdw_df_damax
629        vdw_params%damin = dtset%vdw_df_damin
630        vdw_params%dcut = dtset%vdw_df_dcut
631        vdw_params%dratio = dtset%vdw_df_dratio
632        vdw_params%dsoft = dtset%vdw_df_dsoft
633        vdw_params%gcut = dtset%vdw_df_gcut
634        vdw_params%ndpts = dtset%vdw_df_ndpts
635        vdw_params%ngpts = dtset%vdw_df_ngpts
636        vdw_params%nqpts = dtset%vdw_df_nqpts
637        vdw_params%nrpts = dtset%vdw_df_nrpts
638        vdw_params%nsmooth = dtset%vdw_df_nsmooth
639        vdw_params%phisoft = dtset%vdw_df_phisoft
640        vdw_params%qcut = dtset%vdw_df_qcut
641        vdw_params%qratio = dtset%vdw_df_qratio
642        vdw_params%rcut = dtset%vdw_df_rcut
643        vdw_params%rsoft = dtset%vdw_df_rsoft
644        vdw_params%tolerance = dtset%vdw_df_tolerance
645        vdw_params%tweaks = dtset%vdw_df_tweaks
646        vdw_params%zab = dtset%vdw_df_zab
647        write(message,'(a,1x,a)') ch10,'[vdW-DF] *** Before init ***'
648        call wrtout(std_out,message,'COLL')
649        call xc_vdw_show(std_out,vdw_params)
650        if ( dtset%irdvdw == 1 ) then
651          write(vdw_filnam,'(a,a)') trim(filnam(3)),'_VDW.nc'
652          call xc_vdw_read(vdw_filnam)
653        else
654          call xc_vdw_init(vdw_params)
655        end if
656        call xc_vdw_libxc_init(vdw_params%functional)
657        write(message,'(a,1x,a)') ch10,'[vdW-DF] *** After init ***'
658        call wrtout(std_out,message,'COLL')
659 !      call xc_vdw_get_params(vdw_params)
660        call xc_vdw_memcheck(std_out)
661        call xc_vdw_show(std_out)
662        call xc_vdw_show(ab_out)
663 
664        write (message,'(a,1x,a,e10.3,a)')ch10,&
665 &       '[vdW-DF] activation threshold: vdw_df_threshold=',dtset%vdw_df_threshold,ch10
666        call xc_vdw_trigger(.false.)
667        call wrtout(std_out,message,'COLL')
668      end if
669 #else
670      if ( (dtset%vdw_xc > 0) .and. (dtset%vdw_xc < 3) ) then
671        write(message,'(3a)')&
672 &       'vdW-DF functionals are not fully operational yet.',ch10,&
673 &       'Action : modify vdw_xc'
674        MSG_ERROR(message)
675      end if
676 #endif
677    end if
678 
679 !  FFTW3 threads initialization
680    if (dtset%ngfft(7)/100==FFT_FFTW3)then
681      call fftw3_init_threads()
682    end if
683 
684 !  linalg initialisation:
685    call abi_linalg_init(mpi_enregs(idtset)%comm_bandspinorfft,dtset%np_slk,3*maxval(dtset%nband(:)), mpi_enregs(idtset)%me)
686 
687    call timab(642,2,tsec)
688 
689 !  ****************************************************************************
690 !  Main case selection in driver
691 
692    select case(dtset%optdriver)
693 
694    case(RUNL_GSTATE)
695 
696      ABI_ALLOCATE(fcart_img,(3,dtset%natom,nimage))
697      ABI_ALLOCATE(fred_img,(3,dtset%natom,nimage))
698      ABI_ALLOCATE(etotal_img,(nimage))
699      ABI_ALLOCATE(strten_img,(6,nimage))
700 
701      call status(jdtset_status,filstat,iexit,level,'call gstateimg')
702      call gstateimg(acell_img,amu_img,codvsn,cpui,dtfil,dtset,etotal_img,fcart_img,&
703 &     fred_img,iexit,mixalch_img,mpi_enregs(idtset),nimage,npwtot,occ_img,&
704 &     pawang,pawrad,pawtab,psps,rprim_img,strten_img,vel_cell_img,vel_img,wvl,xred_img,&
705 &     filnam,filstat,idtset,jdtset_,ndtset)
706 
707    case(RUNL_RESPFN)
708 
709      call status(jdtset_status,filstat,iexit,level,'call respfn')
710      call respfn(codvsn,cpui,dtfil,dtset,etotal,iexit,mkmems,mpi_enregs(idtset),&
711 &     npwtot,occ,pawang,pawrad,pawtab,psps,results_respfn,xred)
712 
713    case(RUNL_SCREENING)
714 
715      call status(jdtset_status,filstat,iexit,level,'call screening')
716      call screening(acell,codvsn,dtfil,dtset,pawang,pawrad,pawtab,psps,rprim)
717 
718    case(RUNL_SIGMA)
719 
720      call status(jdtset_status,filstat,iexit,level,'call sigma')
721      call sigma(acell,codvsn,dtfil,dtset,pawang,pawrad,pawtab,psps,rprim,converged)
722 
723    case(RUNL_NONLINEAR)
724 
725      call status(jdtset_status,filstat,iexit,level,'call nonlinear')
726      call nonlinear(codvsn,dtfil,dtset,etotal,iexit,mpi_enregs(idtset),npwtot,occ,&
727 &     pawang,pawrad,pawtab,psps,xred)
728 
729    case(6)
730 
731      write(message,'(3a)')&
732 &     'The optdriver value 6 has been disabled since ABINITv6.0.',ch10,&
733 &     'Action : modify optdriver in the input file.'
734      MSG_ERROR(message)
735 
736    case (RUNL_BSE)
737 
738      call status(jdtset_status,filstat,iexit,level,'call bethe_salpeter')
739      call bethe_salpeter(acell,codvsn,dtfil,dtset,pawang,pawrad,pawtab,psps,rprim,xred)
740 
741    case(RUNL_GWLS) ! For running G0W0 calculations with Lanczos basis for dielectric operator and Sternheimer equation for avoiding the use of conduction states (MC+JJL)
742 
743      ABI_ALLOCATE(fcart_img,(3,dtset%natom,nimage))
744      ABI_ALLOCATE(fred_img,(3,dtset%natom,nimage))
745      ABI_ALLOCATE(etotal_img,(nimage))
746      ABI_ALLOCATE(strten_img,(6,nimage))
747 
748      call status(jdtset_status,filstat,iexit,level,'call gwls_sternheimer')
749      call gwls_sternheimer(acell_img,amu_img,codvsn,cpui,dtfil,dtset,etotal_img,fcart_img,&
750 &     fred_img,iexit,mixalch_img,mpi_enregs(idtset),nimage,npwtot,occ_img,&
751 &     pawang,pawrad,pawtab,psps,rprim_img,strten_img,vel_cell_img,vel_img,xred_img,&
752 &     filnam,filstat,idtset,jdtset_,ndtset)
753 
754      ABI_DEALLOCATE(etotal_img)
755      ABI_DEALLOCATE(fcart_img)
756      ABI_DEALLOCATE(fred_img)
757      ABI_DEALLOCATE(strten_img)
758 
759    case (RUNL_WFK)
760      call status(jdtset_status,filstat,iexit,level,'call wfk_analyze')
761      call wfk_analyze(acell,codvsn,dtfil,dtset,pawang,pawrad,pawtab,psps,rprim,xred)
762 
763    case (RUNL_EPH)
764      call status(jdtset_status,filstat,iexit,level,'call eph      ')
765      call eph(acell,codvsn,dtfil,dtset,pawang,pawrad,pawtab,psps,rprim,xred)
766 
767    case default ! Bad value for optdriver
768      write(message,'(a,i0,4a)')&
769 &     'Unknown value for the variable optdriver: ',dtset%optdriver,ch10,&
770 &     'This is not allowed. ',ch10,&
771 &     'Action: modify optdriver in the input file.'
772      MSG_ERROR(message)
773    end select
774 
775    call timab(643,1,tsec)
776 !  ****************************************************************************
777 
778 !  Transfer of multi dataset outputs from temporaries :
779 !  acell, xred, occ rprim, and vel might be modified from their
780 !  input values
781 !  etotal, fcart, fred, and strten have been computed
782 !  npwtot was already computed before, but is stored only now
783 
784    if(dtset%optdriver==RUNL_GSTATE)then
785      do iimage=1,nimage
786        results_out(idtset)%etotal(iimage)                 =etotal_img(iimage)
787        results_out(idtset)%acell(:,iimage)                =acell_img(:,iimage)
788        results_out(idtset)%amu(1:dtset%ntypat,iimage)     =amu_img(:,iimage)
789        results_out(idtset)%rprim(:,:,iimage)              =rprim_img(:,:,iimage)
790        results_out(idtset)%strten(:,iimage)                =strten_img(:,iimage)
791        results_out(idtset)%fcart(1:3,1:dtset%natom,iimage)=fcart_img(:,:,iimage)
792        results_out(idtset)%fred(1:3,1:dtset%natom,iimage) =fred_img(:,:,iimage)
793        results_out(idtset)%mixalch(1:dtset%npspalch,1:dtset%ntypalch,iimage) &
794 &       =mixalch_img(1:dtset%npspalch,1:dtset%ntypalch,iimage)
795        results_out(idtset)%npwtot(1:dtset%nkpt,iimage)    =npwtot(1:dtset%nkpt)
796        results_out(idtset)%occ(1:dtset%mband*dtset%nkpt*dtset%nsppol,iimage)=&
797 &       occ_img(1:dtset%mband*dtset%nkpt*dtset%nsppol,iimage)
798        results_out(idtset)%vel(:,1:dtset%natom,iimage)    =vel_img(:,:,iimage)
799        results_out(idtset)%vel_cell(:,:,iimage)           =vel_cell_img(:,:,iimage)
800        results_out(idtset)%xred(:,1:dtset%natom,iimage)   =xred_img(:,:,iimage)
801      end do
802      ABI_DEALLOCATE(etotal_img)
803      ABI_DEALLOCATE(fcart_img)
804      ABI_DEALLOCATE(fred_img)
805      ABI_DEALLOCATE(strten_img)
806    else
807      results_out(idtset)%acell(:,1)                =acell(:)
808      results_out(idtset)%etotal(1)                 =etotal
809      results_out(idtset)%rprim(:,:,1)              =rprim(:,:)
810      results_out(idtset)%npwtot(1:dtset%nkpt,1)    =npwtot(1:dtset%nkpt)
811      results_out(idtset)%occ(1:dtset%mband*dtset%nkpt*dtset%nsppol,1)=&
812 &     occ(1:dtset%mband*dtset%nkpt*dtset%nsppol)
813      results_out(idtset)%xred(:,1:dtset%natom,1)   =xred(:,:)
814    end if
815    ABI_DEALLOCATE(acell_img)
816    ABI_DEALLOCATE(amu_img)
817    ABI_DEALLOCATE(mixalch_img)
818    ABI_DEALLOCATE(occ_img)
819    ABI_DEALLOCATE(rprim_img)
820    ABI_DEALLOCATE(vel_img)
821    ABI_DEALLOCATE(vel_cell_img)
822    ABI_DEALLOCATE(xred_img)
823 
824    if (dtset%ngfft(7)/100==FFT_FFTW3) then
825      call fftw3_cleanup()
826    end if
827 
828    if (dtset%ixc<0) then
829      call libxc_functionals_end()
830 
831 #if defined DEV_YP_VDWXC
832      if ( (dtset%vdw_xc > 0) .and. (dtset%vdw_xc < 10) ) then
833        if ( dtset%prtvdw /= 0 ) then
834          write(vdw_filnam,'(a,a)') trim(filnam(4)),'_VDW.nc'
835          call xc_vdw_write(vdw_filnam)
836        end if
837        call xc_vdw_show(std_out,vdw_params)
838        call xc_vdw_done(vdw_params)
839      end if
840 #endif
841    end if
842 
843 !  MG: There are routines such as GW and Berry phase that can change/compute
844 !  entries in the datatype at run-time. These changes won't be visibile
845 !  in the main output file since we are passing a copy of dtsets.
846 !  I tried to update the results with the call below but this creates
847 !  several problems in outvars since one should take into account
848 !  the new dimensions (e.g. nkptgw) and their maximum value.
849 !  For the time being, we continue to pass a copy of dtsets(idtset).
850 #if 0
851    call dtset_free(dtsets(idtset))
852    call dtset_copy(dtsets(idtset),dtset)
853 #endif
854 
855    call dtset_free(dtset)
856 
857    ABI_DEALLOCATE(occ)
858    ABI_DEALLOCATE(xred)
859    ABI_DEALLOCATE(npwtot)
860 
861    call abi_linalg_finalize()
862    call xg_finalize()
863 
864    ! Check whether exiting was required by the user.
865    ! If found then beat a hasty exit from time steps
866    openexit=1; if(dtset%chkexit==0) openexit=0
867    call exit_check(zero,dtfil%filnam_ds(1),iexit,ab_out,mpi_enregs(idtset)%comm_cell,openexit)
868 
869    call timab(643,2,tsec)
870 
871    if (iexit/=0) exit
872 
873  end do ! idtset (allocate statements are present - an exit statement is present)
874 
875 !*********************************************************************
876 
877  call timab(644,1,tsec)
878 
879 !PSP deallocation
880  call psps_free(psps)
881 
882 !XG 121126 : One should not use dtset or idtset in this section, as these might not be defined for all processors.
883 
884 !PAW deallocation
885  if (allocated(pawrad)) then
886    call pawrad_free(pawrad)
887    ABI_DATATYPE_DEALLOCATE(pawrad)
888  end if
889  if (allocated(pawtab)) then
890    call pawtab_free(pawtab)
891    ABI_DATATYPE_DEALLOCATE(pawtab)
892  end if
893  call pawang_free(pawang)
894 
895  ABI_DEALLOCATE(jdtset_)
896 
897 !Results_respfn deallocation
898  call destroy_results_respfn(results_respfn)
899 
900 #if defined HAVE_BIGDFT
901 !XG 121126 : NOTE that the next debugging section was quite problematic : indeed we are
902 !outside the loop over datasets, so the value of dtset%usewvl, that is referring to the last
903 !dtset, might not be initialized, if the dataset is NOT treated by the processor.
904 !See line 253 for the initialisation of dtset, while this might not happen if
905 !if(mpi_enregs(idtset)%me<0) cycle  , as defined by line 216
906 !if (dtset%usewvl == 1) then
907 !!  WVL - debugging
908 !call memocc_abi(0,mpi_enregs(1)%me_wvl,'count','stop')
909 !call wvl_timing(mpi_enregs(1)%me_wvl,'WFN_OPT','PR')
910 !call wvl_timing(mpi_enregs(1)%me_wvl,'              ','RE')
911 !end if
912 #endif
913 
914  call status(0,filstat,iexit,level,'exit')
915 
916  call timab(644,2,tsec)
917  call timab(640,2,tsec)
918 
919  DBG_EXIT("COLL")
920 
921 end subroutine driver

ABINIT/m_driver [ Modules ]

[ Top ] [ Modules ]

NAME

  m_driver

FUNCTION

COPYRIGHT

  Copyright (C) 2008-2018 ABINIT group ()
  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

20 #if defined HAVE_CONFIG_H
21 #include "config.h"
22 #endif
23 
24 #include "abi_common.h"
25 
26 module m_driver
27 
28  implicit none
29 
30  private