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-2018 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

PARENTS

      predictimg

CHILDREN

      convert_coortogen,convert_gentocoor,initialize_perm,metric,mkradim
      mkrdim,randomize_parent,sort_dp,swap,symanal,symfind,symlatt

SOURCE

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