TABLE OF CONTENTS
- ABINIT/gstateimg
- ABINIT/m_gstateimg
- ABINIT/move_1geo
- ABINIT/predict_copy
- ABINIT/predictimg
- ABINIT/prtimg
ABINIT/gstateimg [ Functions ]
NAME
gstateimg
FUNCTION
Routine for conducting DFT calculations for a set of (dynamical) images
INPUTS
codvsn=code version
cpui=initial CPU time
nimage=number of images of the cell (treated by current proc)
=== Optional arguments (needed when nimage>1) ===
filnam(5)=character strings giving file names
filstat=character strings giving name of status file
idtset=index of the dataset
jdtset(0:ndtset)=actual index of the datasets
ndtset=number of datasets
OUTPUT
etotal_img=total energy, for each image fcart_img(3,natom,nimage)=forces, in cartesian coordinates, for each image gred_img(3,natom,nimage)=gradient of E wrt nuclear positions, in reduced coordinates, for each image intgres_img(nspden,natom,nimage)=gradient wrt constraints, for each image npwtot(nkpt) = total number of plane waves at each k point strten_img(6,nimage)=stress tensor, for each image
SIDE EFFECTS
acell_img(3,nimage)=unit cell length scales (bohr), for each image amu_img(ntypat,nimage)=value of mass for each atomic type, for each image dtfil <type(datafiles_type)>=variables related to files dtset <type(dataset_type)>=all input variables in this dataset | mband =maximum number of bands (IN) | mgfft =maximum single fft dimension (IN) | mkmem =maximum number of k points which can fit in core memory (IN) | mpw =maximum number of planewaves in basis sphere (large number) (IN) | natom =number of atoms in unit cell (IN) | nfft =(effective) number of FFT grid points (for this processor) (IN) | nkpt =number of k points (IN) | nspden=number of spin-density components (IN) | nsppol=number of channels for spin-polarization (1 or 2) (IN) | nsym =number of symmetry elements in space group iexit= exit flag mixalch_img(npspalch,ntypalch,nimage)=value of alchemical mixing factors,for each image mpi_enreg=MPI-parallelisation information (some already initialized, some others to be initialized here) occ_img(mband*nkpt*nsppol,nimage) = occupation number for each band and k, for each image pawang <type(pawang_type)>=paw angular mesh and related data pawrad(ntypat*usepaw) <type(pawrad_type)>=paw radial mesh and related data pawtab(ntypat*usepaw) <type(pawtab_type)>=paw tabulated starting data psps <type(pseudopotential_type)>=variables related to pseudopotentials Before entering the first time in gstateimg, a significant part of psps has been initialized : the integers dimekb,lmnmax,lnmax,mpssang,mpssoang,mpsso,mgrid, ntypat,n1xccc,usepaw,useylm, and the arrays dimensioned to npsp All the remaining components of psps are to be initialized in the call to pspini . The next time the code enters gstateimg, psps might be identical to the one of the previous dtset, in which case, no reinitialisation is scheduled in pspini.f . rprim_img(3,3,nimage)=dimensionless real space primitive translations, for each image vel_cell_img(3,3,nimage)=value of cell parameters velocities, for each image vel_img(3,natom,nimage)=value of atomic velocities,for each image xred_img(3,natom,nimage) = reduced atomic coordinates, for each image
NOTES
USE OF FFT GRIDS:
=================
In case of PAW:
---------------
Two FFT grids are used:
- A "coarse" FFT grid (defined by ecut)
for the application of the Hamiltonian on the plane waves basis.
It is defined by nfft, ngfft, mgfft, ...
Hamiltonian, wave-functions, density related to WFs (rhor here), ...
are expressed on this grid.
- A "fine" FFT grid (defined) by ecutdg)
for the computation of the density inside PAW spheres.
It is defined by nfftf, ngfftf, mgfftf, ...
Total density, potentials, ...
are expressed on this grid.
In case of norm-conserving:
---------------------------
- Only the usual FFT grid (defined by ecut) is used.
It is defined by nfft, ngfft, mgfft, ...
For compatibility reasons, (nfftf,ngfftf,mgfftf)
are set equal to (nfft,ngfft,mgfft) in that case.
In case of wavelets:
--------------------
- Only the usual FFT grid (defined by wvl_crmult) is used.
It is defined by nfft, ngfft, mgfft, ... This is strictly not
an FFT grid since its dimensions are not suited for FFTs. They are
defined by wvl_setngfft().
For compatibility reasons, (nfftf,ngfftf,mgfftf)
are set equal to (nfft,ngfft,mgfft) in that case.
TODO
Not yet possible to use restartxf in parallel when localrdwf==0
SOURCE
174 subroutine gstateimg(acell_img,amu_img,codvsn,cpui,dtfil,dtset,etotal_img,fcart_img,& 175 & gred_img,iexit,intgres_img,mixalch_img,mpi_enreg,nimage,npwtot,occ_img,& 176 & pawang,pawrad,pawtab,psps,& 177 & rprim_img,strten_img,vel_cell_img,vel_img,wvl,xred_img,& 178 & filnam,filstat,idtset,jdtset,ndtset) ! optional arguments 179 180 !Arguments ------------------------------------ 181 !scalars 182 integer,intent(in) :: nimage 183 integer,optional,intent(in) :: idtset,ndtset 184 integer,intent(inout) :: iexit 185 real(dp),intent(in) :: cpui 186 character(len=8),intent(in) :: codvsn 187 character(len=fnlen),optional,intent(in) :: filstat 188 type(MPI_type),intent(inout) :: mpi_enreg 189 type(datafiles_type),target,intent(inout) :: dtfil 190 type(dataset_type),target,intent(inout) :: dtset 191 type(pawang_type),intent(inout) :: pawang 192 type(pseudopotential_type),intent(inout) :: psps 193 type(wvl_data),intent(inout) :: wvl 194 !arrays 195 integer,optional,intent(in) :: jdtset(:) 196 integer,intent(out) :: npwtot(dtset%nkpt) 197 character(len=fnlen),optional,intent(in) :: filnam(:) 198 real(dp), intent(out) :: etotal_img(nimage),fcart_img(3,dtset%natom,nimage) 199 real(dp), intent(out) :: gred_img(3,dtset%natom,nimage) 200 real(dp), intent(out) :: intgres_img(dtset%nspden,dtset%natom,nimage) 201 real(dp), intent(out) :: strten_img(6,nimage) 202 real(dp),intent(inout) :: acell_img(3,nimage),amu_img(dtset%ntypat,nimage) 203 real(dp),intent(inout) :: mixalch_img(dtset%npspalch,dtset%ntypalch,nimage) 204 real(dp),intent(inout) :: occ_img(dtset%mband*dtset%nkpt*dtset%nsppol,nimage) 205 real(dp),intent(inout) :: rprim_img(3,3,nimage),vel_cell_img(3,3,nimage),vel_img(3,dtset%natom,nimage) 206 real(dp),intent(inout) :: xred_img(3,dtset%natom,nimage) 207 type(pawrad_type),intent(inout) :: pawrad(psps%ntypat*psps%usepaw) 208 type(pawtab_type),intent(inout) :: pawtab(psps%ntypat*psps%usepaw) 209 210 !Local variables------------------------------- 211 !Define file format for different type of files. Presently, 212 !only one file format is supported for each type of files, but this might 213 !change soon ... 214 !2 for wavefunction file, new format (version 2.0 and after) (fform) NOT USED 215 !52 for density rho(r) (fformr) 216 !102 for potential V(r) file. (fformv) NOT USED 217 !scalars 218 integer,parameter :: formeig=0,level=100,ndtpawuj=0,response=0 219 integer :: history_size,idelta,idynimage,ierr,ifirst 220 integer :: ii,iimage,ih,itimimage,itimimage_eff,itimimage_gstate,itimimage_prev,ndynimage,nocc 221 integer :: ntimimage,ntimimage_stored,ntimimage_max 222 logical :: check_conv,compute_all_images,compute_static_images 223 logical :: isVused,isARused,is_master,is_mep,is_pimd 224 logical :: call_predictor,use_hist,use_hist_prev 225 real(dp) :: delta_energy,dtion 226 character(len=500) :: hist_filename,msg 227 type(args_gs_type) :: args_gs 228 type(mep_type) :: mep_param 229 type(ga_type) :: ga_param 230 type(m1geo_type) :: m1geo_param 231 type(pimd_type) :: pimd_param 232 !arrays 233 integer,allocatable :: list_dynimage(:),scf_initialized(:) 234 character(len=60),parameter :: imagealgo_str(0:13)=(/ & 235 & 'IMAGE COPY ',& ! 0 236 & 'IMAGE STEEPEST DESCENT ',& ! 1 237 & 'STRING METHOD ',& ! 2 238 & 'METADYNAMICS ',& ! 3 239 & 'GENETIC ALGORITHM ',& ! 4 240 & 'NUDGED ELASTIC BAND ',& ! 5 241 & 'LINEAR COMBINATION OF CONSTRAINED DFT ENERGIES ',& ! 6 242 & ' ',& ! 7 243 & ' ',& ! 8 244 & 'PATH-INTEGRAL MOLECULAR DYNAMICS (LANGEVIN) ',& ! 9 245 & 'PATH-INTEGRAL MOLECULAR DYNAMICS (QUANTUM THERMAL BATH) ',& ! 10 246 & ' ',& ! 11 247 & ' ',& ! 12 248 & 'PATH-INTEGRAL MOLECULAR DYNAMICS (CHAIN OF THERMOSTATS) '/) ! 13 249 character(len=24),parameter :: stgalgo_str(0:2)=(/ & 250 & 'ORIGINAL ALGO. ',& ! 0 251 & 'SIMPLIFIED + EQUAL ARC ',& ! 1 252 & 'SIMPLIFIED + ENERGY-WGTH'/) ! 2 253 character(len=20),parameter :: nebalgo_str(0:2)=(/ & 254 & 'ORIGINAL ALGO. ',& ! 0 255 & 'IMPROVED TANGENT ',& ! 1 256 & 'CLIMBING IMAGE '/) ! 2 257 character(len=20),parameter :: mepsolver_str(0:4)=(/ & 258 & 'STEEPEST-DESCENT ',& ! 0 259 & 'QUICK-MIN OPT. ',& ! 1 260 & 'L-BFGS ',& ! 2 261 & 'GL-BFGS ',& ! 3 262 & 'ORDER 4 RUNGE-KUTTA '/) ! 4 263 real(dp) :: acell(3),rprim(3,3),rprimd(3,3),tsec(2),vel_cell(3,3) 264 real(dp),allocatable :: amass(:,:),occ(:),vel(:,:),xred(:,:) 265 type(abihist),allocatable :: hist(:),hist_prev(:) 266 type(results_img_type),pointer :: results_img(:,:),res_img(:) 267 type(scf_history_type),allocatable :: scf_history(:) 268 !type(abiforstr) :: preconforstr ! Preconditioned forces and stress ... Only needed to deallocate an internal matrix in prec_simple 269 270 ! *********************************************************************** 271 272 DBG_ENTER("COLL") 273 274 call timab(1200,1,tsec) 275 call timab(1203,3,tsec) 276 277 !Arguments check 278 if (dtset%nimage>1) then 279 if ((.not.present(filnam)).or.(.not.present(filnam)).or.(.not.present(idtset)).or.& 280 & (.not.present(ndtset)).or.(.not.present(jdtset))) then 281 write(msg,'(3a)') & 282 & 'When nimage>1, all the following argument should be present:',ch10,& 283 & 'filnam, filstat, idtset, ndtset, jdtset !' 284 ABI_BUG(msg) 285 end if 286 end if 287 288 !Set flag for the effective computation (once) of static images 289 !For the time being only set when parallelization is activated 290 !Note: if you modify this flag, do not forget to change it in outvars and outvar1 291 compute_static_images=(dtset%istatimg>0) 292 293 !Prepare the calculation, by computing flags and dimensions 294 is_pimd=(dtset%imgmov==9.or.dtset%imgmov==10.or.dtset%imgmov==13) 295 is_mep =(dtset%imgmov==1.or.dtset%imgmov== 2.or.dtset%imgmov== 5) 296 ntimimage=dtset%ntimimage 297 ntimimage_stored=ntimimage;if(is_pimd)ntimimage_stored=2 298 nocc=dtset%mband*dtset%nkpt*dtset%nsppol 299 is_master=(mpi_enreg%me_cell==0.and.mpi_enreg%me_img==0) 300 delta_energy=zero 301 302 !Management of dynamics/relaxation history (positions, forces, stresses, ...) 303 use_hist=(dtset%imgmov/=0.and.nimage>0) ; use_hist_prev=.false. 304 isVused=is_pimd;isARused=(dtset%optcell/=0) 305 if (use_hist) then 306 !Read history from file (and broadcast if MPI) 307 #if defined HAVE_NETCDF 308 use_hist_prev=(dtset%restartxf==-1.and.nimage>0) 309 #endif 310 hist_filename=trim(dtfil%filnam_ds(4))//'_HIST.nc' 311 if (use_hist_prev)then 312 ABI_MALLOC(hist_prev,(nimage)) 313 if (mpi_enreg%me_cell==0) then 314 call read_md_hist_img(hist_filename,hist_prev,isVused,isARused,& 315 & imgtab=mpi_enreg%my_imgtab) 316 end if 317 call abihist_bcast(hist_prev,0,mpi_enreg%comm_cell) 318 if (nimage>0) then 319 if (any(hist_prev(:)%mxhist/=hist_prev(1)%mxhist)) then 320 msg='History problem: all images should have the same number of time steps!' 321 ABI_ERROR(msg) 322 end if 323 use_hist_prev=(hist_prev(1)%mxhist>0) 324 if (use_hist_prev) ntimimage=ntimimage+hist_prev(1)%mxhist 325 end if 326 if (.not.use_hist_prev) then 327 call abihist_free(hist_prev) 328 ABI_FREE(hist_prev) 329 end if 330 end if 331 !Initialize a variable to write the history 332 ABI_MALLOC(hist,(nimage)) 333 call abihist_init(hist,dtset%natom,ntimimage,isVused,isARused) 334 end if ! imgmov/=0 335 336 !Allocations 337 ABI_MALLOC(occ,(nocc)) 338 ABI_MALLOC(vel,(3,dtset%natom)) 339 ABI_MALLOC(xred,(3,dtset%natom)) 340 ABI_MALLOC(results_img,(nimage,ntimimage_stored)) 341 ABI_MALLOC(list_dynimage,(dtset%ndynimage)) 342 do itimimage=1,ntimimage_stored 343 res_img => results_img(:,itimimage) 344 call init_results_img(dtset%natom,dtset%npspalch,dtset%nspden,dtset%nsppol,dtset%ntypalch,& 345 & dtset%ntypat,res_img) 346 do iimage=1,nimage 347 res_img(iimage)%acell(:) =acell_img(:,iimage) 348 res_img(iimage)%amu(:) =amu_img(:,iimage) 349 res_img(iimage)%mixalch(:,:) =mixalch_img(:,:,iimage) 350 res_img(iimage)%rprim(:,:) =rprim_img(:,:,iimage) 351 res_img(iimage)%xred(:,:) =xred_img(:,:,iimage) 352 res_img(iimage)%vel(:,:) =vel_img(:,:,iimage) 353 res_img(iimage)%vel_cell(:,:)=vel_cell_img(:,:,iimage) 354 end do 355 end do 356 ndynimage=0 357 do iimage=1,nimage 358 ii=mpi_enreg%my_imgtab(iimage) 359 if (dtset%dynimage(ii)==1) then 360 ndynimage=ndynimage+1 361 list_dynimage(ndynimage)=iimage 362 end if 363 end do 364 365 !Management of SCF history (density/WF predictions from one time step to another) 366 ABI_MALLOC(scf_history,(nimage)) 367 ABI_MALLOC(scf_initialized,(nimage)) 368 scf_initialized=0 369 history_size=-1 370 if (dtset%ntimimage<=1) then 371 if (dtset%usewvl==0.and.dtset%ionmov>0.and. & 372 & (abs(dtset%densfor_pred)==5.or.abs(dtset%densfor_pred)==6)) then 373 history_size=2 374 if(dtset%extrapwf==2) history_size=3 375 end if 376 else 377 if (abs(dtset%densfor_pred)==2.or.abs(dtset%densfor_pred)==3) history_size=0 378 if (dtset%imgwfstor==1) history_size=1 379 if (dtset%usewvl==0.and.(abs(dtset%densfor_pred)==5.or.abs(dtset%densfor_pred)==6)) history_size=2 380 end if 381 do iimage=1,nimage 382 call scf_history_nullify(scf_history(iimage)) 383 scf_history(iimage)%history_size=history_size 384 end do 385 386 !In some cases, need amass variable 387 if (use_hist) then 388 ABI_MALLOC(amass,(dtset%natom,nimage)) 389 do iimage=1,nimage 390 if (any(amu_img(:,iimage)/=amu_img(:,1))) then 391 ABI_ERROR('HIST file is not compatible with variable masses!') 392 end if 393 amass(:,iimage)=amu_emass*amu_img(dtset%typat(:),iimage) 394 end do 395 end if 396 397 !In the case of the 4th-order Runge-Kutta solver, 398 !one must have a number of step multiple of 4. 399 ntimimage_max=ntimimage;idelta=1 400 if (dtset%imgmov==2.and.mep_param%mep_solver==4) then 401 ntimimage_max=4*(ntimimage_max/4) 402 idelta=4 403 end if 404 405 !MEP search: fill in eventually the data structure mep_param 406 call mep_init(dtset,mep_param) 407 408 !GA search: fill in eventually the data structure ga_param 409 call ga_init(dtset,ga_param) 410 411 !Move 1GEO approach: fill the data structure m1geo_param 412 call m1geo_init(dtfil,dtset,m1geo_param) 413 414 !PIMD: fill in eventually the data structure pimd_param 415 call pimd_init(dtset,pimd_param,is_master) 416 dtion=one;if (is_pimd) dtion=pimd_param%dtion 417 418 call timab(1203,2,tsec) 419 420 !----------------------------------------------------------------------------------------- 421 !Big loop on the propagation of all images 422 itimimage_eff=1 ; itimimage_gstate=1 423 do itimimage=1,ntimimage 424 425 res_img => results_img(:,itimimage_eff) 426 call_predictor=(ntimimage>1) 427 428 ! If history is activated and if current image is inside it: do not compute anything 429 if (use_hist_prev) then 430 if (all(hist_prev(:)%ihist<=hist_prev(:)%mxhist)) then 431 do iimage=1,nimage 432 ih=hist_prev(iimage)%ihist 433 call abihist_copy(hist_prev(iimage),hist(iimage)) 434 call mkradim(hist_prev(iimage)%acell(:,ih),rprim,hist_prev(iimage)%rprimd(:,:,ih)) 435 res_img(iimage)%acell(:)=hist_prev(iimage)%acell(:,ih) 436 res_img(iimage)%rprim(:,:)=rprim 437 res_img(iimage)%xred(:,:)=hist_prev(iimage)%xred(:,:,ih) 438 res_img(iimage)%vel(:,:)=hist_prev(iimage)%vel(:,:,ih) 439 res_img(iimage)%vel_cell(:,:)=hist_prev(iimage)%vel_cell(:,:,ih) 440 res_img(iimage)%results_gs%fcart(:,:)=hist_prev(iimage)%fcart(:,:,ih) 441 res_img(iimage)%results_gs%strten(:)=hist_prev(iimage)%strten(:,ih) 442 res_img(iimage)%results_gs%etotal=hist_prev(iimage)%etot(ih) 443 res_img(iimage)%results_gs%energies%entropy=hist_prev(iimage)%entropy(ih) 444 call fcart2gred(res_img(iimage)%results_gs%fcart,res_img(iimage)%results_gs%gred,& 445 & hist_prev(iimage)%rprimd(:,:,ih),dtset%natom) 446 hist_prev(iimage)%ihist=hist_prev(iimage)%ihist+1 447 end do 448 !PI-QTB: skip a record in random force file 449 if (pimd_param%use_qtb==1) call pimd_skip_qtb(pimd_param) 450 !call_predictor=.false. 451 goto 110 ! This is temporary 452 end if 453 end if 454 455 call timab(1204,1,tsec) 456 call localfilnam(mpi_enreg%comm_img,mpi_enreg%comm_cell,mpi_enreg%comm_world,filnam,'_IMG',dtset%nimage) 457 compute_all_images=(compute_static_images.and.itimimage==1) 458 459 ! Print title for time step 460 if (dtset%nimage>1.or.dtset%ntimimage>1) then 461 if (dtset%prtvolimg<2) then 462 msg=ch10;if (itimimage >1) write(msg,'(2a)') ch10,ch10 463 write(msg,'(5a)') trim(msg),& 464 & '================================================================================',& 465 & ch10,' ',trim(imagealgo_str(dtset%imgmov)) 466 else 467 msg='';if (itimimage >1) msg=ch10 468 write(msg,'(5a)') trim(msg),& 469 & '--------------------------------------------------------------------------------',& 470 & ch10,' ',trim(imagealgo_str(dtset%imgmov)) 471 end if 472 if (dtset%imgmov==2) then 473 write(msg,'(6a)') trim(msg),' (',trim(stgalgo_str(mep_param%string_algo)),' + ',& 474 & trim(mepsolver_str(mep_param%mep_solver)),')' 475 end if 476 if (dtset%imgmov==5) then 477 ii=merge(mep_param%neb_algo,1,mep_param%neb_algo/=2.or.itimimage>=mep_param%cineb_start) 478 write(msg,'(6a)') trim(msg),' (',trim(nebalgo_str(ii)),' + ',& 479 & trim(mepsolver_str(mep_param%mep_solver)),')' 480 end if 481 if (dtset%ntimimage==1) write(msg,'(2a)') trim(msg),' FOR 1 TIME STEP' 482 if (dtset%ntimimage >1) write(msg,'(2a,i5)') trim(msg),' - TIME STEP ',itimimage 483 if (dtset%prtvolimg<2) then 484 write(msg,'(3a)') trim(msg),ch10,& 485 & '================================================================================' 486 end if 487 call wrtout(ab_out ,msg,'COLL') 488 call wrtout(std_out,msg,'PERS') 489 490 call yaml_iterstart('timimage', itimimage, dev_null, dtset%use_yaml) 491 end if 492 493 if (dtset%use_yaml == 1) call yaml_iterstart('timimage', itimimage, ab_out, dtset%use_yaml) 494 495 call timab(1204,2,tsec) 496 497 ! Loop on the dynamical images 498 idynimage=0 499 do iimage=1,nimage 500 501 ii=mpi_enreg%my_imgtab(iimage) 502 if (dtset%dynimage(ii)==1) idynimage=idynimage+1 503 504 ! Compute static image only at first time step 505 if (dtset%dynimage(ii)==1.or.compute_all_images) then 506 507 call timab(1205,1,tsec) 508 509 ! Change file names according to image index (if nimage>1) 510 if (dtset%nimage>1) then 511 call dtfil_init(dtfil,dtset,filnam,filstat,idtset,jdtset,mpi_enreg,ndtset,& 512 & image_index=ii) 513 if (itimimage>1) then 514 dtfil%ireadwf=0;dtfil%ireadden=0;dtfil%ireadkden=0 515 end if 516 call yaml_iterstart('image', iimage, dev_null, 0) 517 end if 518 519 if (dtset%use_yaml == 1) call yaml_iterstart('image', iimage, ab_out, dtset%use_yaml) 520 521 ! Redefine output units 522 call localwrfile(mpi_enreg%comm_cell,ii,dtset%nimage,mpi_enreg%paral_img,dtset%prtvolimg) 523 524 ! Print title for image 525 if (dtset%nimage>1.and.(dtset%prtvolimg==0.or.do_write_log)) then 526 if (ii==1) write(msg,'(a)' ) ch10 527 if (ii >1) write(msg,'(2a)') ch10,ch10 528 write(msg,'(6a,i4,a,i4,3a)') trim(msg),& 529 '--------------------------------------------------------------------------------',ch10,& 530 ' ',trim(imagealgo_str(dtset%imgmov)),' - CELL # ',ii,'/',dtset%nimage,ch10,& 531 '--------------------------------------------------------------------------------',ch10 532 if (dtset%prtvolimg==0) call wrtout(ab_out ,msg,'COLL') 533 if (do_write_log) call wrtout(std_out,msg,'PERS') 534 end if 535 536 acell(:) =res_img(iimage)%acell(:) 537 rprim(:,:) =res_img(iimage)%rprim(:,:) 538 vel(:,:) =res_img(iimage)%vel(:,:) 539 vel_cell(:,:)=res_img(iimage)%vel_cell(:,:) 540 xred(:,:) =res_img(iimage)%xred(:,:) 541 occ(:) =occ_img(:,iimage) 542 543 call args_gs_init(args_gs, & 544 & res_img(iimage)%amu(:),dtset%cellcharge(ii),res_img(iimage)%mixalch(:,:),& 545 & dtset%dmatpawu(:,:,:,:,ii),dtset%upawu(:,ii),dtset%jpawu(:,ii),& 546 & dtset%rprimd_orig(:,:,ii)) 547 548 call timab(1205,2,tsec) 549 550 call gstate(args_gs,acell,codvsn,cpui,dtfil,dtset,iexit,scf_initialized(iimage),itimimage_gstate,& 551 & mpi_enreg,npwtot,occ,pawang,pawrad,pawtab,psps,& 552 & res_img(iimage)%results_gs,& 553 & rprim,scf_history(iimage),vel,vel_cell,wvl,xred) 554 itimimage_gstate=itimimage_gstate+1 555 556 call timab(1206,1,tsec) 557 558 call args_gs_free(args_gs) 559 560 if (dtset%dynimage(ii)==1) then 561 res_img(iimage)%acell(:) =acell(:) 562 res_img(iimage)%rprim(:,:) =rprim(:,:) 563 res_img(iimage)%vel(:,:) =vel(:,:) 564 res_img(iimage)%vel_cell(:,:)=vel_cell(:,:) 565 res_img(iimage)%xred(:,:) =xred(:,:) 566 occ_img(:,iimage) =occ(:) 567 end if 568 569 ! check change of rprim and rewriting in hist 570 ! check change of xred and rewriting in hist 571 572 ! Close output units ; restore defaults 573 call localredirect(mpi_enreg%comm_cell,mpi_enreg%comm_world,dtset%nimage,mpi_enreg%paral_img,dtset%prtvolimg) 574 call timab(1206,2,tsec) 575 576 else if (itimimage>1) then ! For static images, simply copy one time step to the other 577 itimimage_prev=itimimage_eff-1 578 if (itimimage_prev<1) itimimage_prev=ntimimage_stored 579 call copy_results_img(results_img(iimage,itimimage_prev), & 580 & results_img(iimage,itimimage_eff )) 581 end if 582 583 ! Store results in hist datastructure 584 if (use_hist) then 585 ih=hist(iimage)%ihist 586 call mkrdim(res_img(iimage)%acell(:),res_img(iimage)%rprim(:,:),rprimd) 587 call var2hist(res_img(iimage)%acell(:),hist(iimage),dtset%natom,& 588 & rprimd,res_img(iimage)%xred(:,:),.FALSE.) 589 call vel2hist(amass(:,iimage),hist(iimage),res_img(iimage)%vel(:,:),& 590 & res_img(iimage)%vel_cell(:,:)) 591 hist(iimage)%fcart(:,:,ih)=res_img(iimage)%results_gs%fcart(:,:) 592 hist(iimage)%strten(:,ih)=res_img(iimage)%results_gs%strten(:) 593 hist(iimage)%etot(ih)=res_img(iimage)%results_gs%etotal 594 hist(iimage)%entropy(ih)=res_img(iimage)%results_gs%energies%entropy 595 hist(iimage)%time(ih)=real(itimimage,kind=dp)*dtion 596 end if 597 598 end do ! iimage 599 600 if(mpi_enreg%paral_img==1)then 601 call timab(1208,1,tsec) 602 call xmpi_barrier(mpi_enreg%comm_img) 603 call timab(1208,2,tsec) 604 end if 605 606 call timab(1209,1,tsec) 607 608 ! Output when images are used 609 if (dtset%nimage>1) then 610 ! === 1st option: reduced outputs === 611 if (dtset%prtvolimg>0) then 612 call prtimg(dtset%dynimage,imagealgo_str(dtset%imgmov),dtset%imgmov,ab_out,& 613 & mpi_enreg,nimage,dtset%nimage,compute_all_images,dtset%prtvolimg,res_img) 614 end if 615 end if 616 617 ! Manage log files when images are used 618 call localrdfile(mpi_enreg%comm_img,mpi_enreg%comm_world,compute_all_images,& 619 & dtset%nimage,mpi_enreg%paral_img,dtset%prtvolimg,dyn=dtset%dynimage) 620 621 ! Write hist datastructure in HIST file 622 #if defined HAVE_NETCDF 623 if (use_hist.and.mpi_enreg%me_cell==0) then 624 ifirst=merge(0,1,itimimage>1) 625 call write_md_hist_img(hist,hist_filename,ifirst,itimimage,dtset%natom,dtset%ntypat,& 626 & dtset%typat,amu_img(:,1),dtset%znucl,dtion,& 627 & nimage=dtset%nimage,imgmov=dtset%imgmov,mdtemp=dtset%mdtemp,comm_img=mpi_enreg%comm_img,& 628 & imgtab=mpi_enreg%my_imgtab) 629 end if 630 #endif 631 632 ! TESTS WHETHER ONE CONTINUES THE LOOP 633 ! Here we calculate the change in energy, and exit if delta_energy < tolimg 634 delta_energy=zero 635 ! Doesn't check convergence in case of PIMD 636 check_conv=((.not.is_pimd).and.itimimage>1) 637 ! In case of 4th-order Runge-Kutta, does check convergence every 4 steps 638 if (dtset%imgmov==2.and.mep_param%mep_solver==4) then 639 check_conv=(mod(itimimage,4)==0.and.itimimage>4) 640 end if 641 if (check_conv) then 642 do idynimage=1,ndynimage 643 iimage=list_dynimage(idynimage) 644 delta_energy=delta_energy & 645 & +abs(results_img(iimage,itimimage)%results_gs%etotal & 646 & -results_img(iimage,itimimage-idelta)%results_gs%etotal) 647 end do 648 if (mpi_enreg%paral_img==1) then 649 call xmpi_sum(delta_energy,mpi_enreg%comm_img,ierr) 650 end if 651 delta_energy=delta_energy/dtset%ndynimage 652 if (delta_energy<dtset%tolimg) then 653 if (dtset%prtvolimg<2) then 654 write(msg,'(5a,i5,6a,es11.3,a,es11.3,2a)') ch10,ch10,& 655 & '================================================================================',ch10,& 656 & ' At time step ',itimimage,ch10,& 657 & ' ',trim(imagealgo_str(dtset%imgmov)),' has reached energy convergence',ch10,& 658 & ' with Average[Abs(Etotal(t)-Etotal(t-dt))]=',delta_energy,'<tolimg=',dtset%tolimg,ch10,& 659 & '================================================================================' 660 else 661 write(msg,'(4a,i5,6a,es11.3,a,es11.3)') ch10,& 662 & '--------------------------------------------------------------------------------',ch10,& 663 & ' At time step ',itimimage,ch10,& 664 & ' ',trim(imagealgo_str(dtset%imgmov)),' has reached energy convergence',ch10,& 665 & ' with Average[Abs(Etotal(t)-Etotal(t-dt))]=',delta_energy,'<tolimg=',dtset%tolimg 666 end if 667 call wrtout([std_out, ab_out] ,msg,'COLL') 668 call timab(1209,2,tsec) ! This is the first place where counter 1209 is stopped. 669 exit ! exit itimimage 670 end if 671 end if 672 673 !Temporary statement 674 110 continue 675 676 ! Dont call the predictor at last time step 677 if (itimimage>=ntimimage_max) call_predictor=(call_predictor.and.is_pimd) 678 679 ! Predict the next value of the images 680 if (call_predictor) then 681 call predictimg(delta_energy,imagealgo_str(dtset%imgmov),dtset%imgmov,itimimage,& 682 & itimimage_eff,list_dynimage,ga_param,mep_param,mpi_enreg,m1geo_param,dtset%natom,ndynimage,& 683 & nimage,dtset%nimage,ntimimage_stored,pimd_param,dtset%prtvolimg,results_img) 684 end if 685 686 ! Increment indexes 687 if (itimimage>=ntimimage_max) exit 688 itimimage_eff=itimimage_eff+1;if (itimimage_eff>ntimimage_stored) itimimage_eff=1 689 if (use_hist) then 690 do iimage=1,nimage 691 hist(iimage)%ihist=hist(iimage)%ihist+1 692 end do 693 end if 694 695 call timab(1209,2,tsec) ! This is the second place where counter 1209 is stopped. 696 697 end do ! itimimage 698 !----------------------------------------------------------------------------------------- 699 700 call timab(1210,1,tsec) 701 702 !Copy the results of the computation in the appropriate arguments of the routine 703 do iimage=1,nimage 704 ii=mpi_enreg%my_imgtab(iimage) 705 if (dtset%dynimage(ii)==1) then 706 acell_img(:,iimage) =results_img(iimage,itimimage_eff)%acell(:) 707 amu_img(:,iimage) =results_img(iimage,itimimage_eff)%amu(:) 708 mixalch_img(:,:,iimage) =results_img(iimage,itimimage_eff)%mixalch(:,:) 709 rprim_img(:,:,iimage) =results_img(iimage,itimimage_eff)%rprim(:,:) 710 vel_img(:,:,iimage) =results_img(iimage,itimimage_eff)%vel(:,:) 711 vel_cell_img(:,:,iimage)=results_img(iimage,itimimage_eff)%vel_cell(:,:) 712 xred_img(:,:,iimage) =results_img(iimage,itimimage_eff)%xred(:,:) 713 etotal_img(iimage) =results_img(iimage,itimimage_eff)%results_gs%etotal 714 fcart_img(:,:,iimage) =results_img(iimage,itimimage_eff)%results_gs%fcart(:,:) 715 gred_img(:,:,iimage) =results_img(iimage,itimimage_eff)%results_gs%gred(:,:) 716 intgres_img(:,:,iimage) =results_img(iimage,itimimage_eff)%results_gs%intgres(:,:) 717 strten_img(:,iimage) =results_img(iimage,itimimage_eff)%results_gs%strten(:) 718 else if (compute_static_images) then 719 etotal_img(iimage) =results_img(iimage,1)%results_gs%etotal 720 fcart_img(:,:,iimage) =results_img(iimage,1)%results_gs%fcart(:,:) 721 gred_img(:,:,iimage) =results_img(iimage,1)%results_gs%gred(:,:) 722 intgres_img(:,:,iimage)=results_img(iimage,1)%results_gs%intgres(:,:) 723 strten_img(:,iimage) =results_img(iimage,1)%results_gs%strten(:) 724 end if 725 end do 726 727 !When parallelizattion over images is activated, has to sum number of warnings 728 !and comments written in log file 729 if (mpi_enreg%paral_img==1.and.mpi_enreg%comm_cell==0) then 730 call specialmsg_mpisum(mpi_enreg%comm_img) 731 call libpaw_spmsg_mpisum(mpi_enreg%comm_img) 732 end if 733 734 735 !Final deallocations 736 737 !This call is needed to free internal storages in different routines (prec_simple, pred_bfgs ...) 738 if(dtset%imgmov==6)then 739 m1geo_param%iexit=1 740 call predictimg(delta_energy,imagealgo_str(dtset%imgmov),dtset%imgmov,itimimage,& 741 & itimimage_eff,list_dynimage,ga_param,mep_param,mpi_enreg,m1geo_param,dtset%natom,ndynimage,& 742 & nimage,dtset%nimage,ntimimage_stored,pimd_param,dtset%prtvolimg,results_img) 743 endif 744 745 ABI_FREE(occ) 746 ABI_FREE(vel) 747 ABI_FREE(xred) 748 ABI_FREE(list_dynimage) 749 750 if (allocated(amass)) then 751 ABI_FREE(amass) 752 end if 753 754 do itimimage=1,ntimimage_stored 755 call destroy_results_img(results_img(:,itimimage)) 756 end do 757 ABI_FREE(results_img) 758 do iimage=1,nimage 759 call scf_history_free(scf_history(iimage)) 760 end do 761 ABI_FREE(scf_history) 762 ABI_FREE(scf_initialized) 763 if (allocated(hist_prev)) then 764 call abihist_free(hist_prev) 765 ABI_FREE(hist_prev) 766 end if 767 if (allocated(hist)) then 768 call abihist_free(hist) 769 ABI_FREE(hist) 770 end if 771 772 call mep_destroy(mep_param) 773 call ga_destroy(ga_param) 774 call m1geo_destroy(m1geo_param) 775 call pimd_destroy(pimd_param) 776 777 call timab(1210,2,tsec) 778 call timab(1200,2,tsec) 779 780 DBG_EXIT("COLL") 781 782 end subroutine gstateimg
ABINIT/m_gstateimg [ Modules ]
NAME
m_gstateimg
FUNCTION
COPYRIGHT
Copyright (C) 1998-2024 ABINIT group (XG, AR, GG, MT) 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_gstateimg 22 23 use defs_basis 24 use defs_wvltypes 25 use defs_rectypes 26 use m_abicore 27 use m_abihist 28 use m_mep 29 use m_ga 30 use m_use_ga 31 use m_pimd 32 use m_xmpi 33 use m_errors 34 use m_rec 35 use m_args_gs 36 use m_results_gs 37 use m_results_img 38 use m_scf_history 39 use m_io_redirect 40 use m_m1geo 41 use m_abimover 42 use m_yaml 43 use m_dtfil 44 45 use defs_datatypes, only : pseudopotential_type 46 use defs_abitypes, only : MPI_type 47 use m_time, only : timab 48 use m_geometry, only : mkradim, mkrdim, fcart2gred, xred2xcart, metric 49 use m_specialmsg, only : specialmsg_mpisum 50 use m_libpaw_tools, only : libpaw_spmsg_mpisum 51 use m_pawang, only : pawang_type 52 use m_pawrad, only : pawrad_type 53 use m_pawtab, only : pawtab_type 54 use m_gstate, only : gstate 55 use m_predtk, only : prtxvf 56 use m_precpred_1geo, only : precpred_1geo 57 use m_pred_simple, only : prec_simple 58 59 #if defined HAVE_BIGDFT 60 use BigDFT_API, only: mpi_environment_set 61 #endif 62 63 implicit none 64 65 private
ABINIT/move_1geo [ Functions ]
NAME
move_1geo
FUNCTION
This subroutine uses the forces, stresses and other results obtained for several images with one, common, geometry, weight them to deliver averaged forces, stresses, etc, and uses these to predict the next common geometry. All images must be dynamical. WARNING : at present, only forces are used, to change atomic positions. No change of cell geometry. Since this is not the PIMD algorithm, suppose ntimimage_stored=ntimimage, and itimimage=itimimage_eff.
INPUTS
itimimage_eff=time index in the history nimage=number of images ntimimage_stored=number of time steps stored in the history mpi_enreg=MPI-parallelisation information m1geo_param=parameters for the 1geo algorithms
OUTPUT
SIDE EFFECTS
results_img(ntimimage_stored,nimage)=datastructure that holds the history of previous computations.
results_img(:,:)%acell(3)
at input, history of the values of acell for all images
at output, the predicted values of acell for all images
results_img(:,:)%results_gs
at input, history of the values of energies and forces for all images
results_img(:,:)%rprim(3,3)
at input, history of the values of rprim for all images
at output, the predicted values of rprim for all images
results_img(:,:)%vel(3,natom)
at input, history of the values of vel for all images
at output, the predicted values of vel for all images
results_img(:,:)%vel_cell(3,3)
at input, history of the values of vel_cell for all images
at output, the predicted values of vel_cell for all images
results_img(:,:)%xred(3,natom)
at input, history of the values of xred for all images
at output, the predicted values of xred for all images
SOURCE
1229 subroutine move_1geo(itimimage_eff,m1geo_param,mpi_enreg,nimage,nimage_tot,ntimimage_stored,results_img) 1230 1231 !Arguments ------------------------------------ 1232 !scalars 1233 integer,intent(in) :: itimimage_eff,nimage,nimage_tot,ntimimage_stored 1234 type(MPI_type),intent(in) :: mpi_enreg 1235 type(m1geo_type),intent(inout) :: m1geo_param 1236 !arrays 1237 type(results_img_type),target,intent(inout) :: results_img(nimage,ntimimage_stored) 1238 1239 !Local variables------------------------------- 1240 !scalars 1241 integer :: ihist,iimage,natom,next_itimimage,nspden,nsppol 1242 !integer :: iatom 1243 real(dp) :: deltae,diffor,etotal,entropy,fermie,res2,residm 1244 logical :: test_img 1245 type(results_gs_type) :: results_gs_lincomb 1246 !arrays 1247 real(dp) :: acell(3),rprim(3,3),rprimd(3,3),strten(6),vel_cell(3,3) 1248 real(dp),allocatable :: fcart(:,:),vel(:,:),xred(:,:) 1249 logical :: DEBUG=.FALSE. 1250 type(results_img_type),pointer :: resimg_all(:) 1251 1252 ! ************************************************************************* 1253 1254 natom=m1geo_param%ab_mover%natom 1255 ihist=m1geo_param%hist_1geo%ihist 1256 1257 ABI_MALLOC(fcart,(3,natom)) 1258 ABI_MALLOC(vel,(3,natom)) 1259 ABI_MALLOC(xred,(3,natom)) 1260 1261 !Of course, assume that the geometry parameters are the same for all images, so take them from the first one. 1262 xred(:,:) =results_img(1,itimimage_eff)%xred(:,:) 1263 acell(:) =results_img(1,itimimage_eff)%acell(:) 1264 rprim(:,:) =results_img(1,itimimage_eff)%rprim(:,:) 1265 vel(:,:) =results_img(1,itimimage_eff)%vel(:,:) 1266 vel_cell(:,:)=results_img(1,itimimage_eff)%vel_cell(:,:) 1267 1268 call mkrdim(acell,rprim,rprimd) 1269 1270 !Fill history with the values of xred, acell and rprimd 1271 call var2hist(acell,m1geo_param%hist_1geo,natom,rprimd,xred,DEBUG) 1272 1273 !Fill history with velocities and ionic kinetic energy 1274 call vel2hist(m1geo_param%ab_mover%amass,m1geo_param%hist_1geo,vel,vel_cell) 1275 m1geo_param%hist_1geo%time(ihist)=zero 1276 1277 !In case of image parallelism, collect results accross processors 1278 test_img=(nimage_tot/=1.and.mpi_enreg%paral_img==1) 1279 if (test_img) then 1280 ABI_MALLOC(resimg_all,(nimage_tot)) 1281 call gather_results_img(mpi_enreg,results_img(1:nimage,itimimage_eff),resimg_all,& 1282 & allgather=.true.,only_one_per_img=.false.) 1283 else 1284 resimg_all => results_img(:,itimimage_eff) 1285 end if 1286 1287 !Compute energy, entropy, fermie, forces and stresses for the 1geo : take the weighted average. 1288 !Compute maximum of deltae,diffor,res2,residm 1289 etotal=zero 1290 entropy=zero 1291 fermie=zero 1292 fcart(:,:)=zero 1293 strten(:)=zero 1294 deltae=zero 1295 diffor=zero 1296 res2=zero 1297 residm=zero 1298 1299 do iimage=1,nimage_tot 1300 etotal=etotal+resimg_all(iimage)%results_gs%etotal*m1geo_param%mixesimgf(iimage) 1301 entropy=entropy+resimg_all(iimage)%results_gs%entropy*m1geo_param%mixesimgf(iimage) 1302 fermie=fermie+resimg_all(iimage)%results_gs%fermie*m1geo_param%mixesimgf(iimage) 1303 fcart(:,:)=fcart(:,:)+resimg_all(iimage)%results_gs%fcart(:,:)*m1geo_param%mixesimgf(iimage) 1304 strten(:) =strten(:) +resimg_all(iimage)%results_gs%strten(:)*m1geo_param%mixesimgf(iimage) 1305 if( deltae<resimg_all(iimage)%results_gs%deltae ) deltae=resimg_all(iimage)%results_gs%deltae 1306 if( diffor<resimg_all(iimage)%results_gs%diffor ) diffor=resimg_all(iimage)%results_gs%diffor 1307 if( res2<resimg_all(iimage)%results_gs%res2 ) res2=resimg_all(iimage)%results_gs%res2 1308 if( residm<resimg_all(iimage)%results_gs%residm ) residm=resimg_all(iimage)%results_gs%residm 1309 enddo 1310 1311 !Set up a results_gs datastructure with the linear combination of images 1312 nspden=resimg_all(1)%results_gs%nspden 1313 nsppol=resimg_all(1)%results_gs%nsppol 1314 call init_results_gs(natom,nspden,nsppol,results_gs_lincomb) 1315 call copy_results_gs(resimg_all(1)%results_gs,results_gs_lincomb) 1316 results_gs_lincomb%etotal=etotal 1317 results_gs_lincomb%entropy=entropy 1318 results_gs_lincomb%fermie=fermie 1319 results_gs_lincomb%fcart=fcart 1320 results_gs_lincomb%strten=strten 1321 results_gs_lincomb%deltae=deltae 1322 results_gs_lincomb%diffor=diffor 1323 results_gs_lincomb%res2=res2 1324 results_gs_lincomb%residm=residm 1325 1326 !From now on, all procs contain the same information about the geometry, etotal, forces, stress, etc. 1327 !Nothing more needs to be transmitted, and resimg_all is not needed anymore. 1328 if (test_img) then 1329 call destroy_results_img(resimg_all) 1330 ABI_FREE(resimg_all) 1331 end if 1332 nullify(resimg_all) 1333 1334 !Echo result_gs_lincomb 1335 call results_gs_lincomb%yaml_write(ab_out, info="Linear combination of ground state results") 1336 1337 !Destroy result_gs_lincomb 1338 call destroy_results_gs(results_gs_lincomb) 1339 1340 !Store fcart and strten in hist_1geo 1341 m1geo_param%hist_1geo%fcart(:,:,ihist)=fcart(:,:) 1342 m1geo_param%hist_1geo%strten(:,ihist) =strten(:) 1343 1344 !Store them in ab_xfh 1345 !THIS IS TO BE DONE ! 1346 1347 !Compute new atomic positions and cell characteristics in the single geometry 1348 call precpred_1geo(m1geo_param%ab_mover,& 1349 & m1geo_param%ab_xfh_1geo,& 1350 & m1geo_param%ab_mover%amu_curr,& 1351 & m1geo_param%deloc,& 1352 & m1geo_param%dt_chkdilatmx,& 1353 & mpi_enreg%comm_cell,& 1354 & m1geo_param%dilatmx,& 1355 & m1geo_param%filnam_ds4,& 1356 & m1geo_param%hist_1geo,& 1357 & m1geo_param%hmctt,& 1358 & m1geo_param%icycle,& 1359 & m1geo_param%iexit,& 1360 & itimimage_eff,& 1361 & m1geo_param%mttk_vars,& 1362 & m1geo_param%nctime,& 1363 & m1geo_param%ncycle,& 1364 & m1geo_param%nerr_dilatmx,& 1365 & m1geo_param%npsp,& 1366 & m1geo_param%ntime,& 1367 & m1geo_param%rprimd_orig,& 1368 & m1geo_param%skipcycle,& 1369 & m1geo_param%usewvl) 1370 1371 !Retrieve the new positions, cell parameters [and velocities ?!] 1372 call hist2var(acell,m1geo_param%hist_1geo,natom,rprimd,xred,DEBUG) 1373 1374 !Store acell, rprim, xred and vel for the new iteration if relevant 1375 if(m1geo_param%iexit==0)then 1376 next_itimimage=itimimage_eff+1 1377 if (next_itimimage>ntimimage_stored)then 1378 ABI_ERROR('next_itimimage>ntimimage_stored') 1379 endif 1380 do iimage=1,nimage 1381 results_img(iimage,next_itimimage)%xred(:,:) =xred(:,:) 1382 results_img(iimage,next_itimimage)%acell(:) =acell(:) 1383 results_img(iimage,next_itimimage)%rprim(:,:) =rprim(:,:) 1384 ! WARNING : Should also store vel and vel_cell of course ... 1385 ! results_img(iimage,next_itimimage)%vel(:,:) =vel(:,:) 1386 ! results_img(iimage,next_itimimage)%vel_cell(:,:)=vel_cell(:,:) 1387 end do 1388 endif 1389 ABI_FREE(fcart) 1390 ABI_FREE(vel) 1391 ABI_FREE(xred) 1392 1393 end subroutine move_1geo
ABINIT/predict_copy [ Functions ]
NAME
predict_copy
FUNCTION
Given the past history of images, predict the new set of images. Here, simple copy of the previous image.
INPUTS
itimimage_eff=time index in the history
list_dynimage(nimage)=list of dynamical images. The non-dynamical ones will not change.
Example : in the NEB of string method, one expect the two end images to be fixed.
ndynimage=number of dynamical images
nimage=number of images
ntimimage_stored=number of time steps stored in the history
OUTPUT
SIDE EFFECTS
results_img(ntimimage_stored,nimage)=datastructure that holds the history of previous computations.
results_img(:,:)%acell(3)
at input, history of the values of acell for all images
at output, the predicted values of acell for all images
results_img(:,:)%results_gs
at input, history of the values of energies and forces for all images
results_img(:,:)%rprim(3,3)
at input, history of the values of rprim for all images
at output, the predicted values of rprim for all images
results_img(:,:)%vel(3,natom)
at input, history of the values of vel for all images
at output, the predicted values of vel for all images
results_img(:,:)%vel_cell(3,3)
at input, history of the values of vel_cell for all images
at output, the predicted values of vel_cell for all images
results_img(:,:)%xred(3,natom)
at input, history of the values of xred for all images
at output, the predicted values of xred for all images
SOURCE
1153 subroutine predict_copy(itimimage_eff,list_dynimage,ndynimage,nimage,& 1154 & ntimimage_stored,results_img) 1155 1156 !Arguments ------------------------------------ 1157 !scalars 1158 integer,intent(in) :: itimimage_eff,ndynimage,nimage,ntimimage_stored 1159 !arrays 1160 integer,intent(in) :: list_dynimage(ndynimage) 1161 type(results_img_type),intent(inout) :: results_img(nimage,ntimimage_stored) 1162 1163 !Local variables------------------------------- 1164 !scalars 1165 integer :: idynimage,iimage,next_itimimage 1166 1167 ! ************************************************************************* 1168 1169 next_itimimage=itimimage_eff+1 1170 if (next_itimimage>ntimimage_stored) next_itimimage=1 1171 1172 do idynimage=1,ndynimage 1173 1174 iimage=list_dynimage(idynimage) 1175 1176 results_img(iimage,next_itimimage)%acell(:) =results_img(iimage,itimimage_eff)%acell(:) 1177 results_img(iimage,next_itimimage)%rprim(:,:) =results_img(iimage,itimimage_eff)%rprim(:,:) 1178 results_img(iimage,next_itimimage)%vel(:,:) =results_img(iimage,itimimage_eff)%vel(:,:) 1179 results_img(iimage,next_itimimage)%vel_cell(:,:)=results_img(iimage,itimimage_eff)%vel_cell(:,:) 1180 results_img(iimage,next_itimimage)%xred(:,:) =results_img(iimage,itimimage_eff)%xred(:,:) 1181 1182 end do ! idynimage 1183 1184 end subroutine predict_copy
ABINIT/predictimg [ Functions ]
NAME
predictimg
FUNCTION
Given the past history of images, predict the new set of images
INPUTS
deltae=averaged energy difference used to control convergence over images imagealgo_str=name of the algorithm (with images) used imgmov=gives the algorithm to be used for prediction of new set of images itimimage=time index for image propagation (itimimage+1 is to be predicted here) itimimage_eff=time index in the history list_dynimage(nimage)=list of dynamical images. The non-dynamical ones will not change. Example : in the NEB method, or in the string method, one expect the two end images to be fixed. mep_param=several parameters for Minimal Energy Path (MEP) search mpi_enreg=MPI-parallelisation information natom= number of atoms ndynimage=number of dynamical images nimage=number of images (treated by current proc) nimage_tot=total number of images ntimimage_stored=number of time steps stored in the history pimd_param=several parameters for Path-Integral MD prtvolimg=printing volume
OUTPUT
SIDE EFFECTS
results_img(ntimimage_stored,nimage)=datastructure that holds the history of previous computations.
results_img(:,:)%acell(3)
at input, history of the values of acell for all images
at output, the predicted values of acell for all images
results_img(:,:)%results_gs
at input, history of the values of energies and forces for all images
results_img(:,:)%rprim(3,3)
at input, history of the values of rprim for all images
at output, the predicted values of rprim for all images
results_img(:,:)%vel(3,natom)
at input, history of the values of vel for all images
at output, the predicted values of vel for all images
results_img(:,:)%vel_cell(3,3)
at input, history of the values of vel_cell for all images
at output, the predicted values of vel_cell for all images
results_img(:,:)%xred(3,natom)
at input, history of the values of xred for all images
at output, the predicted values of xred for all images
SOURCE
1006 subroutine predictimg(deltae,imagealgo_str,imgmov,itimimage,itimimage_eff,list_dynimage,& 1007 & ga_param,mep_param,mpi_enreg,m1geo_param,natom,ndynimage,nimage,nimage_tot,& 1008 & ntimimage_stored,pimd_param,prtvolimg,results_img) 1009 1010 use m_results_gs , only : results_gs_type 1011 use m_predict_neb, only : predict_neb 1012 use m_predict_steepest, only : predict_steepest 1013 use m_predict_pimd, only : predict_pimd 1014 use m_predict_string, only : predict_string 1015 1016 !Arguments ------------------------------------ 1017 !scalars 1018 integer,intent(in) :: imgmov,itimimage,itimimage_eff,natom,ndynimage 1019 integer,intent(in) :: nimage,nimage_tot,ntimimage_stored,prtvolimg 1020 character(len=60),intent(in) :: imagealgo_str 1021 real(dp),intent(in) :: deltae 1022 type(mep_type),intent(inout) :: mep_param 1023 type(m1geo_type),intent(inout) :: m1geo_param 1024 type(ga_type),intent(inout) :: ga_param 1025 type(pimd_type),intent(inout) :: pimd_param 1026 type(MPI_type),intent(in) :: mpi_enreg 1027 !arrays 1028 integer,intent(in) :: list_dynimage(ndynimage) 1029 type(results_img_type) :: results_img(nimage,ntimimage_stored) 1030 1031 !Local variables------------------------------- 1032 !scalars 1033 integer,save :: idum=5 1034 logical :: is_pimd 1035 character(len=500) :: msg 1036 1037 ! ************************************************************************* 1038 1039 is_pimd=(imgmov==9.or.imgmov==10.or.imgmov==13) 1040 1041 !Write convergence info 1042 write(msg,'(3a)') ch10,& 1043 & '------------------------------------------------------------',ch10 1044 if (prtvolimg<2) write(msg,'(5a)') trim(msg),' ',trim(imagealgo_str),':',ch10 1045 1046 !Specific case of 4th-order RK algorithm 1047 if (mep_param%mep_solver==4) then 1048 if (mod(itimimage,4)==0) then 1049 write(msg,'(4a)') trim(msg),& 1050 & ' Fourth-order Runge-Kutta algorithm - final step',ch10 1051 if (itimimage>4) write(msg,'(2a,es11.3,2a)') trim(msg),& 1052 & ' Average[Abs(Etotal(t)-Etotal(t-dt))]=',deltae,' Hartree',ch10 1053 write(msg,'(2a)') trim(msg),' Moving images of the cell...' 1054 else 1055 write(msg,'(2a,i1,2a)') trim(msg),& 1056 & ' Fourth-order Runge-Kutta algorithm - intermediate step ',mod(itimimage,4),ch10 1057 write(msg,'(2a)') trim(msg),' Computing new intermediate positions...' 1058 end if 1059 else if (is_pimd) then 1060 1061 ! PIMD 1062 write(msg,'(2a)') trim(msg),' Moving images of the cell...' 1063 else 1064 1065 ! Other cases 1066 if (itimimage>1) write(msg,'(2a,es11.3,2a)') trim(msg),& 1067 & ' Average[Abs(Etotal(t)-Etotal(t-dt))]=',deltae,' Hartree',ch10 1068 write(msg,'(2a)') trim(msg),' Moving images of the cell...' 1069 1070 end if 1071 1072 !Write the msg 1073 !Prevent writing if iexit==1, which at present only happens for imgmov==6 algo 1074 if(imgmov/=6 .or. m1geo_param%iexit==0) call wrtout([std_out, ab_out] ,msg) 1075 1076 select case(imgmov) 1077 1078 case(0) 1079 call predict_copy(itimimage_eff,list_dynimage,ndynimage,nimage,& 1080 & ntimimage_stored,results_img) 1081 1082 case(1) 1083 call predict_steepest(itimimage,itimimage_eff,list_dynimage,mep_param,natom,ndynimage,nimage,& 1084 & ntimimage_stored,results_img) 1085 1086 case(2) 1087 call predict_string(itimimage,itimimage_eff,list_dynimage,mep_param,mpi_enreg,natom,& 1088 & ndynimage,nimage,nimage_tot,ntimimage_stored,results_img) 1089 1090 case(4) 1091 call predict_ga(itimimage_eff,idum,ga_param,natom,nimage,& 1092 & ntimimage_stored,results_img) 1093 1094 case(5) 1095 call predict_neb(itimimage,itimimage_eff,list_dynimage,mep_param,mpi_enreg,natom,& 1096 & ndynimage,nimage,nimage_tot,ntimimage_stored,results_img) 1097 1098 case(6) 1099 call move_1geo(itimimage_eff,m1geo_param,mpi_enreg,nimage,nimage_tot,ntimimage_stored,results_img) 1100 1101 case(9, 10, 13) 1102 ! Path Integral Molecular Dynamics 1103 call predict_pimd(imgmov,itimimage,itimimage_eff,mpi_enreg,natom,nimage,nimage_tot,& 1104 & ntimimage_stored,pimd_param,prtvolimg,results_img) 1105 1106 case default 1107 1108 end select 1109 1110 end subroutine predictimg
ABINIT/prtimg [ Functions ]
NAME
prtimg
FUNCTION
Print out results obtained by as ground-state calculation of an image of the cell. The printing format is condensed in order to facilitate the reading.
INPUTS
dynimage(nimagetot)=flags defining static/dynamic state of images
imagealgo_str=name of the algorithm (with images) used
imgmov=index of algorithm (with images) used
iout=unit number for output
mpi_enreg=MPI-parallelisation information
nimage=number of images stored on current proc
nimage_tot=total number of images (should be dtset%nimage)
prt_all_images=true if all images have to be printed out (ignoring dynimage)
prtvolimg=printing volume for each image
<0 : nothing
0 : only a title
1 : energy, residuals, forces, stresses, velocities, atomic positions
2 : energy, residuals
resimg(nimage) <type(results_img_type)>=results of the ground-state computations
for all images treated by current proc
OUTPUT
(data written to unit iout)
SOURCE
816 subroutine prtimg(dynimage,imagealgo_str,imgmov,iout,mpi_enreg,nimage,nimage_tot,& 817 & prt_all_images,prtvolimg,resimg) 818 819 !Arguments ------------------------------------ 820 !scalars 821 integer,intent(in) :: nimage_tot,dynimage(nimage_tot),imgmov,iout,nimage,prtvolimg !vz_d 822 logical,intent(in) :: prt_all_images 823 character(len=60),intent(in) :: imagealgo_str 824 type(MPI_type),intent(in) :: mpi_enreg 825 !arrays 826 type(results_img_type),target,intent(inout) :: resimg(nimage) 827 828 !Local variables------------------------------- 829 !scalars 830 integer :: ii,prtvel 831 logical :: test_img 832 real(dp) :: ucvol_img 833 character(len=500) :: msg 834 !arrays 835 integer,allocatable :: iatfix_img(:,:) 836 real(dp),allocatable :: gmet_img(:,:),gprimd_img(:,:),rmet_img(:,:),xcart_img(:,:) 837 type(results_img_type),pointer :: resimg_all(:) 838 839 ! **************************************************************** 840 841 DBG_ENTER('COLL') 842 843 if (prtvolimg<=0) return 844 if (mpi_enreg%me_cell/=0) return 845 846 !Gather data 847 if (prtvolimg==1.or.prtvolimg==2) then 848 test_img=(nimage_tot/=1.and.mpi_enreg%paral_img==1) 849 if (test_img) then 850 if (mpi_enreg%me==0) then 851 ABI_MALLOC(resimg_all,(nimage_tot)) 852 end if 853 call gather_results_img(mpi_enreg,resimg,resimg_all,master=0,& 854 & allgather=.false.,only_one_per_img=.true.) 855 else 856 resimg_all => resimg 857 end if 858 end if 859 860 !===== First option for the printing volume === 861 if (prtvolimg==1.and.mpi_enreg%me==0) then 862 863 prtvel=0;if (imgmov==0.or.imgmov==9.or.imgmov==10.or.imgmov==13) prtvel=1 864 865 do ii=1,nimage_tot 866 if (dynimage(ii)==1.or.prt_all_images) then 867 868 ! Title 869 write(msg,'(6a,i4,a,i4,2a)') ch10,& 870 '----------------------------------------------------------------------',ch10,& 871 ' ',trim(imagealgo_str),' - CELL # ',ii,'/',nimage_tot,ch10,& 872 '----------------------------------------------------------------------' 873 call wrtout(iout,msg,'COLL') 874 875 ! Total energy 876 write(msg,'(2a,es20.12)') ch10,' Total energy for the cell [Ha]: ',resimg_all(ii)%results_gs%etotal 877 call wrtout(iout,msg,'COLL') 878 879 ! Residuals of the SCF cycle 880 write(msg,'(3a,4(a,es16.8,a))') ch10,& 881 ' Residuals from SCF cycle: ',ch10,& 882 ' Total energy difference =',resimg_all(ii)%results_gs%deltae,ch10,& 883 ' Maximal forces difference =',resimg_all(ii)%results_gs%diffor,ch10,& 884 ' Max. residual of wave-functions=',resimg_all(ii)%results_gs%residm,ch10,& 885 ' Density/potential residual (^2)=',resimg_all(ii)%results_gs%res2,ch10 886 call wrtout(iout,msg,'COLL') 887 888 ! Cell parameters 889 ABI_MALLOC(rmet_img,(3,3)) 890 ABI_MALLOC(gmet_img,(3,3)) 891 ABI_MALLOC(gprimd_img,(3,3)) 892 call metric(gmet_img,gprimd_img,iout,rmet_img,resimg_all(ii)%rprim,ucvol_img) 893 ABI_FREE(rmet_img) 894 ABI_FREE(gmet_img) 895 ABI_FREE(gprimd_img) 896 897 ! Positions, forces and velocities 898 ABI_MALLOC(iatfix_img,(3,resimg_all(ii)%natom)) 899 ABI_MALLOC(xcart_img,(3,resimg_all(ii)%natom)) 900 iatfix_img=0 901 call xred2xcart(resimg_all(ii)%natom,resimg_all(ii)%rprim,xcart_img,resimg_all(ii)%xred) 902 call prtxvf(resimg_all(ii)%results_gs%fcart,resimg_all(ii)%results_gs%gred,& 903 & iatfix_img,iout,resimg_all(ii)%natom,prtvel,& 904 & resimg_all(ii)%vel,xcart_img,resimg_all(ii)%xred) 905 ABI_FREE(iatfix_img) 906 ABI_FREE(xcart_img) 907 908 ! Stress tensor 909 write(msg, '(a,es12.4,a)' ) & 910 & '-Cartesian components of stress tensor (GPa) [Pressure=',& 911 & -(resimg_all(ii)%results_gs%strten(1)+resimg_all(ii)%results_gs%strten(2) & 912 & +resimg_all(ii)%results_gs%strten(3))*HaBohr3_GPa/three,' GPa]' 913 call wrtout(iout,msg,'COLL') 914 write(msg, '(2(a,1p,e16.8))' ) '- sigma(1 1)=',resimg_all(ii)%results_gs%strten(1)*HaBohr3_GPa,& 915 & ' sigma(3 2)=',resimg_all(ii)%results_gs%strten(4)*HaBohr3_GPa 916 call wrtout(iout,msg,'COLL') 917 write(msg, '(2(a,1p,e16.8))' ) '- sigma(2 2)=',resimg_all(ii)%results_gs%strten(2)*HaBohr3_GPa,& 918 & ' sigma(3 1)=',resimg_all(ii)%results_gs%strten(5)*HaBohr3_GPa 919 call wrtout(iout,msg,'COLL') 920 write(msg, '(2(a,1p,e16.8))' ) '- sigma(3 3)=',resimg_all(ii)%results_gs%strten(3)*HaBohr3_GPa,& 921 & ' sigma(2 1)=',resimg_all(ii)%results_gs%strten(6)*HaBohr3_GPa 922 call wrtout(iout,msg,'COLL') 923 end if 924 end do 925 end if 926 927 928 !===== 2nd option for the printing volume === 929 if (prtvolimg==2.and.mpi_enreg%me==0) then 930 write(msg,'(a,1x,a)') ch10,'Cell Total_energy[Ha] deltae diffor residm res2' 931 call wrtout(iout,msg,'COLL') 932 do ii=1,nimage_tot 933 if (dynimage(ii)==1.or.prt_all_images) then 934 write(msg,'(1x,i4,2x,es16.8,4(1x,es13.5))') & 935 & ii,resimg_all(ii)%results_gs%etotal,resimg_all(ii)%results_gs%deltae,& 936 & resimg_all(ii)%results_gs%diffor,resimg_all(ii)%results_gs%residm,& 937 & resimg_all(ii)%results_gs%res2 938 call wrtout(iout,msg,'COLL') 939 end if 940 end do 941 end if 942 943 !===== 944 if (prtvolimg==1.or.prtvolimg==2) then 945 if (test_img.and.mpi_enreg%me==0) then 946 call destroy_results_img(resimg_all) 947 ABI_FREE(resimg_all) 948 end if 949 nullify(resimg_all) 950 end if 951 952 DBG_EXIT('COLL') 953 954 end subroutine prtimg