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