Line data Source code
1 : !!****m* ABINIT/m_predict_string
2 : !! NAME
3 : !! m_predict_string
4 : !!
5 : !! FUNCTION
6 : !!
7 : !!
8 : !! COPYRIGHT
9 : !! Copyright (C) 2009-2026 ABINIT group (XG,ARom,MT)
10 : !! This file is distributed under the terms of the
11 : !! GNU General Public License, see ~abinit/COPYING
12 : !! or http://www.gnu.org/copyleft/gpl.txt .
13 : !!
14 : !! SOURCE
15 :
16 : #if defined HAVE_CONFIG_H
17 : #include "config.h"
18 : #endif
19 :
20 : #include "abi_common.h"
21 :
22 : module m_predict_string
23 :
24 : use defs_basis
25 : use m_abicore
26 : use m_splines
27 : use m_mep
28 : use m_errors
29 : use m_xmpi
30 :
31 : use defs_abitypes, only : MPI_type
32 : use m_results_img, only : results_img_type, gather_array_img, get_geometry_img
33 :
34 : implicit none
35 :
36 : private
37 : !!***
38 :
39 : public :: predict_string
40 : !!***
41 :
42 : contains
43 : !!***
44 :
45 : !!****f* ABINIT/predict_string
46 : !! NAME
47 : !! predict_string
48 : !!
49 : !! FUNCTION
50 : !! Given the past history of images, predict the new set of images using String Method.
51 : !! The changes on the geometry and others are predicted by rescaling the path
52 : !! No change of acell, rprim and vel at present.
53 : !!
54 : !! INPUTS
55 : !! itimimage=time index for image propagation (itimimage+1 is to be predicted here)
56 : !! itimimage_eff=time index in the history
57 : !! list_dynimage(nimage)=list of dynamical images. The non-dynamical ones will not change.
58 : !! Example : in the NEB of string method, one expect the two end images to be fixed.
59 : !! This is quite useful when ground states of the A and B states is known
60 : !! mpi_enreg=MPI-parallelisation information
61 : !! natom=dimension of vel_timimage and xred_timimage
62 : !! ndynimage=number of dynamical images
63 : !! nimage=number of images (on current proc)
64 : !! nimage_tot=total number of images
65 : !! ntimimage_stored=number of time steps stored in the history
66 : !!
67 : !! OUTPUT
68 : !!
69 : !! SIDE EFFECTS
70 : !! mep_param=several parameters for Minimal Energy Path (MEP) search
71 : !! results_img(ntimimage_stored,nimage)=datastructure that holds the history of previous computations.
72 : !! results_img(:,:)%acell(3)
73 : !! at input, history of the values of acell for all images
74 : !! at output, the predicted values of acell for all images
75 : !! results_img(:,:)%results_gs
76 : !! at input, history of the values of energies and forces for all images
77 : !! results_img(:,:)%rprim(3,3)
78 : !! at input, history of the values of rprim for all images
79 : !! at output, the predicted values of rprim for all images
80 : !! results_img(:,:)%vel(3,natom)
81 : !! at input, history of the values of vel for all images
82 : !! at output, the predicted values of vel for all images
83 : !! results_img(:,:)%vel_cell(3,3)
84 : !! at input, history of the values of vel_cell for all images
85 : !! at output, the predicted values of vel_cell for all images
86 : !! results_img(:,:)%xred(3,natom)
87 : !! at input, history of the values of xred for all images
88 : !! at output, the predicted values of xred for all images
89 : !!
90 : !! SOURCE
91 :
92 182 : subroutine predict_string(itimimage,itimimage_eff,list_dynimage,mep_param,mpi_enreg,natom,&
93 182 : & ndynimage,nimage,nimage_tot,ntimimage_stored,results_img)
94 :
95 : !Arguments ------------------------------------
96 : !scalars
97 : integer,intent(in) :: itimimage,itimimage_eff,natom,ndynimage
98 : integer,intent(in) :: nimage,nimage_tot,ntimimage_stored
99 : type(mep_type),intent(inout) :: mep_param
100 : type(MPI_type),intent(in) :: mpi_enreg
101 : !arrays
102 : integer,intent(in) :: list_dynimage(ndynimage)
103 : type(results_img_type) :: results_img(nimage,ntimimage_stored)
104 :
105 : !Local variables-------------------------------
106 : !scalars
107 : integer :: idynimage,ierr,ii,iimage,iatom,next_itimimage
108 : real(dp) :: emax,emin,step
109 : !arrays
110 182 : real(dp),allocatable :: buffer(:,:),buffer_all(:,:)
111 182 : real(dp),allocatable :: darc(:),dimage(:),fcart(:,:,:),rprimd(:,:,:),vect(:,:),wimage(:)
112 182 : real(dp),allocatable :: x(:),y(:),z(:),x2(:),y2(:),z2(:)
113 182 : real(dp),allocatable :: xout(:),yout(:),zout(:)
114 182 : real(dp),allocatable,target :: etotal(:),xcart(:,:,:),xred(:,:,:),strten(:,:)
115 182 : real(dp),pointer :: etotal_all(:),xcart_all(:,:,:),xred_all(:,:,:)
116 :
117 : ! *************************************************************************
118 :
119 728 : ABI_MALLOC(xred,(3,natom,nimage))
120 :
121 : !Parallelism over images: only one process per image of the cell
122 182 : if (mpi_enreg%me_cell==0) then
123 :
124 : ! Retrieve positions and forces
125 258 : ABI_MALLOC(etotal,(nimage))
126 258 : ABI_MALLOC(xcart,(3,natom,nimage))
127 258 : ABI_MALLOC(fcart,(3,natom,nimage))
128 258 : ABI_MALLOC(rprimd,(3,3,nimage))
129 258 : ABI_MALLOC(strten,(6,nimage))
130 : call get_geometry_img(results_img(:,itimimage_eff),etotal,natom,nimage,&
131 86 : & fcart,rprimd,strten,xcart,xred)
132 :
133 : ! EVOLUTION STEP
134 : ! ===============================================
135 :
136 : ! Compute new atomic positions in each cell
137 86 : if (mep_param%mep_solver==MEP_SOLVER_STEEPEST) then ! Steepest-descent
138 86 : call mep_steepest(fcart,list_dynimage,mep_param,natom,natom,ndynimage,nimage,rprimd,xcart,xred)
139 0 : else if (mep_param%mep_solver==MEP_SOLVER_RK4) then ! 4th-order Runge-Kutta
140 0 : call mep_rk4(fcart,itimimage,list_dynimage,mep_param,natom,ndynimage,nimage,rprimd,xcart,xred)
141 : else
142 0 : ABI_BUG("Inconsistent solver !")
143 : end if
144 :
145 : ! REPARAMETRIZATION STEP
146 : ! ===============================================
147 :
148 : ! No reparametrization step in case of Runge-Kutta and mod(istep,4)>0
149 86 : if (mep_param%mep_solver/=MEP_SOLVER_RK4.or.mod(itimimage,4)==0) then
150 :
151 : ! Parallelism: gather data of all images
152 86 : if (mpi_enreg%paral_img==1) then
153 0 : ABI_MALLOC(buffer,(6*natom+1,nimage))
154 0 : ABI_MALLOC(buffer_all,(6*natom+1,nimage_tot))
155 0 : ABI_MALLOC(xred_all,(3,natom,nimage_tot))
156 0 : ABI_MALLOC(xcart_all,(3,natom,nimage_tot))
157 0 : ABI_MALLOC(etotal_all,(nimage_tot))
158 0 : buffer=zero;ii=0
159 0 : buffer(ii+1:ii+3*natom,1:nimage)=reshape(xred ,(/3*natom,nimage/));ii=3*natom
160 0 : buffer(ii+1:ii+3*natom,1:nimage)=reshape(xcart,(/3*natom,nimage/));ii=6*natom
161 0 : buffer(ii+1 ,1:nimage)=etotal(1:nimage);ii=0
162 0 : call gather_array_img(buffer,buffer_all,mpi_enreg,allgather=.true.)
163 0 : xred_all(:,:,:) =reshape(buffer_all(ii+1:ii+3*natom,1:nimage_tot),(/3,natom,nimage_tot/));ii=3*natom
164 0 : xcart_all(:,:,:)=reshape(buffer_all(ii+1:ii+3*natom,1:nimage_tot),(/3,natom,nimage_tot/));ii=6*natom
165 0 : etotal_all(:)=buffer_all(ii+1,1:nimage_tot);ii=0
166 0 : ABI_FREE(buffer)
167 0 : ABI_FREE(buffer_all)
168 : else
169 86 : xred_all => xred
170 86 : xcart_all => xcart
171 86 : etotal_all => etotal
172 : end if
173 :
174 : ! dimage is the distance between two images
175 : ! darc is the parametrization on the string
176 : ! wimage is the weight
177 258 : ABI_MALLOC(darc,(nimage_tot))
178 172 : ABI_MALLOC(dimage,(nimage_tot))
179 172 : ABI_MALLOC(wimage,(nimage_tot))
180 :
181 : ! === Weights for equal arc length
182 86 : if (mep_param%string_algo/=STRING_ALGO_SIMPLIFIED_ENERGY) then
183 674 : wimage(:)=one
184 :
185 : ! === Weights for energy-weight arc length
186 : else
187 0 : emin=min(etotal_all(1),etotal_all(nimage_tot))
188 0 : emax=maxval(etotal_all(:))
189 0 : wimage(1)=one
190 0 : do iimage=2,nimage_tot
191 0 : wimage(iimage)=exp((half*(etotal_all(iimage)+etotal_all(iimage-1))-emin)/(emax-emin))
192 : end do
193 : end if
194 :
195 : ! The distance between images is calculated
196 : ! and normalized to a string length of 1.0
197 674 : dimage=zero
198 258 : ABI_MALLOC(vect,(3,natom))
199 588 : do iimage=2,nimage_tot
200 502 : dimage(iimage)=dimage(iimage-1)
201 : ! MT april 2012: distance must be computed with cartesian coordinates
202 : ! vect(:,:)=xred_all(:,:,iimage)-xred_all(:,:,iimage-1)
203 4518 : vect(:,:)=xcart_all(:,:,iimage)-xcart_all(:,:,iimage-1)
204 588 : dimage(iimage)=dimage(iimage)+wimage(iimage)*mep_img_norm(vect)
205 : end do
206 674 : dimage(:)=dimage(:)/dimage(nimage_tot)
207 86 : ABI_FREE(vect)
208 :
209 : ! Arc lengths
210 86 : darc(1)=zero
211 86 : step=one/dble(nimage_tot-1)
212 588 : do iimage=2,nimage_tot
213 588 : darc(iimage)=darc(iimage-1)+step
214 : end do
215 :
216 : ! New image coordinates are calculated and such that now the mesh is uniform
217 172 : ABI_MALLOC(x,(nimage_tot))
218 172 : ABI_MALLOC(y,(nimage_tot))
219 172 : ABI_MALLOC(z,(nimage_tot))
220 172 : ABI_MALLOC(x2,(nimage_tot))
221 172 : ABI_MALLOC(y2,(nimage_tot))
222 172 : ABI_MALLOC(z2,(nimage_tot))
223 172 : ABI_MALLOC(xout,(nimage_tot))
224 172 : ABI_MALLOC(yout,(nimage_tot))
225 172 : ABI_MALLOC(zout,(nimage_tot))
226 258 : do iatom=1,natom
227 1348 : do iimage=1,nimage_tot
228 1176 : x(iimage)=xred_all(1,iatom,iimage)
229 1176 : y(iimage)=xred_all(2,iatom,iimage)
230 1348 : z(iimage)=xred_all(3,iatom,iimage)
231 : end do
232 172 : call spline(dimage,x,nimage_tot,greatest_real,greatest_real,x2)
233 172 : call spline(dimage,y,nimage_tot,greatest_real,greatest_real,y2)
234 172 : call spline(dimage,z,nimage_tot,greatest_real,greatest_real,z2)
235 172 : call splint(nimage_tot,dimage,x,x2,nimage_tot,darc,xout)
236 172 : call splint(nimage_tot,dimage,y,y2,nimage_tot,darc,yout)
237 172 : call splint(nimage_tot,dimage,z,z2,nimage_tot,darc,zout)
238 : ! After a spline, new image coordinate for that particular
239 : ! atom are generated only if they are dynamical
240 1090 : do idynimage=1,ndynimage
241 832 : iimage=list_dynimage(idynimage)
242 832 : ii=mpi_enreg%my_imgtab(iimage)
243 832 : xred(1,iatom,iimage)=xout(ii)
244 832 : xred(2,iatom,iimage)=yout(ii)
245 1004 : xred(3,iatom,iimage)=zout(ii)
246 : end do
247 : end do ! iatom
248 :
249 : ! Free memory
250 86 : ABI_FREE(x)
251 86 : ABI_FREE(y)
252 86 : ABI_FREE(z)
253 86 : ABI_FREE(x2)
254 86 : ABI_FREE(y2)
255 86 : ABI_FREE(z2)
256 86 : ABI_FREE(xout)
257 86 : ABI_FREE(yout)
258 86 : ABI_FREE(zout)
259 86 : ABI_FREE(darc)
260 86 : ABI_FREE(dimage)
261 86 : ABI_FREE(wimage)
262 86 : if (mpi_enreg%paral_img==1) then
263 0 : ABI_FREE(xred_all)
264 0 : ABI_FREE(xcart_all)
265 0 : ABI_FREE(etotal_all)
266 : end if
267 :
268 : end if ! Reparametrization
269 :
270 : ! ===============================================
271 :
272 86 : ABI_FREE(etotal)
273 86 : ABI_FREE(xcart)
274 86 : ABI_FREE(fcart)
275 86 : ABI_FREE(rprimd)
276 86 : ABI_FREE(strten)
277 : end if ! mpi_enreg%me_cell==0
278 :
279 : !Store acell, rprim, xred and vel for the new iteration
280 182 : call xmpi_bcast(xred,0,mpi_enreg%comm_cell,ierr)
281 182 : next_itimimage=itimimage_eff+1
282 182 : if (next_itimimage>ntimimage_stored) next_itimimage=1
283 1442 : do iimage=1,nimage
284 11340 : results_img(iimage,next_itimimage)%xred(:,:) =xred(:,:,iimage)
285 5040 : results_img(iimage,next_itimimage)%acell(:) =results_img(iimage,itimimage_eff)%acell(:)
286 16380 : results_img(iimage,next_itimimage)%rprim(:,:) =results_img(iimage,itimimage_eff)%rprim(:,:)
287 11340 : results_img(iimage,next_itimimage)%vel(:,:) =results_img(iimage,itimimage_eff)%vel(:,:)
288 16562 : results_img(iimage,next_itimimage)%vel_cell(:,:)=results_img(iimage,itimimage_eff)%vel_cell(:,:)
289 : end do
290 182 : ABI_FREE(xred)
291 :
292 364 : end subroutine predict_string
293 : !!***
294 :
295 : end module m_predict_string
296 : !!***
|