TABLE OF CONTENTS


ABINIT/predict_ga [ Modules ]

[ Top ] [ 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