TABLE OF CONTENTS
ABINIT/predict_ga [ Modules ]
NAME
predict_ga
FUNCTION
Given a given set of images, which represent a population, it predicts a new set of images. The implementation is based on a Genetic Algorithm idea, where the best fit candidates are passed to the next generation. Those are chosen from ga_opt_percent% best fit and (1-ga_opt_percent)% from Genetic rules
COPYRIGHT
Copyright (C) 2009-2022 ABINIT group (XG, AHR) 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
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. This is quite useful when ground states of the A and B states is known natom=dimension of vel_timimage and xred_timimage ndynimage=number of dynamical images nimage= population size 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
51 #if defined HAVE_CONFIG_H 52 #include "config.h" 53 #endif 54 55 #include "abi_common.h" 56 57 MODULE m_use_ga 58 59 use defs_basis 60 use m_abicore 61 use m_ga 62 use m_sort 63 64 use m_symfind, only : symfind, symanal, symlatt 65 use m_geometry, only : mkradim, mkrdim, metric, dist2 66 use m_results_img, only : results_img_type,gather_array_img 67 use m_numeric_tools, only : uniformrandom 68 69 implicit none 70 71 private 72 73 public :: predict_ga 74 75 CONTAINS 76 77 subroutine predict_ga(itimimage_eff,idum,ga_param,natom,nimage,& 78 & ntimimage_stored,results_img) 79 80 implicit none 81 82 !Arguments ------------------------------------ 83 !scalars 84 integer,intent(in) :: itimimage_eff,natom,nimage,ntimimage_stored 85 integer,intent(inout) :: idum 86 !arrays 87 type(results_img_type) :: results_img(nimage,ntimimage_stored) 88 type(ga_type),intent(inout) :: ga_param 89 90 !Local variables------------------------------ 91 !scalars 92 integer :: ii,jj,kk,itmp,iimage,indiv,oper,iparent1,iparent2,ndimen 93 integer :: next_itimimage,nsurvivor,nspinor 94 !character(len=500) :: message 95 !arrays 96 integer,allocatable :: ieperm(:),ihperm(:),ibgoptperm(:,:),ibgperm(:,:) 97 integer,allocatable :: zperm1(:),zperm2(:) 98 real(dp) :: gprimd(3,3),rmet(3,3),gmet(3,3) 99 real(dp) :: rprimdparent1(3,3),rprimdparent2(3,3) 100 real(dp) :: rprimson1(3,3),rprimson2(3,3) 101 real(dp) :: rprimdson1(3,3),rprimdson2(3,3) 102 real(dp) :: acellson1(3),acellson2(3) 103 real(dp) :: son1(3,natom),son2(3,natom) 104 real(dp) :: parent1(3,natom),parent2(3,natom) 105 ! store the energy and the enthalpy of all elements of the population 106 real(dp),allocatable :: etotal_img(:),enthalpy_img(:),bg_img(:,:),bg_opt_img(:,:) 107 real(dp),allocatable :: acell(:,:),acell_old(:,:),rprim(:,:,:),rprim_old(:,:,:),rprimd(:,:,:) 108 real(dp),allocatable :: fitness(:),zcoor1(:),zcoor2(:) 109 real(dp),allocatable :: coor(:,:),coor_old(:,:) 110 real(dp),allocatable :: vson1(:),vson2(:) 111 112 !real quantities 113 real(dp) :: denom,sumH,Hmin,Hmax,ucvol,rtmp(3,3) 114 115 ! ************************************************************************* 116 117 !DEBUG 118 !write(std_out,*)' MODULE predict_ga : enter ' 119 !ENDDEBUG 120 121 ! gen dimension 122 123 ndimen=3*natom 124 nspinor=results_img(nimage,itimimage_eff)%nsppol 125 126 ABI_MALLOC(coor,(ndimen,nimage)) 127 ABI_MALLOC(coor_old,(ndimen,nimage)) 128 ABI_MALLOC(acell,(3,nimage)) 129 ABI_MALLOC(acell_old,(3,nimage)) 130 ABI_MALLOC(rprim,(3,3,nimage)) 131 ABI_MALLOC(rprimd,(3,3,nimage)) 132 ABI_MALLOC(rprim_old,(3,3,nimage)) 133 134 ABI_MALLOC(zperm1,(natom)) 135 ABI_MALLOC(zperm2,(natom)) 136 ABI_MALLOC(ieperm,(nimage)) 137 ABI_MALLOC(ihperm,(nimage)) 138 ABI_MALLOC(ibgperm,(nspinor,nimage)) 139 ABI_MALLOC(ibgoptperm,(nspinor,nimage)) 140 141 ABI_MALLOC(etotal_img,(nimage)) 142 ABI_MALLOC(enthalpy_img,(nimage)) 143 144 !to store locally the calculated band gaps 145 ABI_MALLOC(bg_img,(nspinor,nimage)) 146 ABI_MALLOC(bg_opt_img,(nspinor,nimage)) 147 148 ABI_MALLOC(fitness,(nimage)) 149 150 ABI_MALLOC(zcoor1,(natom)) 151 ABI_MALLOC(zcoor2,(natom)) 152 ABI_MALLOC(vson1,(ndimen)) 153 ABI_MALLOC(vson2,(ndimen)) 154 155 call initialize_perm(ieperm,nimage) 156 call initialize_perm(ihperm,nimage) 157 do ii=1,nspinor 158 call initialize_perm(ibgperm(ii,:),nimage) 159 call initialize_perm(ibgoptperm(ii,:),nimage) 160 enddo 161 162 !from gs_results_image take energies, rprim and acell and build energy and enthalpy vectors reordered from larger to smaller 163 164 do iimage=1,nimage 165 etotal_img(iimage)=results_img(iimage,itimimage_eff)%results_gs%etotal 166 call convert_coortogen(results_img(iimage,itimimage_eff)%xred(:,:),coor_old(:,iimage),natom) 167 acell_old(:,iimage)=results_img(iimage,itimimage_eff)%acell(:) 168 rprim_old(:,:,iimage)=results_img(iimage,itimimage_eff)%rprim(:,:) 169 if (results_img(iimage,itimimage_eff)%results_gs%gaps(3,1) > 0.0_dp) then 170 bg_img(:,iimage)=results_img(iimage,itimimage_eff)%results_gs%gaps(1,:) 171 bg_opt_img(:,iimage)=results_img(iimage,itimimage_eff)%results_gs%gaps(2,:) 172 else 173 bg_img(:,iimage)=-100.0_dp 174 bg_opt_img(:,iimage)=-100.0_dp 175 endif 176 call mkrdim(acell_old(:,iimage),rprim_old(:,:,iimage),rprimd(:,:,iimage)) 177 call metric(gmet,gprimd,-1,rmet,rprimd(:,:,iimage),ucvol) 178 enthalpy_img(iimage)=etotal_img(iimage) & 179 & +sum(results_img(iimage,itimimage_eff)%results_gs%strten(1:3))*ucvol 180 enddo 181 182 !sort energies 183 184 call sort_dp(nimage,etotal_img,ieperm,tol9) 185 186 !sort enthalpies 187 188 call sort_dp(nimage,enthalpy_img,ihperm,tol9) 189 190 !sort band gaps 191 192 do ii=1,nspinor 193 call sort_dp(nimage,bg_img(ii,:),ibgperm(ii,:),tol9) 194 call sort_dp(nimage,bg_opt_img(ii,:),ibgoptperm(ii,:),tol9) 195 enddo 196 197 ! Fitness is calculated 198 199 sumH=0.d0 200 Hmin=minval(enthalpy_img) 201 Hmax=maxval(enthalpy_img) 202 fitness = zero 203 204 select case (ga_param%ga_fitness) 205 case ( 1 ) 206 207 ! Function weighted over the difference with respect to the minimum 208 209 do iimage=1,nimage 210 sumh=sumh+enthalpy_img(iimage); fitness(iimage)=sumh 211 enddo 212 denom=nimage*enthalpy_img(nimage)-sumh 213 do iimage=1,nimage 214 fitness(iimage)=(iimage*enthalpy_img(nimage)-fitness(iimage))/denom 215 enddo 216 217 case ( 2 ) 218 219 ! Botzmann like, using the enthalpy difference. 220 221 do iimage=1,nimage 222 fitness(iimage)=exp(-one*(enthalpy_img(iimage)-enthalpy_img(1))) 223 sumH = sumH + fitness(iimage) 224 enddo 225 fitness = fitness/sumH 226 do iimage=2,nimage 227 fitness(iimage)=fitness(iimage)+fitness(iimage-1) 228 enddo 229 230 case ( 3 ) 231 232 ! weighted ove the position in the ordered list, with probability 1/i (non uniform) 233 234 do iimage=1,nimage 235 fitness(iimage)=one/float(iimage) 236 sumH = sumH + fitness(iimage) 237 enddo 238 fitness = fitness/sumH 239 do iimage=2,nimage 240 fitness(iimage)=fitness(iimage)+fitness(iimage-1) 241 enddo 242 243 end select 244 245 ! do a single GA boocle 246 247 indiv=0 248 249 ! Selection over the best ga_opt_percent of the population 250 251 nsurvivor=int(ga_param%ga_opt_percent*nimage) 252 253 if (nsurvivor < one) nsurvivor=1 254 255 ! pass coordinates,rprim, acell of survivors to next generation 256 257 do iimage=1,nsurvivor 258 indiv=indiv+1 259 coor(:,iimage)=coor_old(:,ihperm(iimage)) 260 acell(:,iimage)=acell_old(:,ihperm(iimage)) 261 rprim(:,:,iimage)=rprim_old(:,:,ihperm(iimage)) 262 enddo 263 264 ! complete the number of individuals of the generation by choosing them through GA 265 266 do while(indiv<nimage) 267 268 ! ga_n_rules corresponds to the number of chosen Genetic rules 269 270 oper=ga_param%ga_rules(int(ga_param%ga_n_rules*uniformrandom(idum)+1)) 271 272 select case(oper) 273 case(1) ! cut and splice elements after a randomization 274 iparent1=choosefather(fitness,nimage,idum) 275 iparent2=choosefather(fitness,nimage,idum) 276 call convert_gentocoor(parent1,coor_old(:,ihperm(iparent1)),natom) 277 call convert_gentocoor(parent2,coor_old(:,ihperm(iparent2)),natom) 278 ! randomize the cell: random rotation and translation 279 call randomize_parent(parent1,natom,idum) 280 call randomize_parent(parent2,natom,idum) 281 ! choose direction of cutting plane 282 itmp = int(3*uniformrandom(idum)+1) 283 ! order coordinates from small to large along that random direction axis of both parents 284 zcoor1(:)=parent1(itmp,:) 285 zcoor2(:)=parent2(itmp,:) 286 call initialize_perm(zperm1,natom) 287 call initialize_perm(zperm2,natom) 288 call sort_dp(natom,zcoor1,zperm1,tol9) 289 call sort_dp(natom,zcoor2,zperm2,tol9) 290 ! choose the atom position to take the cut and take atoms below from one parent and above from the other parent 291 itmp=int(natom*uniformrandom(idum)+1) 292 do ii=1,itmp 293 son1(:,ii)=parent1(:,zperm1(ii)) 294 son2(:,ii)=parent2(:,zperm2(ii)) 295 enddo 296 do ii=itmp+1,natom 297 son1(:,ii)=parent2(:,zperm2(ii)) 298 son2(:,ii)=parent1(:,zperm1(ii)) 299 enddo 300 ! random combination of rprimd from parents 301 call mkrdim(acell_old(:,ihperm(iparent1)),rprim_old(:,:,ihperm(iparent1)),rprimdparent1) 302 call mkrdim(acell_old(:,ihperm(iparent2)),rprim_old(:,:,ihperm(iparent2)),rprimdparent2) 303 do ii=1,3 304 do jj=1,3 305 rtmp(ii,jj)=uniformrandom(idum) 306 enddo 307 enddo 308 rprimdson1=(one-rtmp)*rprimdparent1+rtmp*rprimdparent2 309 do ii=1,3 310 do jj=1,3 311 rtmp(ii,jj)=uniformrandom(idum) 312 enddo 313 enddo 314 rprimdson2=(one-rtmp)*rprimdparent1+rtmp*rprimdparent2 315 !create acell and rprim from rprimd of springoffs 316 call mkradim(acellson1,rprimson1,rprimdson1) 317 call mkradim(acellson2,rprimson2,rprimdson2) 318 ! check distances of atoms of every springoff 319 if (checkatomicdist(natom,son1,rprimdson1)==0 .and. indiv<nimage) then 320 indiv=indiv+1 321 !if any fix coordinate restore the parent coordinate without modifications 322 do ii=1,natom 323 do jj=1,3 324 if (ga_param%ga_iatfix(jj,ii) == 1) son1(jj,ii)=parent1(jj,ii) 325 enddo 326 enddo 327 call convert_coortogen(son1,coor(:,indiv),natom) 328 acell(:,indiv)=acellson1 329 rprim(:,:,indiv)=rprimson1 330 endif 331 if (checkatomicdist(natom,son2,rprimdson2)==0 .and. indiv<nimage) then 332 indiv=indiv+1 333 !if any fix coordinate restore the parent coordinate without modifications 334 do ii=1,natom 335 do jj=1,3 336 if (ga_param%ga_iatfix(jj,ii) == 1) son2(jj,ii)=parent1(jj,ii) 337 enddo 338 enddo 339 call convert_coortogen(son2,coor(:,indiv),natom) 340 acell(:,indiv)=acellson2 341 rprim(:,:,indiv)=rprimson2 342 endif 343 case(2)! vector flip mutation 344 iparent1=choosefather(fitness,nimage,idum) 345 ii=int(ndimen*uniformrandom(idum)+1) 346 jj=int(ndimen*uniformrandom(idum)+1) 347 if (ii>jj) then 348 call swap(ii,jj) 349 end if 350 vson1(1:ii)=coor_old(1:ii,ihperm(iparent1)) 351 if (jj<ndimen) vson1(jj+1:ndimen)=coor_old(jj+1:ndimen,ihperm(iparent1)) 352 do kk=ii,jj 353 vson1(kk)=coor_old(jj+1-kk,ihperm(iparent1)) 354 enddo 355 rprimson1=rprim_old(:,:,ihperm(iparent1)) 356 acellson1=acell_old(:,ihperm(iparent1)) 357 call convert_gentocoor(son1,vson1,natom) 358 call mkrdim(acellson1,rprimson1,rprimdson1) 359 !if any fix coordinate restore the parent coordinate without modifications 360 do ii=1,natom 361 do jj=1,3 362 if (ga_param%ga_iatfix(jj,ii) == 1) son1(jj,ii)=parent1(jj,ii) 363 enddo 364 enddo 365 call convert_coortogen(son1,vson1,natom) 366 if (checkatomicdist(natom,son1,rprimdson1)==0 .and. indiv<nimage) then 367 indiv=indiv+1 368 coor(:,indiv)=vson1 369 acell(:,indiv)=acellson1 370 rprim(:,:,indiv)=rprimson1 371 endif 372 case(3) ! random strain - nondiagonal 373 iparent1=choosefather(fitness,nimage,idum) 374 rtmp(1,1)=one+gaussian_random(idum,0.1_dp) 375 rtmp(2,2)=one+gaussian_random(idum,0.1_dp) 376 rtmp(3,3)=one+gaussian_random(idum,0.1_dp) 377 rtmp(1,2)=gaussian_random(idum,0.1_dp)*half 378 rtmp(1,3)=gaussian_random(idum,0.1_dp)*half 379 rtmp(2,3)=gaussian_random(idum,0.1_dp)*half 380 rtmp(2,1)=rtmp(1,2) 381 rtmp(3,1)=rtmp(1,3) 382 rtmp(3,2)=rtmp(2,3) 383 rprimson1=matmul(rtmp,rprim_old(:,:,ihperm(iparent1))) 384 vson1=coor_old(:,ihperm(iparent1)) 385 acellson1=acell_old(:,ihperm(iparent1)) 386 call convert_gentocoor(son1,vson1,natom) 387 call mkrdim(acellson1,rprimson1,rprimdson1) 388 !if any fix coordinate restore the parent coordinate without modifications 389 do ii=1,natom 390 do jj=1,3 391 if (ga_param%ga_iatfix(jj,ii) == 1) son1(jj,ii)=parent1(jj,ii) 392 enddo 393 enddo 394 call convert_coortogen(son1,vson1,natom) 395 if (checkatomicdist(natom,son1,rprimdson1)==0 .and. indiv<nimage) then 396 indiv=indiv+1 397 coor(:,indiv)=vson1 398 acell(:,indiv)=acellson1 399 rprim(:,:,indiv)=rprimson1 400 endif 401 case(4) ! coordinates mutation 402 iparent1=choosefather(fitness,nimage,idum) 403 vson1=coor_old(:,ihperm(iparent1)) 404 itmp=ndimen/4 405 if (itmp<1) itmp=1 406 do jj=1,itmp 407 ii=int(ndimen*uniformrandom(idum)+1) 408 vson1(ii)=vson1(ii)+0.15*uniformrandom(idum) 409 if (vson1(ii)>one) vson1(ii)=vson1(ii)-one 410 if (vson1(ii)<zero) vson1(ii)=vson1(ii)+one 411 enddo 412 rprimson1=rprim_old(:,:,ihperm(iparent1)) 413 acellson1=acell_old(:,ihperm(iparent1)) 414 call convert_gentocoor(son1,vson1,natom) 415 call mkrdim(acellson1,rprimson1,rprimdson1) 416 !if any fix coordinate restore the parent coordinate without modifications 417 do ii=1,natom 418 do jj=1,3 419 if (ga_param%ga_iatfix(jj,ii) == 1) son1(jj,ii)=parent1(jj,ii) 420 enddo 421 enddo 422 call convert_coortogen(son1,vson1,natom) 423 if (checkatomicdist(natom,son1,rprimdson1)==0 .and. indiv<nimage) then 424 indiv=indiv+1 425 coor(:,indiv)=vson1 426 acell(:,indiv)=acellson1 427 rprim(:,:,indiv)=rprimson1 428 endif 429 end select 430 enddo 431 432 next_itimimage=itimimage_eff+1 433 if (next_itimimage>ntimimage_stored) next_itimimage=1 434 435 do iimage=1,nimage 436 results_img(iimage,next_itimimage)%acell = acell(:,iimage) 437 results_img(iimage,next_itimimage)%rprim = rprim(:,:,iimage) 438 results_img(iimage,next_itimimage)%vel = results_img(iimage,itimimage_eff)%vel 439 results_img(iimage,next_itimimage)%vel_cell = results_img(iimage,itimimage_eff)%vel_cell 440 call convert_gentocoor(results_img(iimage,next_itimimage)%xred,coor(:,iimage),natom) 441 enddo 442 443 ABI_FREE(coor) 444 ABI_FREE(acell) 445 ABI_FREE(acell_old) 446 ABI_FREE(rprim) 447 ABI_FREE(rprimd) 448 ABI_FREE(rprim_old) 449 ABI_FREE(zcoor1) 450 ABI_FREE(zcoor2) 451 ABI_FREE(coor_old) 452 ABI_FREE(zperm1) 453 ABI_FREE(zperm2) 454 ABI_FREE(ieperm) 455 ABI_FREE(ihperm) 456 ABI_FREE(ibgperm) 457 ABI_FREE(ibgoptperm) 458 ABI_FREE(etotal_img) 459 ABI_FREE(bg_img) 460 ABI_FREE(bg_opt_img) 461 ABI_FREE(enthalpy_img) 462 ABI_FREE(fitness) 463 ABI_FREE(vson1) 464 ABI_FREE(vson2) 465 466 end subroutine predict_ga