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