Branch data Line data Source code
1 : : !!****m* ABINIT/m_gstateimg
2 : : !! NAME
3 : : !! m_gstateimg
4 : : !!
5 : : !! FUNCTION
6 : : !!
7 : : !! COPYRIGHT
8 : : !! Copyright (C) 1998-2026 ABINIT group (XG, AR, GG, MT)
9 : : !! This file is distributed under the terms of the
10 : : !! GNU General Public License, see ~abinit/COPYING
11 : : !! or http://www.gnu.org/copyleft/gpl.txt .
12 : : !!
13 : : !! SOURCE
14 : :
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
66 : : !!***
67 : :
68 : : public :: gstateimg
69 : : !!***
70 : :
71 : : contains
72 : : !!***
73 : :
74 : : !!****f* ABINIT/gstateimg
75 : : !! NAME
76 : : !! gstateimg
77 : : !!
78 : : !! FUNCTION
79 : : !! Routine for conducting DFT calculations for a set of (dynamical) images
80 : : !!
81 : : !! INPUTS
82 : : !! codvsn=code version
83 : : !! cpui=initial CPU time
84 : : !! nimage=number of images of the cell (treated by current proc)
85 : : !! === Optional arguments (needed when nimage>1) ===
86 : : !! filnam(5)=character strings giving file names
87 : : !! filstat=character strings giving name of status file
88 : : !! idtset=index of the dataset
89 : : !! jdtset(0:ndtset)=actual index of the datasets
90 : : !! ndtset=number of datasets
91 : : !!
92 : : !! OUTPUT
93 : : !! etotal_img=total energy, for each image
94 : : !! fcart_img(3,natom,nimage)=forces, in cartesian coordinates, for each image
95 : : !! gred_img(3,natom,nimage)=gradient of E wrt nuclear positions, in reduced coordinates, for each image
96 : : !! intgres_img(nspden,natom,nimage)=gradient wrt constraints, for each image
97 : : !! npwtot(nkpt) = total number of plane waves at each k point
98 : : !! strten_img(6,nimage)=stress tensor, for each image
99 : : !!
100 : : !! SIDE EFFECTS
101 : : !! acell_img(3,nimage)=unit cell length scales (bohr), for each image
102 : : !! amu_img(ntypat,nimage)=value of mass for each atomic type, for each image
103 : : !! dtfil <type(datafiles_type)>=variables related to files
104 : : !! dtset <type(dataset_type)>=all input variables in this dataset
105 : : !! | mband =maximum number of bands (IN)
106 : : !! | mgfft =maximum single fft dimension (IN)
107 : : !! | mkmem =maximum number of k points which can fit in core memory (IN)
108 : : !! | mpw =maximum number of planewaves in basis sphere (large number) (IN)
109 : : !! | natom =number of atoms in unit cell (IN)
110 : : !! | nfft =(effective) number of FFT grid points (for this processor) (IN)
111 : : !! | nkpt =number of k points (IN)
112 : : !! | nspden=number of spin-density components (IN)
113 : : !! | nsppol=number of channels for spin-polarization (1 or 2) (IN)
114 : : !! | nsym =number of symmetry elements in space group
115 : : !! iexit= exit flag
116 : : !! mixalch_img(npspalch,ntypalch,nimage)=value of alchemical mixing factors,for each image
117 : : !! mpi_enreg=MPI-parallelisation information (some already initialized,
118 : : !! some others to be initialized here)
119 : : !! occ_img(mband*nkpt*nsppol,nimage) = occupation number for each band and k, for each image
120 : : !! pawang <type(pawang_type)>=paw angular mesh and related data
121 : : !! pawrad(ntypat*usepaw) <type(pawrad_type)>=paw radial mesh and related data
122 : : !! pawtab(ntypat*usepaw) <type(pawtab_type)>=paw tabulated starting data
123 : : !! psps <type(pseudopotential_type)>=variables related to pseudopotentials
124 : : !! Before entering the first time in gstateimg, a significant part of
125 : : !! psps has been initialized :
126 : : !! the integers dimekb,lmnmax,lnmax,mpssang,mpssoang,mpsso,mgrid,
127 : : !! ntypat,n1xccc,usepaw,useylm, and the arrays dimensioned to npsp
128 : : !! All the remaining components of psps are to be initialized in the call
129 : : !! to pspini .
130 : : !! The next time the code enters gstateimg, psps might be identical to the
131 : : !! one of the previous dtset, in which case, no reinitialisation is scheduled
132 : : !! in pspini.f .
133 : : !! rprim_img(3,3,nimage)=dimensionless real space primitive translations, for each image
134 : : !! vel_cell_img(3,3,nimage)=value of cell parameters velocities, for each image
135 : : !! vel_img(3,natom,nimage)=value of atomic velocities,for each image
136 : : !! xred_img(3,natom,nimage) = reduced atomic coordinates, for each image
137 : : !!
138 : : !! NOTES
139 : : !! USE OF FFT GRIDS:
140 : : !! =================
141 : : !! In case of PAW:
142 : : !! ---------------
143 : : !! Two FFT grids are used:
144 : : !! - A "coarse" FFT grid (defined by ecut)
145 : : !! for the application of the Hamiltonian on the plane waves basis.
146 : : !! It is defined by nfft, ngfft, mgfft, ...
147 : : !! Hamiltonian, wave-functions, density related to WFs (rhor here), ...
148 : : !! are expressed on this grid.
149 : : !! - A "fine" FFT grid (defined) by ecutdg)
150 : : !! for the computation of the density inside PAW spheres.
151 : : !! It is defined by nfftf, ngfftf, mgfftf, ...
152 : : !! Total density, potentials, ...
153 : : !! are expressed on this grid.
154 : : !! In case of norm-conserving:
155 : : !! ---------------------------
156 : : !! - Only the usual FFT grid (defined by ecut) is used.
157 : : !! It is defined by nfft, ngfft, mgfft, ...
158 : : !! For compatibility reasons, (nfftf,ngfftf,mgfftf)
159 : : !! are set equal to (nfft,ngfft,mgfft) in that case.
160 : : !! In case of wavelets:
161 : : !! --------------------
162 : : !! - Only the usual FFT grid (defined by wvl_crmult) is used.
163 : : !! It is defined by nfft, ngfft, mgfft, ... This is strictly not
164 : : !! an FFT grid since its dimensions are not suited for FFTs. They are
165 : : !! defined by wvl_setngfft().
166 : : !! For compatibility reasons, (nfftf,ngfftf,mgfftf)
167 : : !! are set equal to (nfft,ngfft,mgfft) in that case.
168 : : !!
169 : : !! TODO
170 : : !! Not yet possible to use restartxf in parallel when localrdwf==0
171 : : !!
172 : : !! SOURCE
173 : :
174 : 8020 : subroutine gstateimg(acell_img,amu_img,codvsn,cpui,dtfil,dtset,etotal_img,fcart_img,&
175 : 4010 : & gred_img,iexit,intgres_img,mixalch_img,mpi_enreg,nimage,npwtot,occ_img,&
176 : 4010 : & pawang,pawrad,pawtab,psps,&
177 : 4010 : & rprim_img,strten_img,vel_cell_img,vel_img,wvl,xred_img,&
178 [ + - + - : 4010 : & 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 : 4010 : type(mep_type) :: mep_param
229 : 4010 : type(ga_type) :: ga_param
230 : 4010 : type(m1geo_type) :: m1geo_param
231 : 4010 : type(pimd_type) :: pimd_param
232 : : !arrays
233 : 4010 : 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 : 4010 : real(dp),allocatable :: amass(:,:),occ(:),vel(:,:),xred(:,:)
265 : : !real(dp),pointer :: vel_prev(:,:)
266 : 4010 : type(abihist),allocatable :: hist(:),hist_prev(:)
267 : 4010 : type(results_img_type),pointer :: results_img(:,:),res_img(:)
268 : 4010 : type(scf_history_type),allocatable :: scf_history(:)
269 : : !type(abiforstr) :: preconforstr ! Preconditioned forces and stress ... Only needed to deallocate an internal matrix in prec_simple
270 : :
271 : : ! ***********************************************************************
272 : :
273 : : DBG_ENTER("COLL")
274 : :
275 : 4010 : call timab(1200,1,tsec)
276 : 4010 : call timab(1203,3,tsec)
277 : :
278 : : !Arguments check
279 [ + + ]: 4010 : if (dtset%nimage>1) then
280 : : if ((.not.present(filnam)).or.(.not.present(filnam)).or.(.not.present(idtset)).or.&
281 [ + - - + ]: 50 : & (.not.present(ndtset)).or.(.not.present(jdtset))) then
282 : : write(msg,'(3a)') &
283 : 0 : & 'When nimage>1, all the following argument should be present:',ch10,&
284 : 0 : & 'filnam, filstat, idtset, ndtset, jdtset !'
285 : 0 : ABI_BUG(msg)
286 : : end if
287 : : end if
288 : :
289 : : !Set flag for the effective computation (once) of static images
290 : : !For the time being only set when parallelization is activated
291 : : !Note: if you modify this flag, do not forget to change it in outvars and outvar1
292 : 4010 : compute_static_images=(dtset%istatimg>0)
293 : :
294 : : !Prepare the calculation, by computing flags and dimensions
295 : 4010 : is_pimd=(dtset%imgmov==9.or.dtset%imgmov==10.or.dtset%imgmov==13)
296 : 4010 : is_mep =(dtset%imgmov==1.or.dtset%imgmov== 2.or.dtset%imgmov== 5)
297 : 4010 : ntimimage=dtset%ntimimage
298 [ + + ]: 4010 : ntimimage_stored=ntimimage;if(is_pimd)ntimimage_stored=2
299 : 4010 : nocc=dtset%mband*dtset%nkpt*dtset%nsppol
300 [ + + - + ]: 4010 : is_master=(mpi_enreg%me_cell==0.and.mpi_enreg%me_img==0)
301 : 4010 : delta_energy=zero
302 : :
303 : : !Management of dynamics/relaxation history (positions, forces, stresses, ...)
304 [ + + - + ]: 4010 : use_hist=(dtset%imgmov/=0.and.nimage>0) ; use_hist_prev=.false.
305 : 4010 : isVused=is_pimd;isARused=(dtset%optcell/=0)
306 [ + + ]: 4010 : if (use_hist) then
307 : : !Read history from file (and broadcast if MPI)
308 [ + + - + ]: 43 : use_hist_prev=(dtset%restartxf==-1.and.nimage>0)
309 [ + - + - ]: 43 : hist_filename=trim(dtfil%filnam_ds(4))//'_HIST.nc'
310 [ + + ]: 43 : if (use_hist_prev)then
311 [ - + + - : 8 : ABI_MALLOC(hist_prev,(nimage))
+ - + - -
+ - + - +
+ + + - -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ ]
312 [ + - ]: 1 : if (mpi_enreg%me_cell==0) then
313 : : call read_md_hist_img(hist_filename,hist_prev,isVused,isARused,&
314 : 1 : & imgtab=mpi_enreg%my_imgtab)
315 : : end if
316 : 1 : call abihist_bcast(hist_prev,0,mpi_enreg%comm_cell)
317 [ - + ]: 1 : if (nimage>0) then
318 [ + + + - : 6 : if (any(hist_prev(:)%mxhist/=hist_prev(1)%mxhist)) then
- + ]
319 : 0 : msg='History problem: all images should have the same number of time steps!'
320 : 0 : ABI_ERROR(msg)
321 : : end if
322 : 1 : use_hist_prev=(hist_prev(1)%mxhist>0)
323 [ + - ]: 1 : if (use_hist_prev) ntimimage=ntimimage+hist_prev(1)%mxhist
324 : : end if
325 [ - + ]: 1 : if (.not.use_hist_prev) then
326 : 0 : call abihist_free(hist_prev)
327 [ # # # # : 0 : ABI_FREE(hist_prev)
# # # # #
# # # # #
# # # # #
# # # # #
# # # # ]
328 : : end if
329 : : end if
330 : : !Initialize a variable to write the history
331 [ - + + - : 395 : ABI_MALLOC(hist,(nimage))
+ - + - -
+ - + - +
+ + + - -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ ]
332 : 43 : call abihist_init(hist,dtset%natom,ntimimage,isVused,isARused)
333 : : end if ! imgmov/=0
334 : :
335 : : !Various allocations
336 [ - + + - : 12030 : ABI_MALLOC(occ,(nocc))
+ - + - -
+ - + ]
337 [ - + + - : 12030 : ABI_MALLOC(vel,(3,dtset%natom))
+ - + - -
+ - + ]
338 [ - + + - : 8020 : ABI_MALLOC(xred,(3,dtset%natom))
- + - + ]
339 : :
340 : : !Select dynamical images
341 [ - + + - : 12030 : ABI_MALLOC(list_dynimage,(dtset%ndynimage))
+ - + - -
+ - + ]
342 : 4010 : ndynimage=0
343 [ + + ]: 8280 : do iimage=1,nimage
344 : 4270 : ii=mpi_enreg%my_imgtab(iimage)
345 [ + + ]: 8280 : if (dtset%dynimage(ii)==1) then
346 : 4210 : ndynimage=ndynimage+1
347 : 4210 : list_dynimage(ndynimage)=iimage
348 : : end if
349 : : end do
350 : :
351 : : !Management of SCF history (density/WF predictions from one time step to another)
352 [ - + + - : 16300 : ABI_MALLOC(scf_history,(nimage))
+ - + - -
+ - + + +
+ - - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - - - -
- - - - -
- - - - -
- + - - -
- - - - -
- - - - -
- - + - -
- - - - ]
353 [ - + + - : 12030 : ABI_MALLOC(scf_initialized,(nimage))
+ - + - -
+ - + ]
354 [ + + ]: 8280 : scf_initialized=0
355 : 4010 : history_size=-1
356 [ + + ]: 4010 : if (dtset%ntimimage<=1) then
357 [ + - + + : 3967 : if (dtset%usewvl==0.and.dtset%ionmov>0.and. (abs(dtset%densfor_pred)==5.or.abs(dtset%densfor_pred)==6)) then
+ + ]
358 : 98 : history_size=2
359 [ + - ]: 98 : if(dtset%extrapwf==2) history_size=3
360 : : end if
361 : : else
362 [ + + ]: 43 : if (abs(dtset%densfor_pred)==2.or.abs(dtset%densfor_pred)==3) history_size=0
363 [ + + ]: 43 : if (dtset%imgwfstor==1) history_size=1
364 [ + - - + ]: 141 : if (dtset%usewvl==0.and.(abs(dtset%densfor_pred)==5.or.abs(dtset%densfor_pred)==6)) history_size=2
365 : : end if
366 [ + + ]: 8280 : do iimage=1,nimage
367 : 4270 : call scf_history_nullify(scf_history(iimage))
368 : 8280 : scf_history(iimage)%history_size=history_size
369 : : end do
370 : :
371 : : !In some cases, need amass variable
372 [ + + ]: 4010 : if (use_hist) then
373 [ - + + - : 172 : ABI_MALLOC(amass,(dtset%natom,nimage))
- + + - +
- + - - +
- + ]
374 [ + + ]: 309 : do iimage=1,nimage
375 [ + + + - : 737 : if (any(amu_img(:,iimage)/=amu_img(:,1))) then
- + ]
376 : 0 : ABI_ERROR('Problem with variable masses!')
377 : : end if
378 [ + - + + : 1713 : amass(:,iimage)=amu_emass*amu_img(dtset%typat(:),iimage)
+ + ]
379 : : end do
380 : : end if
381 : :
382 : : !In the case of the 4th-order Runge-Kutta solver,
383 : : !one must have a number of step multiple of 4.
384 : 4010 : ntimimage_max=ntimimage;idelta=1
385 [ + + - + ]: 4010 : if (dtset%imgmov==2.and.mep_param%mep_solver==4) then
386 : 0 : ntimimage_max=4*(ntimimage_max/4)
387 : 0 : idelta=4
388 : : end if
389 : :
390 : : !MEP search: fill in eventually the data structure mep_param
391 : 4010 : call mep_init(dtset,mep_param)
392 : :
393 : : !GA search: fill in eventually the data structure ga_param
394 : 4010 : call ga_init(dtset,ga_param)
395 : :
396 : : !Move 1GEO approach: fill the data structure m1geo_param
397 : 4010 : call m1geo_init(dtfil,dtset,m1geo_param)
398 : :
399 : : !PIMD: fill in the data structure pimd_param and init velocities
400 : 4010 : call pimd_init(dtset,pimd_param,is_master)
401 [ + + ]: 4010 : dtion=one;if (is_pimd) dtion=pimd_param%dtion
402 : :
403 : : !Initialization of main variables
404 [ - + + - : 26735 : ABI_MALLOC(results_img,(nimage,ntimimage_stored))
- + + - +
- + - - +
- + + + +
+ + - - +
- + - + -
+ - + - +
- + ]
405 [ + + ]: 8361 : do itimimage=1,ntimimage_stored
406 : 4351 : res_img => results_img(:,itimimage)
407 : : call init_results_img(dtset%natom,dtset%npspalch,dtset%nspden,dtset%nsppol,dtset%ntypalch,&
408 : 4351 : & dtset%ntypat,res_img)
409 [ + + ]: 14705 : do iimage=1,nimage
410 [ + + ]: 25376 : res_img(iimage)%acell(:) =acell_img(:,iimage)
411 [ + + ]: 15793 : res_img(iimage)%amu(:) =amu_img(:,iimage)
412 [ + + + + ]: 6419 : res_img(iimage)%mixalch(:,:) =mixalch_img(:,:,iimage)
413 [ + + + + ]: 82472 : res_img(iimage)%rprim(:,:) =rprim_img(:,:,iimage)
414 [ + + + + ]: 64736 : res_img(iimage)%xred(:,:) =xred_img(:,:,iimage)
415 [ + + + + ]: 64736 : res_img(iimage)%vel(:,:) =vel_img(:,:,iimage)
416 [ + + + + ]: 86823 : res_img(iimage)%vel_cell(:,:)=vel_cell_img(:,:,iimage)
417 : : end do
418 : : end do
419 : :
420 : : !Set Number of degrees Of Freedom for PIMD algorithms
421 [ + + ]: 4010 : if(use_hist) then
422 [ + + ]: 309 : hist(:)%ndof=3*dtset%natom ! Init ndof to 3N (no iatfix for PIMD)
423 [ + + ]: 43 : if(dtset%imgmov==9.or.dtset%imgmov==10) then
424 : : if(pimd_param%pitransform==1.or.pimd_param%pitransform==2.or.&
425 [ + + + + : 4 : & pimd_param%constraint==1.or.pimd_param%optcell==2) then
- + ]
426 [ + + ]: 21 : hist(:)%ndof=hist(:)%ndof-3
427 : : end if
428 : : end if
429 : : end if
430 : :
431 : 4010 : call timab(1203,2,tsec)
432 : :
433 : : !-----------------------------------------------------------------------------------------
434 : : !Big loop on the propagation of all images
435 : 4010 : itimimage_eff=1 ; itimimage_gstate=1
436 [ + - ]: 4306 : do itimimage=1,ntimimage
437 : :
438 : 4306 : res_img => results_img(:,itimimage_eff)
439 : 4306 : call_predictor=(ntimimage>1)
440 : :
441 : : ! If history is activated and if current image is inside it: do not compute anything
442 [ + + ]: 4306 : if (use_hist_prev) then
443 [ + + + + : 35 : if (all(hist_prev(:)%ihist<=hist_prev(:)%mxhist)) then
+ + ]
444 [ + + ]: 30 : do iimage=1,nimage
445 : 25 : ih=hist_prev(iimage)%ihist
446 : 25 : call abihist_copy(hist_prev(iimage),hist(iimage))
447 : 25 : call mkradim(hist_prev(iimage)%acell(:,ih),rprim,hist_prev(iimage)%rprimd(:,:,ih))
448 [ + + ]: 100 : res_img(iimage)%acell(:)=hist_prev(iimage)%acell(:,ih)
449 [ + + + + ]: 325 : res_img(iimage)%rprim(:,:)=rprim
450 [ + + + + ]: 625 : res_img(iimage)%xred(:,:)=hist_prev(iimage)%xred(:,:,ih)
451 [ + + + + ]: 625 : res_img(iimage)%vel(:,:)=hist_prev(iimage)%vel(:,:,ih)
452 [ + + + + ]: 325 : res_img(iimage)%vel_cell(:,:)=hist_prev(iimage)%vel_cell(:,:,ih)
453 [ + + + + ]: 625 : res_img(iimage)%results_gs%fcart(:,:)=hist_prev(iimage)%fcart(:,:,ih)
454 [ + + ]: 175 : res_img(iimage)%results_gs%strten(:)=hist_prev(iimage)%strten(:,ih)
455 : 25 : res_img(iimage)%results_gs%etotal=hist_prev(iimage)%etot(ih)
456 : 25 : res_img(iimage)%results_gs%energies%entropy=hist_prev(iimage)%entropy(ih)
457 : : call fcart2gred(res_img(iimage)%results_gs%fcart,res_img(iimage)%results_gs%gred,&
458 : 25 : hist_prev(iimage)%rprimd(:,:,ih),dtset%natom)
459 : 30 : hist_prev(iimage)%ihist=hist_prev(iimage)%ihist+1
460 : : end do
461 : : !PI-QTB: skip a record in random force file
462 [ - + ]: 5 : if (pimd_param%use_qtb==1) call pimd_skip_qtb(pimd_param)
463 : : !call_predictor=.false.
464 : : goto 110 ! This is temporary
465 : : end if
466 : : end if
467 : :
468 : 4301 : call timab(1204,1,tsec)
469 : 4301 : call localfilnam(mpi_enreg%comm_img,mpi_enreg%comm_cell,mpi_enreg%comm_world,filnam,'_IMG',dtset%nimage)
470 : 4301 : compute_all_images=(compute_static_images.and.itimimage==1)
471 : :
472 : : ! Print title for time step
473 [ + + + + ]: 4301 : if (dtset%nimage>1.or.dtset%ntimimage>1) then
474 [ + + ]: 342 : if (dtset%prtvolimg<2) then
475 [ + + ]: 274 : msg=ch10;if (itimimage >1) write(msg,'(2a)') ch10,ch10
476 : 274 : write(msg,'(5a)') trim(msg),&
477 : 274 : '================================================================================',&
478 [ + - ]: 548 : ch10,' ',trim(imagealgo_str(dtset%imgmov))
479 : : else
480 [ + + ]: 68 : msg='';if (itimimage >1) msg=ch10
481 : 68 : write(msg,'(5a)') trim(msg),&
482 : 68 : '--------------------------------------------------------------------------------',&
483 [ + - ]: 136 : ch10,' ',trim(imagealgo_str(dtset%imgmov))
484 : : end if
485 [ + + ]: 342 : if (dtset%imgmov==2) then
486 [ + - ]: 206 : write(msg,'(6a)') trim(msg),' (',trim(stgalgo_str(mep_param%string_algo)),' + ',&
487 [ + - ]: 412 : trim(mepsolver_str(mep_param%mep_solver)),')'
488 : : end if
489 [ + + ]: 342 : if (dtset%imgmov==5) then
490 [ + + + + ]: 46 : ii=merge(mep_param%neb_algo,1,mep_param%neb_algo/=2.or.itimimage>=mep_param%cineb_start)
491 [ + - ]: 46 : write(msg,'(6a)') trim(msg),' (',trim(nebalgo_str(ii)),' + ',&
492 [ + - ]: 92 : trim(mepsolver_str(mep_param%mep_solver)),')'
493 : : end if
494 [ + + ]: 342 : if (dtset%ntimimage==1) write(msg,'(2a)') trim(msg),' FOR 1 TIME STEP'
495 [ + + ]: 342 : if (dtset%ntimimage >1) write(msg,'(2a,i5)') trim(msg),' - TIME STEP ',itimimage
496 [ + + ]: 342 : if (dtset%prtvolimg<2) then
497 : 274 : write(msg,'(3a)') trim(msg),ch10,&
498 : 548 : & '================================================================================'
499 : : end if
500 : 342 : call wrtout(ab_out ,msg,'COLL')
501 : 342 : call wrtout(std_out,msg,'PERS')
502 : :
503 : 342 : call yaml_iterstart('timimage', itimimage, dev_null, dtset%use_yaml)
504 : : end if
505 : :
506 [ + + ]: 4301 : if (dtset%use_yaml == 1) call yaml_iterstart('timimage', itimimage, ab_out, dtset%use_yaml)
507 : :
508 : 4301 : call timab(1204,2,tsec)
509 : :
510 : : ! Loop on the dynamical images
511 : 4301 : idynimage=0
512 [ + + ]: 10417 : do iimage=1,nimage
513 : :
514 : 6116 : ii=mpi_enreg%my_imgtab(iimage)
515 : 6116 : if (dtset%dynimage(ii)==1) idynimage=idynimage+1
516 : :
517 : : ! Compute static image only at first time step
518 [ + + + + ]: 6116 : if (dtset%dynimage(ii)==1.or.compute_all_images) then
519 : :
520 : 5664 : call timab(1205,1,tsec)
521 : :
522 : : ! Change file names according to image index (if nimage>1)
523 [ + + ]: 5664 : if (dtset%nimage>1) then
524 [ - + - - : 1700 : call dtfil%init(dtset,filnam,filstat,idtset,jdtset,mpi_enreg,ndtset,image_index=ii)
- - - + -
- - - - +
- + - - ]
525 [ + + ]: 1700 : if (itimimage>1) then
526 : 1399 : dtfil%ireadwf=0;dtfil%ireadden=0;dtfil%ireadkden=0
527 : : end if
528 : 1700 : call yaml_iterstart('image', iimage, dev_null, 0)
529 : : end if
530 : :
531 [ + + ]: 5664 : if (dtset%use_yaml == 1) call yaml_iterstart('image', iimage, ab_out, dtset%use_yaml)
532 : :
533 : : ! Redefine output units
534 : 5664 : call localwrfile(mpi_enreg%comm_cell,ii,dtset%nimage,mpi_enreg%paral_img,dtset%prtvolimg)
535 : :
536 : : ! Print title for image
537 [ + + + + : 5664 : if (dtset%nimage>1.and.(dtset%prtvolimg==0.or.do_write_log)) then
+ + ]
538 [ + + ]: 1324 : if (ii==1) write(msg,'(a)' ) ch10
539 [ + + ]: 1324 : if (ii >1) write(msg,'(2a)') ch10,ch10
540 : 1324 : write(msg,'(6a,i4,a,i4,3a)') trim(msg),&
541 : 1324 : '--------------------------------------------------------------------------------',ch10,&
542 [ + - ]: 1324 : ' ',trim(imagealgo_str(dtset%imgmov)),' - CELL # ',ii,'/',dtset%nimage,ch10,&
543 : 2648 : '--------------------------------------------------------------------------------',ch10
544 [ + + ]: 1324 : if (dtset%prtvolimg==0) call wrtout(ab_out ,msg,'COLL')
545 [ + + ]: 1324 : if (do_write_log) call wrtout(std_out,msg,'PERS')
546 : : end if
547 : :
548 [ + + ]: 22656 : acell(:) =res_img(iimage)%acell(:)
549 [ + + + + ]: 73632 : rprim(:,:) =res_img(iimage)%rprim(:,:)
550 [ + + + + ]: 58856 : vel(:,:) =res_img(iimage)%vel(:,:)
551 [ + + + + ]: 73632 : vel_cell(:,:)=res_img(iimage)%vel_cell(:,:)
552 [ + + + + ]: 58856 : xred(:,:) =res_img(iimage)%xred(:,:)
553 [ + + ]: 453298 : occ(:) =occ_img(:,iimage)
554 : :
555 : : call args_gs_init(args_gs, &
556 : : res_img(iimage)%amu(:),dtset%cellcharge(ii),res_img(iimage)%mixalch(:,:),&
557 : 5664 : dtset%dmatpawu(:,:,:,:,ii),dtset%upawu(:,ii),dtset%jpawu(:,ii),dtset%rprimd_orig(:,:,ii))
558 : :
559 : 5664 : call timab(1205,2,tsec)
560 : :
561 : : call gstate(args_gs,acell,codvsn,cpui,dtfil,dtset,iexit,scf_initialized(iimage),itimimage_gstate,&
562 : : mpi_enreg,npwtot,occ,pawang,pawrad,pawtab,psps,&
563 : : res_img(iimage)%results_gs,&
564 : 5664 : rprim,scf_history(iimage),vel,vel_cell,wvl,xred)
565 : :
566 : 5664 : itimimage_gstate=itimimage_gstate+1
567 : :
568 : 5664 : call timab(1206,1,tsec)
569 : :
570 : 5664 : call args_gs_free(args_gs)
571 : :
572 [ + + ]: 5664 : if (dtset%dynimage(ii)==1) then
573 [ + + ]: 22432 : res_img(iimage)%acell(:) =acell(:)
574 [ + + + + ]: 72904 : res_img(iimage)%rprim(:,:) =rprim(:,:)
575 [ + + + + ]: 58304 : res_img(iimage)%vel(:,:) =vel(:,:)
576 [ + + + + ]: 72904 : res_img(iimage)%vel_cell(:,:)=vel_cell(:,:)
577 [ + + + + ]: 58304 : res_img(iimage)%xred(:,:) =xred(:,:)
578 [ + + ]: 452126 : occ_img(:,iimage) =occ(:)
579 : : end if
580 : :
581 : : ! check change of rprim and rewriting in hist
582 : : ! check change of xred and rewriting in hist
583 : :
584 : : ! Close output units ; restore defaults
585 : 5664 : call localredirect(mpi_enreg%comm_cell,mpi_enreg%comm_world,dtset%nimage,mpi_enreg%paral_img,dtset%prtvolimg)
586 : 5664 : call timab(1206,2,tsec)
587 : :
588 [ + + ]: 452 : else if (itimimage>1) then ! For static images, simply copy one time step to the other
589 [ - + ]: 448 : itimimage_prev=itimimage_eff-1;if (itimimage_prev<1) itimimage_prev=ntimimage_stored
590 : : call copy_results_img(results_img(iimage,itimimage_prev), &
591 : 448 : & results_img(iimage,itimimage_eff ))
592 : : end if
593 : :
594 : : ! Store results in hist datastructure
595 [ + + ]: 10417 : if (use_hist) then
596 : 2110 : ih=hist(iimage)%ihist
597 : 2110 : call mkrdim(res_img(iimage)%acell(:),res_img(iimage)%rprim(:,:),rprimd)
598 : : call var2hist(res_img(iimage)%acell(:),hist(iimage),dtset%natom,&
599 : 2110 : & rprimd,res_img(iimage)%xred(:,:),.FALSE.)
600 : 2110 : call vel2hist(amass(:,iimage),hist(iimage),res_img(iimage)%vel(:,:),res_img(iimage)%vel_cell(:,:))
601 [ + + + + ]: 21094 : hist(iimage)%fcart(:,:,ih)=res_img(iimage)%results_gs%fcart(:,:)
602 [ + + ]: 14770 : hist(iimage)%strten(:,ih)=res_img(iimage)%results_gs%strten(:)
603 : 2110 : hist(iimage)%etot(ih)=res_img(iimage)%results_gs%etotal
604 : 2110 : hist(iimage)%entropy(ih)=res_img(iimage)%results_gs%energies%entropy
605 : 2110 : hist(iimage)%time(ih)=real(itimimage,kind=dp)*dtion
606 : : end if
607 : :
608 : : end do ! iimage
609 : :
610 [ - + ]: 4301 : if(mpi_enreg%paral_img==1)then
611 : 0 : call timab(1208,1,tsec)
612 : 0 : call xmpi_barrier(mpi_enreg%comm_img)
613 : 0 : call timab(1208,2,tsec)
614 : : end if
615 : :
616 : 4301 : call timab(1209,1,tsec)
617 : :
618 : : ! Output when images are used
619 [ + + ]: 4301 : if (dtset%nimage>1) then
620 : : ! === 1st option: reduced outputs ===
621 [ + + ]: 337 : if (dtset%prtvolimg>0) then
622 : : call prtimg(dtset%dynimage,imagealgo_str(dtset%imgmov),dtset%imgmov,ab_out,&
623 [ - + ]: 131 : & mpi_enreg,nimage,dtset%nimage,compute_all_images,dtset%prtvolimg,res_img)
624 : : end if
625 : : end if
626 : :
627 : : ! Manage log files when images are used
628 : : call localrdfile(mpi_enreg%comm_img,mpi_enreg%comm_world,compute_all_images,&
629 : 4301 : & dtset%nimage,mpi_enreg%paral_img,dtset%prtvolimg,dyn=dtset%dynimage)
630 : :
631 : : ! Write hist datastructure in HIST file
632 : : ! Note: for PIMD, writing is done later
633 : : ! if (use_hist.and.mpi_enreg%me_cell==0) then
634 [ + + + + : 4301 : if (use_hist.and.mpi_enreg%me_cell==0.and.(.not.is_pimd)) then
+ + ]
635 [ + + ]: 190 : ifirst=merge(0,1,itimimage>1)
636 : : call write_md_hist_img(hist,hist_filename,ifirst,itimimage,dtset%natom,dtset%ntypat,&
637 : : & dtset%typat,amu_img(:,1),dtset%znucl,dtion,&
638 : : & nimage=dtset%nimage,imgmov=dtset%imgmov,mdtemp=dtset%mdtemp,comm_img=mpi_enreg%comm_img,&
639 : 190 : & imgtab=mpi_enreg%my_imgtab)
640 : : end if
641 : :
642 : : ! TESTS WHETHER ONE CONTINUES THE LOOP
643 : : ! Here we calculate the change in energy, and exit if delta_energy < tolimg
644 : 35 : delta_energy=zero
645 : : ! Doesn't check convergence in case of PIMD
646 [ + - + + ]: 4266 : check_conv=((.not.is_pimd).and.itimimage>1)
647 : : ! In case of 4th-order Runge-Kutta, does check convergence every 4 steps
648 [ + + - + ]: 4301 : if (dtset%imgmov==2.and.mep_param%mep_solver==4) then
649 [ # # # # ]: 0 : check_conv=(mod(itimimage,4)==0.and.itimimage>4)
650 : : end if
651 [ + + ]: 4301 : if (check_conv) then
652 [ + + ]: 1497 : do idynimage=1,ndynimage
653 : 1234 : iimage=list_dynimage(idynimage)
654 : : delta_energy=delta_energy &
655 : : & +abs(results_img(iimage,itimimage)%results_gs%etotal &
656 : 1497 : & -results_img(iimage,itimimage-idelta)%results_gs%etotal)
657 : : end do
658 [ - + ]: 263 : if (mpi_enreg%paral_img==1) then
659 : 0 : call xmpi_sum(delta_energy,mpi_enreg%comm_img,ierr)
660 : : end if
661 : 263 : delta_energy=delta_energy/dtset%ndynimage
662 [ + + ]: 263 : if (delta_energy<dtset%tolimg) then
663 [ + + ]: 33 : if (dtset%prtvolimg<2) then
664 : 26 : write(msg,'(5a,i5,6a,es11.3,a,es11.3,2a)') ch10,ch10,&
665 : 26 : & '================================================================================',ch10,&
666 : 26 : & ' At time step ',itimimage,ch10,&
667 [ + - ]: 26 : & ' ',trim(imagealgo_str(dtset%imgmov)),' has reached energy convergence',ch10,&
668 : 26 : & ' with Average[Abs(Etotal(t)-Etotal(t-dt))]=',delta_energy,'<tolimg=',dtset%tolimg,ch10,&
669 : 52 : & '================================================================================'
670 : : else
671 : 7 : write(msg,'(4a,i5,6a,es11.3,a,es11.3)') ch10,&
672 : 7 : & '--------------------------------------------------------------------------------',ch10,&
673 : 7 : & ' At time step ',itimimage,ch10,&
674 [ + - ]: 7 : & ' ',trim(imagealgo_str(dtset%imgmov)),' has reached energy convergence',ch10,&
675 : 14 : & ' with Average[Abs(Etotal(t)-Etotal(t-dt))]=',delta_energy,'<tolimg=',dtset%tolimg
676 : : end if
677 [ + + ]: 99 : call wrtout([std_out, ab_out] ,msg,'COLL')
678 : 33 : call timab(1209,2,tsec) ! This is the first place where counter 1209 is stopped.
679 : 33 : exit ! exit itimimage
680 : : end if
681 : : end if
682 : :
683 : : !Temporary statement
684 : : ! 110 continue
685 : :
686 : : ! Dont call the predictor at last time step (except for PIMD)
687 [ + + ]: 4268 : if (itimimage>=ntimimage_max) call_predictor=(call_predictor.and.is_pimd)
688 : :
689 : : ! Predict the next value of the images
690 [ + + ]: 4268 : if (call_predictor) then
691 : : call predictimg(delta_energy,imagealgo_str(dtset%imgmov),dtset%imgmov,itimimage,&
692 : : & itimimage_eff,list_dynimage,ga_param,mep_param,mpi_enreg,m1geo_param,dtset%natom,ndynimage,&
693 [ - + ]: 298 : & nimage,dtset%nimage,ntimimage_stored,pimd_param,dtset%prtvolimg,results_img)
694 : : end if
695 : :
696 : : ! Write hist datastructure in HIST file in case of PIMD
697 : : ! Note : velocities have been updated
698 : : #if defined HAVE_NETCDF
699 [ + + ]: 4268 : if (use_hist.and.is_pimd) then
700 [ + + ]: 240 : do iimage=1,nimage
701 : 205 : ih=hist(iimage)%ihist
702 : : call vel2hist(amass(:,iimage),hist(iimage),results_img(iimage,itimimage_eff)%vel(:,:),&
703 : 240 : & results_img(iimage,itimimage_eff)%vel_cell(:,:))
704 : : end do
705 [ - + ]: 35 : if (mpi_enreg%me_cell==0) then
706 [ + + ]: 35 : ifirst=merge(0,1,itimimage>1)
707 : : call write_md_hist_img(hist,hist_filename,ifirst,itimimage,dtset%natom,dtset%ntypat,&
708 : : & dtset%typat,amu_img(:,1),dtset%znucl,dtion,&
709 : : & nimage=dtset%nimage,imgmov=dtset%imgmov,mdtemp=dtset%mdtemp,comm_img=mpi_enreg%comm_img,&
710 : 35 : & imgtab=mpi_enreg%my_imgtab)
711 : : end if
712 : : end if
713 : : #endif
714 : :
715 : : !Temporary statement
716 : : 110 continue
717 : :
718 : : ! Increment indexes
719 [ + + ]: 4273 : if (itimimage>=ntimimage_max) exit
720 [ + + ]: 296 : itimimage_eff=itimimage_eff+1;if (itimimage_eff>ntimimage_stored) itimimage_eff=1
721 [ + + ]: 296 : if (use_hist) then
722 [ + + ]: 2164 : do iimage=1,nimage
723 : 2164 : hist(iimage)%ihist=hist(iimage)%ihist+1
724 : : end do
725 : : end if
726 : :
727 : 4273 : call timab(1209,2,tsec) ! This is the second place where counter 1209 is stopped.
728 : :
729 : : end do ! itimimage
730 : : !-----------------------------------------------------------------------------------------
731 : :
732 : 4010 : call timab(1210,1,tsec)
733 : :
734 : : !Copy the results of the computation in the appropriate arguments of the routine
735 [ + + ]: 8280 : do iimage=1,nimage
736 : 4270 : ii=mpi_enreg%my_imgtab(iimage)
737 [ + + ]: 8280 : if (dtset%dynimage(ii)==1) then
738 [ + + ]: 16840 : acell_img(:,iimage) =results_img(iimage,itimimage_eff)%acell(:)
739 [ + + ]: 9702 : amu_img(:,iimage) =results_img(iimage,itimimage_eff)%amu(:)
740 [ + + + + ]: 4285 : mixalch_img(:,:,iimage) =results_img(iimage,itimimage_eff)%mixalch(:,:)
741 [ + + + + ]: 54730 : rprim_img(:,:,iimage) =results_img(iimage,itimimage_eff)%rprim(:,:)
742 [ + + + + ]: 44298 : vel_img(:,:,iimage) =results_img(iimage,itimimage_eff)%vel(:,:)
743 [ + + + + ]: 54730 : vel_cell_img(:,:,iimage)=results_img(iimage,itimimage_eff)%vel_cell(:,:)
744 [ + + + + ]: 44298 : xred_img(:,:,iimage) =results_img(iimage,itimimage_eff)%xred(:,:)
745 : 4210 : etotal_img(iimage) =results_img(iimage,itimimage_eff)%results_gs%etotal
746 [ + + + + ]: 44298 : fcart_img(:,:,iimage) =results_img(iimage,itimimage_eff)%results_gs%fcart(:,:)
747 [ + + + + ]: 44298 : gred_img(:,:,iimage) =results_img(iimage,itimimage_eff)%results_gs%gred(:,:)
748 [ + + + + ]: 26212 : intgres_img(:,:,iimage) =results_img(iimage,itimimage_eff)%results_gs%intgres(:,:)
749 [ + + ]: 29470 : strten_img(:,iimage) =results_img(iimage,itimimage_eff)%results_gs%strten(:)
750 [ + + ]: 60 : else if (compute_static_images) then
751 : 56 : etotal_img(iimage) =results_img(iimage,1)%results_gs%etotal
752 [ + + + + ]: 552 : fcart_img(:,:,iimage) =results_img(iimage,1)%results_gs%fcart(:,:)
753 [ + + + + ]: 552 : gred_img(:,:,iimage) =results_img(iimage,1)%results_gs%gred(:,:)
754 [ + + + + ]: 304 : intgres_img(:,:,iimage)=results_img(iimage,1)%results_gs%intgres(:,:)
755 [ + + ]: 392 : strten_img(:,iimage) =results_img(iimage,1)%results_gs%strten(:)
756 : : end if
757 : : end do
758 : :
759 : : !When parallelizattion over images is activated, has to sum number of warnings
760 : : !and comments written in log file
761 [ - + - - ]: 4010 : if (mpi_enreg%paral_img==1.and.mpi_enreg%comm_cell==0) then
762 : 0 : call specialmsg_mpisum(mpi_enreg%comm_img)
763 : 0 : call libpaw_spmsg_mpisum(mpi_enreg%comm_img)
764 : : end if
765 : :
766 : :
767 : : !Final deallocations
768 : :
769 : : !This call is needed to free internal storages in different routines (prec_simple, pred_bfgs ...)
770 [ + + ]: 4010 : if(dtset%imgmov==6)then
771 : 6 : m1geo_param%iexit=1
772 : : call predictimg(delta_energy,imagealgo_str(dtset%imgmov),dtset%imgmov,itimimage,&
773 : : & itimimage_eff,list_dynimage,ga_param,mep_param,mpi_enreg,m1geo_param,dtset%natom,ndynimage,&
774 [ - + ]: 6 : & nimage,dtset%nimage,ntimimage_stored,pimd_param,dtset%prtvolimg,results_img)
775 : : endif
776 : :
777 [ - + ]: 4010 : ABI_FREE(occ)
778 [ - + ]: 4010 : ABI_FREE(vel)
779 [ - + ]: 4010 : ABI_FREE(xred)
780 [ - + ]: 4010 : ABI_FREE(list_dynimage)
781 : :
782 [ + + ]: 4010 : ABI_SFREE(amass)
783 : :
784 [ + + ]: 8361 : do itimimage=1,ntimimage_stored
785 : 8361 : call destroy_results_img(results_img(:,itimimage))
786 : : end do
787 [ - + ]: 4010 : ABI_FREE(results_img)
788 [ + + ]: 8280 : do iimage=1,nimage
789 : 8280 : call scf_history_free(scf_history(iimage))
790 : : end do
791 [ + - + + : 8280 : ABI_FREE(scf_history)
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
- - - - -
- - - - -
- - - - +
- - - - -
- - - - -
- - - - -
+ - - - -
- - - + ]
792 [ - + ]: 4010 : ABI_FREE(scf_initialized)
793 [ + + ]: 4010 : if (allocated(hist_prev)) then
794 : 1 : call abihist_free(hist_prev)
795 [ + - + + : 6 : ABI_FREE(hist_prev)
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + ]
796 : : end if
797 [ + + ]: 4010 : if (allocated(hist)) then
798 : 43 : call abihist_free(hist)
799 [ + - + + : 309 : ABI_FREE(hist)
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + ]
800 : : end if
801 : :
802 : 4010 : call mep_destroy(mep_param)
803 : 4010 : call ga_destroy(ga_param)
804 : 4010 : call m1geo_destroy(m1geo_param)
805 : 4010 : call pimd_destroy(pimd_param)
806 : :
807 : 4010 : call timab(1210,2,tsec)
808 : 4010 : call timab(1200,2,tsec)
809 : :
810 : : DBG_EXIT("COLL")
811 : :
812 [ - + - + : 8020 : end subroutine gstateimg
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- + - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - + ]
813 : : !!***
814 : :
815 : : !!****f* ABINIT/prtimg
816 : : !! NAME
817 : : !! prtimg
818 : : !!
819 : : !! FUNCTION
820 : : !! Print out results obtained by as ground-state calculation of
821 : : !! an image of the cell. The printing format is condensed in order
822 : : !! to facilitate the reading.
823 : : !!
824 : : !! INPUTS
825 : : !! dynimage(nimagetot)=flags defining static/dynamic state of images
826 : : !! imagealgo_str=name of the algorithm (with images) used
827 : : !! imgmov=index of algorithm (with images) used
828 : : !! iout=unit number for output
829 : : !! mpi_enreg=MPI-parallelisation information
830 : : !! nimage=number of images stored on current proc
831 : : !! nimage_tot=total number of images (should be dtset%nimage)
832 : : !! prt_all_images=true if all images have to be printed out (ignoring dynimage)
833 : : !! prtvolimg=printing volume for each image
834 : : !! <0 : nothing
835 : : !! 0 : only a title
836 : : !! 1 : energy, residuals, forces, stresses, velocities, atomic positions
837 : : !! 2 : energy, residuals
838 : : !! resimg(nimage) <type(results_img_type)>=results of the ground-state computations
839 : : !! for all images treated by current proc
840 : : !!
841 : : !! OUTPUT
842 : : !! (data written to unit iout)
843 : : !!
844 : : !! SOURCE
845 : :
846 : 131 : subroutine prtimg(dynimage,imagealgo_str,imgmov,iout,mpi_enreg,nimage,nimage_tot,&
847 : 131 : & prt_all_images,prtvolimg,resimg)
848 : :
849 : : !Arguments ------------------------------------
850 : : !scalars
851 : : integer,intent(in) :: nimage_tot,dynimage(nimage_tot),imgmov,iout,nimage,prtvolimg !vz_d
852 : : logical,intent(in) :: prt_all_images
853 : : character(len=60),intent(in) :: imagealgo_str
854 : : type(MPI_type),intent(in) :: mpi_enreg
855 : : !arrays
856 : : type(results_img_type),target,intent(inout) :: resimg(nimage)
857 : :
858 : : !Local variables-------------------------------
859 : : !scalars
860 : : integer :: ii,prtvel
861 : : logical :: test_img
862 : : real(dp) :: ucvol_img
863 : : character(len=500) :: msg
864 : : !arrays
865 : 131 : integer,allocatable :: iatfix_img(:,:)
866 : 131 : real(dp),allocatable :: gmet_img(:,:),gprimd_img(:,:),rmet_img(:,:),xcart_img(:,:)
867 : 131 : type(results_img_type),pointer :: resimg_all(:)
868 : :
869 : : ! ****************************************************************
870 : :
871 : : DBG_ENTER('COLL')
872 : :
873 : 131 : if (prtvolimg<=0) return
874 [ + + ]: 131 : if (mpi_enreg%me_cell/=0) return
875 : :
876 : : !Gather data
877 [ + - ]: 59 : if (prtvolimg==1.or.prtvolimg==2) then
878 [ + - - + ]: 59 : test_img=(nimage_tot/=1.and.mpi_enreg%paral_img==1)
879 : : if (test_img) then
880 [ # # ]: 0 : if (mpi_enreg%me==0) then
881 [ # # # # : 0 : ABI_MALLOC(resimg_all,(nimage_tot))
# # # # #
# # # # #
# # ]
882 : : end if
883 : : call gather_results_img(mpi_enreg,resimg,resimg_all,master=0,&
884 : 0 : & allgather=.false.,only_one_per_img=.true.)
885 : : else
886 : 59 : resimg_all => resimg
887 : : end if
888 : : end if
889 : :
890 : : !===== First option for the printing volume ===
891 [ + + + - ]: 59 : if (prtvolimg==1.and.mpi_enreg%me==0) then
892 : :
893 [ - + ]: 27 : prtvel=0;if (imgmov==0.or.imgmov==9.or.imgmov==10.or.imgmov==13) prtvel=1
894 : :
895 [ + + ]: 216 : do ii=1,nimage_tot
896 [ + + + + ]: 216 : if (dynimage(ii)==1.or.prt_all_images) then
897 : :
898 : : ! Title
899 : 141 : write(msg,'(6a,i4,a,i4,2a)') ch10,&
900 : 141 : '----------------------------------------------------------------------',ch10,&
901 : 141 : ' ',trim(imagealgo_str),' - CELL # ',ii,'/',nimage_tot,ch10,&
902 : 282 : '----------------------------------------------------------------------'
903 : 141 : call wrtout(iout,msg,'COLL')
904 : :
905 : : ! Total energy
906 : 141 : write(msg,'(2a,es20.12)') ch10,' Total energy for the cell [Ha]: ',resimg_all(ii)%results_gs%etotal
907 : 141 : call wrtout(iout,msg,'COLL')
908 : :
909 : : ! Residuals of the SCF cycle
910 : 141 : write(msg,'(3a,4(a,es16.8,a))') ch10,&
911 : 141 : ' Residuals from SCF cycle: ',ch10,&
912 : 141 : ' Total energy difference =',resimg_all(ii)%results_gs%deltae,ch10,&
913 : 141 : ' Maximal forces difference =',resimg_all(ii)%results_gs%diffor,ch10,&
914 : 141 : ' Max. residual of wave-functions=',resimg_all(ii)%results_gs%residm,ch10,&
915 : 282 : ' Density/potential residual (^2)=',resimg_all(ii)%results_gs%res2,ch10
916 : 141 : call wrtout(iout,msg,'COLL')
917 : :
918 : : ! Cell parameters
919 [ - + ]: 141 : ABI_MALLOC(rmet_img,(3,3))
920 [ - + ]: 141 : ABI_MALLOC(gmet_img,(3,3))
921 [ - + ]: 141 : ABI_MALLOC(gprimd_img,(3,3))
922 : 141 : call metric(gmet_img,gprimd_img,iout,rmet_img,resimg_all(ii)%rprim,ucvol_img)
923 [ - + ]: 141 : ABI_FREE(rmet_img)
924 [ - + ]: 141 : ABI_FREE(gmet_img)
925 [ - + ]: 141 : ABI_FREE(gprimd_img)
926 : :
927 : : ! Positions, forces and velocities
928 [ - + + - : 423 : ABI_MALLOC(iatfix_img,(3,resimg_all(ii)%natom))
+ - + - -
+ - + ]
929 [ - + + - : 423 : ABI_MALLOC(xcart_img,(3,resimg_all(ii)%natom))
+ - + - -
+ - + ]
930 [ + + + + ]: 1269 : iatfix_img=0
931 : 141 : call xred2xcart(resimg_all(ii)%natom,resimg_all(ii)%rprim,xcart_img,resimg_all(ii)%xred)
932 : : call prtxvf(resimg_all(ii)%results_gs%fcart,resimg_all(ii)%results_gs%gred,&
933 : : & iatfix_img,iout,resimg_all(ii)%natom,prtvel,&
934 : 141 : & resimg_all(ii)%vel,xcart_img,resimg_all(ii)%xred)
935 [ - + ]: 141 : ABI_FREE(iatfix_img)
936 [ - + ]: 141 : ABI_FREE(xcart_img)
937 : :
938 : : ! Stress tensor
939 : : write(msg, '(a,es12.4,a)' ) &
940 : 141 : & '-Cartesian components of stress tensor (GPa) [Pressure=',&
941 : : & -(resimg_all(ii)%results_gs%strten(1)+resimg_all(ii)%results_gs%strten(2) &
942 : 282 : & +resimg_all(ii)%results_gs%strten(3))*HaBohr3_GPa/three,' GPa]'
943 : 141 : call wrtout(iout,msg,'COLL')
944 : 141 : write(msg, '(2(a,1p,e16.8))' ) '- sigma(1 1)=',resimg_all(ii)%results_gs%strten(1)*HaBohr3_GPa,&
945 : 282 : & ' sigma(3 2)=',resimg_all(ii)%results_gs%strten(4)*HaBohr3_GPa
946 : 141 : call wrtout(iout,msg,'COLL')
947 : 141 : write(msg, '(2(a,1p,e16.8))' ) '- sigma(2 2)=',resimg_all(ii)%results_gs%strten(2)*HaBohr3_GPa,&
948 : 282 : & ' sigma(3 1)=',resimg_all(ii)%results_gs%strten(5)*HaBohr3_GPa
949 : 141 : call wrtout(iout,msg,'COLL')
950 : 141 : write(msg, '(2(a,1p,e16.8))' ) '- sigma(3 3)=',resimg_all(ii)%results_gs%strten(3)*HaBohr3_GPa,&
951 : 282 : & ' sigma(2 1)=',resimg_all(ii)%results_gs%strten(6)*HaBohr3_GPa
952 : 282 : call wrtout(iout,msg,'COLL')
953 : : end if
954 : : end do
955 : : end if
956 : :
957 : :
958 : : !===== 2nd option for the printing volume ===
959 [ + + + - ]: 59 : if (prtvolimg==2.and.mpi_enreg%me==0) then
960 : 32 : write(msg,'(a,1x,a)') ch10,'Cell Total_energy[Ha] deltae diffor residm res2'
961 : 32 : call wrtout(iout,msg,'COLL')
962 [ + + ]: 246 : do ii=1,nimage_tot
963 [ + + + + ]: 246 : if (dynimage(ii)==1.or.prt_all_images) then
964 : : write(msg,'(1x,i4,2x,es16.8,4(1x,es13.5))') &
965 : 166 : & ii,resimg_all(ii)%results_gs%etotal,resimg_all(ii)%results_gs%deltae,&
966 : 166 : & resimg_all(ii)%results_gs%diffor,resimg_all(ii)%results_gs%residm,&
967 : 332 : & resimg_all(ii)%results_gs%res2
968 : 166 : call wrtout(iout,msg,'COLL')
969 : : end if
970 : : end do
971 : : end if
972 : :
973 : : !=====
974 [ + - ]: 59 : if (prtvolimg==1.or.prtvolimg==2) then
975 [ - + - - ]: 59 : if (test_img.and.mpi_enreg%me==0) then
976 : 0 : call destroy_results_img(resimg_all)
977 [ # # ]: 0 : ABI_FREE(resimg_all)
978 : : end if
979 : : nullify(resimg_all)
980 : : end if
981 : :
982 : : DBG_EXIT('COLL')
983 : :
984 [ + - ]: 131 : end subroutine prtimg
985 : : !!***
986 : :
987 : : !!****f* ABINIT/predictimg
988 : : !! NAME
989 : : !! predictimg
990 : : !!
991 : : !! FUNCTION
992 : : !! Given the past history of images, predict the new set of images
993 : : !!
994 : : !! INPUTS
995 : : !! deltae=averaged energy difference used to control convergence over images
996 : : !! imagealgo_str=name of the algorithm (with images) used
997 : : !! imgmov=gives the algorithm to be used for prediction of new set of images
998 : : !! itimimage=time index for image propagation (itimimage+1 is to be predicted here)
999 : : !! itimimage_eff=time index in the history
1000 : : !! list_dynimage(nimage)=list of dynamical images. The non-dynamical ones will not change.
1001 : : !! Example : in the NEB method, or in the string method, one expect the two end images to be fixed.
1002 : : !! mep_param=several parameters for Minimal Energy Path (MEP) search
1003 : : !! mpi_enreg=MPI-parallelisation information
1004 : : !! natom= number of atoms
1005 : : !! ndynimage=number of dynamical images
1006 : : !! nimage=number of images (treated by current proc)
1007 : : !! nimage_tot=total number of images
1008 : : !! ntimimage_stored=number of time steps stored in the history
1009 : : !! pimd_param=several parameters for Path-Integral MD
1010 : : !! prtvolimg=printing volume
1011 : : !!
1012 : : !! OUTPUT
1013 : : !!
1014 : : !! SIDE EFFECTS
1015 : : !! results_img(ntimimage_stored,nimage)=datastructure that holds the history of previous computations.
1016 : : !! results_img(:,:)%acell(3)
1017 : : !! at input, history of the values of acell for all images
1018 : : !! at output, the predicted values of acell for all images
1019 : : !! results_img(:,:)%results_gs
1020 : : !! at input, history of the values of energies and forces for all images
1021 : : !! results_img(:,:)%rprim(3,3)
1022 : : !! at input, history of the values of rprim for all images
1023 : : !! at output, the predicted values of rprim for all images
1024 : : !! results_img(:,:)%vel(3,natom)
1025 : : !! at input, history of the values of vel for all images
1026 : : !! at output, the predicted values of vel for all images
1027 : : !! results_img(:,:)%vel_cell(3,3)
1028 : : !! at input, history of the values of vel_cell for all images
1029 : : !! at output, the predicted values of vel_cell for all images
1030 : : !! results_img(:,:)%xred(3,natom)
1031 : : !! at input, history of the values of xred for all images
1032 : : !! at output, the predicted values of xred for all images
1033 : : !!
1034 : : !! SOURCE
1035 : :
1036 : 608 : subroutine predictimg(deltae,imagealgo_str,imgmov,itimimage,itimimage_eff,list_dynimage,&
1037 : : & ga_param,mep_param,mpi_enreg,m1geo_param,natom,ndynimage,nimage,nimage_tot,&
1038 : 304 : & ntimimage_stored,pimd_param,prtvolimg,results_img)
1039 : :
1040 : : use m_results_gs , only : results_gs_type
1041 : : use m_predict_neb, only : predict_neb
1042 : : use m_predict_steepest, only : predict_steepest
1043 : : use m_predict_pimd, only : predict_pimd
1044 : : use m_predict_string, only : predict_string
1045 : :
1046 : : !Arguments ------------------------------------
1047 : : !scalars
1048 : : integer,intent(in) :: imgmov,itimimage,itimimage_eff,natom,ndynimage
1049 : : integer,intent(in) :: nimage,nimage_tot,ntimimage_stored,prtvolimg
1050 : : character(len=60),intent(in) :: imagealgo_str
1051 : : real(dp),intent(in) :: deltae
1052 : : type(mep_type),intent(inout) :: mep_param
1053 : : type(m1geo_type),intent(inout) :: m1geo_param
1054 : : type(ga_type),intent(inout) :: ga_param
1055 : : type(pimd_type),intent(inout) :: pimd_param
1056 : : type(MPI_type),intent(in) :: mpi_enreg
1057 : : !arrays
1058 : : integer,intent(in) :: list_dynimage(ndynimage)
1059 : : type(results_img_type) :: results_img(nimage,ntimimage_stored)
1060 : :
1061 : : !Local variables-------------------------------
1062 : : !scalars
1063 : : integer,save :: idum=5
1064 : : logical :: is_pimd
1065 : : character(len=500) :: msg
1066 : : ! *************************************************************************
1067 : :
1068 : 304 : is_pimd=(imgmov==9.or.imgmov==10.or.imgmov==13)
1069 : :
1070 : : !Write convergence info
1071 : 304 : write(msg,'(3a)') ch10,&
1072 : 608 : & '------------------------------------------------------------',ch10
1073 [ + + ]: 304 : if (prtvolimg<2) write(msg,'(5a)') trim(msg),' ',trim(imagealgo_str),':',ch10
1074 : :
1075 : : !Specific case of 4th-order RK algorithm
1076 [ - + ]: 304 : if (mep_param%mep_solver==4) then
1077 [ # # ]: 0 : if (mod(itimimage,4)==0) then
1078 : 0 : write(msg,'(4a)') trim(msg),&
1079 : 0 : & ' Fourth-order Runge-Kutta algorithm - final step',ch10
1080 [ # # ]: 0 : if (itimimage>4) write(msg,'(2a,es11.3,2a)') trim(msg),&
1081 : 0 : & ' Average[Abs(Etotal(t)-Etotal(t-dt))]=',deltae,' Hartree',ch10
1082 : 0 : write(msg,'(2a)') trim(msg),' Moving images of the cell...'
1083 : : else
1084 : 0 : write(msg,'(2a,i1,2a)') trim(msg),&
1085 : 0 : & ' Fourth-order Runge-Kutta algorithm - intermediate step ',mod(itimimage,4),ch10
1086 : 0 : write(msg,'(2a)') trim(msg),' Computing new intermediate positions...'
1087 : : end if
1088 [ + + ]: 304 : else if (is_pimd) then
1089 : :
1090 : : ! PIMD
1091 : 35 : write(msg,'(2a)') trim(msg),' Moving images of the cell...'
1092 : : else
1093 : :
1094 : : ! Other cases
1095 [ + + ]: 502 : if (itimimage>1) write(msg,'(2a,es11.3,2a)') trim(msg),&
1096 : 466 : & ' Average[Abs(Etotal(t)-Etotal(t-dt))]=',deltae,' Hartree',ch10
1097 : 269 : write(msg,'(2a)') trim(msg),' Moving images of the cell...'
1098 : :
1099 : : end if
1100 : :
1101 : : !Write the msg
1102 : : !Prevent writing if iexit==1, which at present only happens for imgmov==6 algo
1103 [ + + + + : 900 : if(imgmov/=6 .or. m1geo_param%iexit==0) call wrtout([std_out, ab_out] ,msg)
+ + ]
1104 : :
1105 : 1 : select case(imgmov)
1106 : :
1107 : : case(0)
1108 : : call predict_copy(itimimage_eff,list_dynimage,ndynimage,nimage,&
1109 : 1 : & ntimimage_stored,results_img)
1110 : :
1111 : : case(1)
1112 : : call predict_steepest(itimimage,itimimage_eff,list_dynimage,mep_param,natom,ndynimage,nimage,&
1113 : 10 : & ntimimage_stored,results_img)
1114 : :
1115 : : case(2)
1116 : : call predict_string(itimimage,itimimage_eff,list_dynimage,mep_param,mpi_enreg,natom,&
1117 : 182 : & ndynimage,nimage,nimage_tot,ntimimage_stored,results_img)
1118 : :
1119 : : case(4)
1120 : : call predict_ga(itimimage_eff,idum,ga_param,natom,nimage,&
1121 : 4 : & ntimimage_stored,results_img)
1122 : :
1123 : : case(5)
1124 : : call predict_neb(itimimage,itimimage_eff,list_dynimage,mep_param,mpi_enreg,natom,&
1125 : 42 : & ndynimage,nimage,nimage_tot,ntimimage_stored,results_img)
1126 : :
1127 : : case(6)
1128 : 30 : call move_1geo(itimimage_eff,m1geo_param,mpi_enreg,nimage,nimage_tot,ntimimage_stored,pimd_param,results_img)
1129 : :
1130 : : case(9, 10, 13)
1131 : : ! Path Integral Molecular Dynamics
1132 : : call predict_pimd(imgmov,itimimage,itimimage_eff,mpi_enreg,natom,nimage,nimage_tot,&
1133 [ + + + + : 304 : & ntimimage_stored,pimd_param,prtvolimg,results_img)
+ + + - ]
1134 : :
1135 : : case default
1136 : :
1137 : : end select
1138 : :
1139 : 304 : end subroutine predictimg
1140 : : !!***
1141 : :
1142 : : !!****f* ABINIT/predict_copy
1143 : : !! NAME
1144 : : !! predict_copy
1145 : : !!
1146 : : !! FUNCTION
1147 : : !! Given the past history of images, predict the new set of images.
1148 : : !! Here, simple copy of the previous image.
1149 : : !!
1150 : : !! INPUTS
1151 : : !! itimimage_eff=time index in the history
1152 : : !! list_dynimage(nimage)=list of dynamical images. The non-dynamical ones will not change.
1153 : : !! Example : in the NEB of string method, one expect the two end images to be fixed.
1154 : : !! ndynimage=number of dynamical images
1155 : : !! nimage=number of images
1156 : : !! ntimimage_stored=number of time steps stored in the history
1157 : : !!
1158 : : !! OUTPUT
1159 : : !!
1160 : : !! SIDE EFFECTS
1161 : : !! results_img(ntimimage_stored,nimage)=datastructure that holds the history of previous computations.
1162 : : !! results_img(:,:)%acell(3)
1163 : : !! at input, history of the values of acell for all images
1164 : : !! at output, the predicted values of acell for all images
1165 : : !! results_img(:,:)%results_gs
1166 : : !! at input, history of the values of energies and forces for all images
1167 : : !! results_img(:,:)%rprim(3,3)
1168 : : !! at input, history of the values of rprim for all images
1169 : : !! at output, the predicted values of rprim for all images
1170 : : !! results_img(:,:)%vel(3,natom)
1171 : : !! at input, history of the values of vel for all images
1172 : : !! at output, the predicted values of vel for all images
1173 : : !! results_img(:,:)%vel_cell(3,3)
1174 : : !! at input, history of the values of vel_cell for all images
1175 : : !! at output, the predicted values of vel_cell for all images
1176 : : !! results_img(:,:)%xred(3,natom)
1177 : : !! at input, history of the values of xred for all images
1178 : : !! at output, the predicted values of xred for all images
1179 : : !!
1180 : : !! SOURCE
1181 : :
1182 : 1 : subroutine predict_copy(itimimage_eff,list_dynimage,ndynimage,nimage,&
1183 : 1 : & ntimimage_stored,results_img)
1184 : :
1185 : : !Arguments ------------------------------------
1186 : : !scalars
1187 : : integer,intent(in) :: itimimage_eff,ndynimage,nimage,ntimimage_stored
1188 : : !arrays
1189 : : integer,intent(in) :: list_dynimage(ndynimage)
1190 : : type(results_img_type),intent(inout) :: results_img(nimage,ntimimage_stored)
1191 : :
1192 : : !Local variables-------------------------------
1193 : : !scalars
1194 : : integer :: idynimage,iimage,next_itimimage
1195 : : ! *************************************************************************
1196 : :
1197 : 1 : next_itimimage=itimimage_eff+1
1198 [ - + ]: 1 : if (next_itimimage>ntimimage_stored) next_itimimage=1
1199 : :
1200 [ + + ]: 3 : do idynimage=1,ndynimage
1201 : :
1202 : 2 : iimage=list_dynimage(idynimage)
1203 : :
1204 [ + + ]: 8 : results_img(iimage,next_itimimage)%acell(:) =results_img(iimage,itimimage_eff)%acell(:)
1205 [ + + + + ]: 26 : results_img(iimage,next_itimimage)%rprim(:,:) =results_img(iimage,itimimage_eff)%rprim(:,:)
1206 [ + + + + ]: 18 : results_img(iimage,next_itimimage)%vel(:,:) =results_img(iimage,itimimage_eff)%vel(:,:)
1207 [ + + + + ]: 26 : results_img(iimage,next_itimimage)%vel_cell(:,:)=results_img(iimage,itimimage_eff)%vel_cell(:,:)
1208 [ + + + + ]: 19 : results_img(iimage,next_itimimage)%xred(:,:) =results_img(iimage,itimimage_eff)%xred(:,:)
1209 : :
1210 : : end do ! idynimage
1211 : :
1212 : 304 : end subroutine predict_copy
1213 : : !!***
1214 : :
1215 : : !!****f* ABINIT/move_1geo
1216 : : !! NAME
1217 : : !! move_1geo
1218 : : !!
1219 : : !! FUNCTION
1220 : : !! This subroutine uses the forces, stresses and other results obtained for several images with one, common, geometry,
1221 : : !! weight them to deliver averaged forces, stresses, etc, and uses these to predict the next common geometry.
1222 : : !! All images must be dynamical.
1223 : : !! WARNING : at present, only forces are used, to change atomic positions. No change of cell geometry.
1224 : : !! Since this is not the PIMD algorithm, suppose ntimimage_stored=ntimimage, and itimimage=itimimage_eff.
1225 : : !!
1226 : : !! INPUTS
1227 : : !! itimimage_eff=time index in the history
1228 : : !! nimage=number of images
1229 : : !! ntimimage_stored=number of time steps stored in the history
1230 : : !! mpi_enreg=MPI-parallelisation information
1231 : : !! m1geo_param=parameters for the 1geo algorithms
1232 : : !! pimd_param=datastructure that contains all the parameters necessary to Path-Integral MD
1233 : : !!
1234 : : !! OUTPUT
1235 : : !!
1236 : : !! SIDE EFFECTS
1237 : : !! results_img(ntimimage_stored,nimage)=datastructure that holds the history of previous computations.
1238 : : !! results_img(:,:)%acell(3)
1239 : : !! at input, history of the values of acell for all images
1240 : : !! at output, the predicted values of acell for all images
1241 : : !! results_img(:,:)%results_gs
1242 : : !! at input, history of the values of energies and forces for all images
1243 : : !! results_img(:,:)%rprim(3,3)
1244 : : !! at input, history of the values of rprim for all images
1245 : : !! at output, the predicted values of rprim for all images
1246 : : !! results_img(:,:)%vel(3,natom)
1247 : : !! at input, history of the values of vel for all images
1248 : : !! at output, the predicted values of vel for all images
1249 : : !! results_img(:,:)%vel_cell(3,3)
1250 : : !! at input, history of the values of vel_cell for all images
1251 : : !! at output, the predicted values of vel_cell for all images
1252 : : !! results_img(:,:)%xred(3,natom)
1253 : : !! at input, history of the values of xred for all images
1254 : : !! at output, the predicted values of xred for all images
1255 : : !!
1256 : : !! SOURCE
1257 : :
1258 : 30 : subroutine move_1geo(itimimage_eff,m1geo_param,mpi_enreg,nimage,nimage_tot,ntimimage_stored,pimd_param,results_img)
1259 : :
1260 : : !Arguments ------------------------------------
1261 : : !scalars
1262 : : integer,intent(in) :: itimimage_eff,nimage,nimage_tot,ntimimage_stored
1263 : : type(MPI_type),intent(in) :: mpi_enreg
1264 : : type(m1geo_type),intent(inout) :: m1geo_param
1265 : : type(pimd_type),intent(in) :: pimd_param
1266 : : !arrays
1267 : : type(results_img_type),target,intent(inout) :: results_img(nimage,ntimimage_stored)
1268 : :
1269 : : !Local variables-------------------------------
1270 : : !scalars
1271 : : integer :: ihist,iimage,natom,next_itimimage,nspden,nsppol
1272 : : !integer :: iatom
1273 : : real(dp) :: deltae,diffor,etotal,entropy,fermie,res2,residm
1274 : : logical :: test_img
1275 : 30 : type(results_gs_type) :: results_gs_lincomb
1276 : : !arrays
1277 : : real(dp) :: acell(3),rprim(3,3),rprimd(3,3),strten(6),vel_cell(3,3)
1278 : 30 : real(dp),allocatable :: fcart(:,:),vel(:,:),xred(:,:)
1279 : : logical :: DEBUG=.FALSE.
1280 : 30 : type(results_img_type),pointer :: resimg_all(:)
1281 : : ! *************************************************************************
1282 : :
1283 : 30 : natom=m1geo_param%ab_mover%natom
1284 : 30 : ihist=m1geo_param%hist_1geo%ihist
1285 : :
1286 [ + - + - : 90 : ABI_MALLOC(fcart,(3,natom))
+ - - + -
+ ]
1287 [ - + + - : 60 : ABI_MALLOC(vel,(3,natom))
- + - + ]
1288 [ - + + - : 60 : ABI_MALLOC(xred,(3,natom))
- + - + ]
1289 : :
1290 : : !Of course, assume that the geometry parameters are the same for all images, so take them from the first one.
1291 [ + + + + ]: 342 : xred(:,:) =results_img(1,itimimage_eff)%xred(:,:)
1292 [ + + ]: 120 : acell(:) =results_img(1,itimimage_eff)%acell(:)
1293 [ + + + + ]: 390 : rprim(:,:) =results_img(1,itimimage_eff)%rprim(:,:)
1294 [ + + + + ]: 342 : vel(:,:) =results_img(1,itimimage_eff)%vel(:,:)
1295 [ + + + + ]: 390 : vel_cell(:,:)=results_img(1,itimimage_eff)%vel_cell(:,:)
1296 : :
1297 : 30 : call mkrdim(acell,rprim,rprimd)
1298 : :
1299 : : !Fill history with the values of xred, acell and rprimd
1300 : 30 : call var2hist(acell,m1geo_param%hist_1geo,natom,rprimd,xred,DEBUG)
1301 : :
1302 : : !Fill history with velocities and ionic kinetic energy
1303 : 30 : call vel2hist(m1geo_param%ab_mover%amass,m1geo_param%hist_1geo,vel,vel_cell)
1304 : 30 : m1geo_param%hist_1geo%time(ihist)=zero
1305 : :
1306 : : !In case of image parallelism, collect results accross processors
1307 [ + + - + ]: 30 : test_img=(nimage_tot/=1.and.mpi_enreg%paral_img==1)
1308 : : if (test_img) then
1309 [ # # # # : 0 : ABI_MALLOC(resimg_all,(nimage_tot))
# # # # #
# # # # #
# # ]
1310 : : call gather_results_img(mpi_enreg,results_img(1:nimage,itimimage_eff),resimg_all,&
1311 : 0 : & allgather=.true.,only_one_per_img=.false.)
1312 : : else
1313 : 30 : resimg_all => results_img(:,itimimage_eff)
1314 : : end if
1315 : :
1316 : : !Compute energy, entropy, fermie, forces and stresses for the 1geo : take the weighted average.
1317 : : !Compute maximum of deltae,diffor,res2,residm
1318 : 30 : etotal=zero
1319 : 30 : entropy=zero
1320 : 30 : fermie=zero
1321 [ + + + + ]: 342 : fcart(:,:)=zero
1322 : 30 : strten(:)=zero
1323 : 30 : deltae=zero
1324 : 30 : diffor=zero
1325 : 30 : res2=zero
1326 : 30 : residm=zero
1327 : :
1328 [ + + ]: 88 : do iimage=1,nimage_tot
1329 : 58 : etotal=etotal+resimg_all(iimage)%results_gs%etotal*m1geo_param%mixesimgf(iimage)
1330 : 58 : entropy=entropy+resimg_all(iimage)%results_gs%entropy*m1geo_param%mixesimgf(iimage)
1331 : 58 : fermie=fermie+resimg_all(iimage)%results_gs%fermie*m1geo_param%mixesimgf(iimage)
1332 [ + + + + ]: 738 : fcart(:,:)=fcart(:,:)+resimg_all(iimage)%results_gs%fcart(:,:)*m1geo_param%mixesimgf(iimage)
1333 [ + + ]: 406 : strten(:) =strten(:) +resimg_all(iimage)%results_gs%strten(:)*m1geo_param%mixesimgf(iimage)
1334 : 58 : if( deltae<resimg_all(iimage)%results_gs%deltae ) deltae=resimg_all(iimage)%results_gs%deltae
1335 : 58 : if( diffor<resimg_all(iimage)%results_gs%diffor ) diffor=resimg_all(iimage)%results_gs%diffor
1336 : 58 : if( res2<resimg_all(iimage)%results_gs%res2 ) res2=resimg_all(iimage)%results_gs%res2
1337 : 88 : if( residm<resimg_all(iimage)%results_gs%residm ) residm=resimg_all(iimage)%results_gs%residm
1338 : : enddo
1339 : :
1340 : : !Set up a results_gs datastructure with the linear combination of images
1341 : 30 : nspden=resimg_all(1)%results_gs%nspden
1342 : 30 : nsppol=resimg_all(1)%results_gs%nsppol
1343 : 30 : call results_gs_lincomb%init(natom,nspden,nsppol)
1344 : 30 : call resimg_all(1)%results_gs%copy(results_gs_lincomb)
1345 : 30 : results_gs_lincomb%etotal=etotal
1346 : 30 : results_gs_lincomb%entropy=entropy
1347 : 30 : results_gs_lincomb%fermie=fermie
1348 [ + - + - : 372 : results_gs_lincomb%fcart=fcart
- + - - +
+ + + ]
1349 [ + + ]: 210 : results_gs_lincomb%strten=strten
1350 : 30 : results_gs_lincomb%deltae=deltae
1351 : 30 : results_gs_lincomb%diffor=diffor
1352 : 30 : results_gs_lincomb%res2=res2
1353 : 30 : results_gs_lincomb%residm=residm
1354 : :
1355 : : !From now on, all procs contain the same information about the geometry, etotal, forces, stress, etc.
1356 : : !Nothing more needs to be transmitted, and resimg_all is not needed anymore.
1357 [ - + ]: 30 : if (test_img) then
1358 : 0 : call destroy_results_img(resimg_all)
1359 [ # # ]: 0 : ABI_FREE(resimg_all)
1360 : : end if
1361 : 30 : nullify(resimg_all)
1362 : :
1363 : : !Echo result_gs_lincomb
1364 : 30 : call results_gs_lincomb%yaml_write(ab_out, info="Linear combination of ground state results")
1365 : :
1366 : : !Destroy result_gs_lincomb
1367 : 30 : call results_gs_lincomb%free()
1368 : :
1369 : : !Store fcart and strten in hist_1geo
1370 [ + + + + ]: 342 : m1geo_param%hist_1geo%fcart(:,:,ihist)=fcart(:,:)
1371 [ + + ]: 210 : m1geo_param%hist_1geo%strten(:,ihist) =strten(:)
1372 : :
1373 : : !Store them in ab_xfh
1374 : : !THIS IS TO BE DONE !
1375 : :
1376 : : !Compute new atomic positions and cell characteristics in the single geometry
1377 : : call precpred_1geo(m1geo_param%ab_mover,&
1378 : : & m1geo_param%ab_xfh_1geo,&
1379 : : & m1geo_param%ab_mover%amu_curr,&
1380 : : & m1geo_param%deloc,&
1381 : : & m1geo_param%dt_chkdilatmx,&
1382 : : & mpi_enreg%comm_cell,&
1383 : : & m1geo_param%dilatmx,&
1384 : : & m1geo_param%filnam_ds4,&
1385 : : & m1geo_param%hist_1geo,&
1386 : : & m1geo_param%hmctt,&
1387 : : & m1geo_param%icycle,&
1388 : : & m1geo_param%iexit,&
1389 : : & itimimage_eff,&
1390 : : & m1geo_param%mttk_vars,&
1391 : : & m1geo_param%nctime,&
1392 : : & m1geo_param%ncycle,&
1393 : : & m1geo_param%nerr_dilatmx,&
1394 : : & m1geo_param%npsp,&
1395 : : & m1geo_param%ntime,&
1396 : : & pimd_param,&
1397 : : & m1geo_param%rprimd_orig,&
1398 : : & m1geo_param%skipcycle,&
1399 [ + - - - : 30 : & m1geo_param%usewvl)
- - - + ]
1400 : :
1401 : : !Retrieve the new positions, cell parameters [and velocities ?!]
1402 : 30 : call hist2var(acell,m1geo_param%hist_1geo,natom,rprimd,xred,DEBUG)
1403 : :
1404 : : !Store acell, rprim, xred and vel for the new iteration if relevant
1405 [ + + ]: 30 : if(m1geo_param%iexit==0)then
1406 : 24 : next_itimimage=itimimage_eff+1
1407 [ - + ]: 24 : if (next_itimimage>ntimimage_stored)then
1408 : 0 : ABI_ERROR('next_itimimage>ntimimage_stored')
1409 : : endif
1410 [ + + ]: 70 : do iimage=1,nimage
1411 [ + + + + ]: 558 : results_img(iimage,next_itimimage)%xred(:,:) =xred(:,:)
1412 [ + + ]: 184 : results_img(iimage,next_itimimage)%acell(:) =acell(:)
1413 [ + + + + ]: 622 : results_img(iimage,next_itimimage)%rprim(:,:) =rprim(:,:)
1414 : : ! WARNING : Should also store vel and vel_cell of course ...
1415 : : ! results_img(iimage,next_itimimage)%vel(:,:) =vel(:,:)
1416 : : ! results_img(iimage,next_itimimage)%vel_cell(:,:)=vel_cell(:,:)
1417 : : end do
1418 : : endif
1419 [ - + ]: 30 : ABI_FREE(fcart)
1420 [ - + ]: 30 : ABI_FREE(vel)
1421 [ - + ]: 30 : ABI_FREE(xred)
1422 : :
1423 [ - + - + : 60 : end subroutine move_1geo
- + - + -
+ - + - +
- + - + -
+ - + -
+ ]
1424 : : !!***
1425 : :
1426 : : end module m_gstateimg
1427 : : !!***
|
