Line data Source code
1 : !!****m* ABINIT/m_fit_data
2 : !!
3 : !! NAME
4 : !! m_fit_data
5 : !!
6 : !! FUNCTION
7 : !!
8 : !! COPYRIGHT
9 : !! Copyright (C) 2010-2026 ABINIT group (AM)
10 : !! This file is distributed under the terms of the
11 : !! GNU General Public Licence, see ~abinit/COPYING
12 : !! or http://www.gnu.org/copyleft/gpl.txt .
13 : !! For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
14 : !!
15 : !! SOURCE
16 :
17 :
18 : #if defined HAVE_CONFIG_H
19 : #include "config.h"
20 : #endif
21 :
22 : #include "abi_common.h"
23 :
24 : module m_fit_data
25 :
26 : use defs_basis
27 : use m_errors
28 : use m_abicore
29 :
30 : use m_geometry, only : metric
31 :
32 : implicit none
33 : !!***
34 :
35 : !!****t* m_fit_data/training_set_type
36 : !! NAME
37 : !! training_set_type
38 : !!
39 : !! FUNCTION
40 : !! datatype for with all the information for the fit process
41 : !! This structure contains all the information of the training set:
42 : !! - ntime
43 : !! - displacement
44 : !! - du_delta
45 : !! - strain
46 : !! - sqomega
47 : !! - ucvol
48 : !!
49 : !! SOURCE
50 :
51 : type, public :: training_set_type
52 :
53 : integer :: ntime
54 : ! Number of time in the training set
55 :
56 : integer :: natom
57 : ! Number of atoms in the training set
58 :
59 : real(dp), allocatable :: displacement(:,:,:)
60 : ! displacement(3,natom,ntime)
61 : ! displacement array, difference of the position in cartisian coordinates
62 : ! between training set and reference
63 :
64 : real(dp), allocatable :: du_delta(:,:,:,:)
65 : ! du_delta(6,3,natom,ntime)
66 : ! du_delta array, variation of displacement wrt to the strain
67 :
68 : real(dp), allocatable :: strain(:,:)
69 : ! strain(6,ntime)
70 : ! strain array, strain in the training set
71 :
72 : real(dp), allocatable :: sqomega(:)
73 : ! sqomega(ntime)
74 : ! sqomega(itime) = (((ucvol(itime)**(-2.))* ((natom_sc)**(0.5)))**(-1.0/3.0))**2
75 :
76 : real(dp), allocatable :: ucvol(:)
77 : ! ucvol(ntime)
78 : ! ucvol array, volume of the cell in the training set
79 :
80 : end type training_set_type
81 :
82 : !routine for training_set
83 : public :: training_set_init
84 : public :: training_set_free
85 : !!***
86 :
87 : !----------------------------------------------------------------------
88 :
89 : !!****t* m_fit_data/fit_data_type
90 : !! NAME
91 : !! fit_data_type
92 : !!
93 : !! FUNCTION
94 : !!
95 : !! SOURCE
96 :
97 : type, public :: fit_data_type
98 :
99 : integer :: ntime
100 : ! Number of time in the training set
101 :
102 : integer :: natom
103 : ! Number of atoms in the training set
104 :
105 : real(dp),allocatable :: energy_diff(:)
106 : ! energy(ntime)
107 : ! Array with the diffence between energy from training set and energy from initial model.
108 : ! The model constains only harmonic part
109 :
110 : real(dp),allocatable :: fcart_diff(:,:,:)
111 : ! fcart_diff(3,natom,ntime)
112 : ! Array with the diffence between cartesian forces from training set and forces from initial model.
113 : ! The model constains only harmonic part
114 :
115 : real(dp),allocatable :: strten_diff(:,:)
116 : ! strten_diff(6,ntime)
117 : ! Array with the diffence between strain from training set and strain from initial model.
118 : ! The model constains only harmonic part
119 :
120 : type(training_set_type) :: training_set
121 : ! datatype with the information of the training set
122 :
123 : end type fit_data_type
124 :
125 : !routine for fit_data
126 : public :: fit_data_compute
127 : public :: fit_data_init
128 : public :: fit_data_free
129 : !!***
130 :
131 : CONTAINS !===========================================================================================
132 :
133 : !!****f* m_fit_data/fit_data_init
134 : !!
135 : !! NAME
136 : !! fit_data_init
137 : !!
138 : !! FUNCTION
139 : !! Initialize fit_data datatype
140 : !!
141 : !! INPUTS
142 : !! energy_diff(3,natom,ntime) = Difference of energy between DFT calculation and
143 : !! fixed part of the model (more often harmonic part)
144 : !! fcart_diff(3,natom,ntime) = Difference of cartesian forces between DFT calculation and
145 : !! fixed part of the model (more often harmonic part)
146 : !! natom = Number of atoms
147 : !! ntime = Number of time (number of snapshot, number of md step...)
148 : !! strten_diff(6,natom) = Difference of stress tensor between DFT calculation and
149 : !! fixed part of the model (more often harmonic part)
150 : !! sqomega(ntime) = Sheppard and al Factors \Omega^{2} see J.Chem Phys 136, 074103 (2012) [[cite:Sheppard2012]]
151 : !! ucvol(ntime) = Volume of the system for each time
152 : !! ts<training_set_type> = datatype with the information about the training set
153 : !!
154 : !! OUTPUT
155 : !! fit_data<fit_data_type> = fit_data datatype to be initialized
156 : !!
157 : !! SOURCE
158 :
159 42 : subroutine fit_data_init(fit_data,energy_diff,fcart_diff,natom,ntime,strten_diff,ts)
160 :
161 : !Arguments ------------------------------------
162 : !scalars
163 : integer,intent(in) :: natom,ntime
164 : !arrays
165 : real(dp),intent(in) :: energy_diff(ntime),fcart_diff(3,natom,ntime)
166 : real(dp),intent(in) :: strten_diff(6,ntime)
167 : type(training_set_type),intent(in) :: ts
168 : type(fit_data_type),intent(inout) :: fit_data
169 : !Local variables-------------------------------
170 : !scalar
171 : !arrays
172 : character(len=500) :: message
173 : ! *************************************************************************
174 :
175 42 : if(natom /= ts%natom)then
176 : write(message, '(a)')&
177 0 : & ' The number of atoms does not correspond to the training set'
178 0 : ABI_BUG(message)
179 : end if
180 :
181 42 : if(ntime /= ts%ntime)then
182 : write(message, '(a)')&
183 0 : & ' The number of time does not correspond to the training set'
184 0 : ABI_BUG(message)
185 : end if
186 :
187 : !Free the output
188 42 : call fit_data_free(fit_data)
189 :
190 : !Set integer values
191 42 : fit_data%ntime = ntime
192 42 : fit_data%natom = natom
193 :
194 : !allocate arrays
195 168 : ABI_MALLOC(fit_data%fcart_diff,(3,natom,ntime))
196 154763 : fit_data%fcart_diff(:,:,:) = fcart_diff(:,:,:)
197 :
198 126 : ABI_MALLOC(fit_data%strten_diff,(6,ntime))
199 6769 : fit_data%strten_diff(:,:) = strten_diff(:,:)
200 :
201 126 : ABI_MALLOC(fit_data%energy_diff,(ntime))
202 1003 : fit_data%energy_diff(:) = energy_diff
203 :
204 : call training_set_init(fit_data%training_set,ts%displacement,ts%du_delta,&
205 42 : & natom,ntime,ts%strain,ts%sqomega,ts%ucvol)
206 :
207 42 : end subroutine fit_data_init
208 : !!***
209 :
210 : !!****f* m_fit_data/fit_data_free
211 : !!
212 : !! NAME
213 : !! fit_data_free
214 : !!
215 : !! FUNCTION
216 : !! Free the fit_data datatype
217 : !!
218 : !! INPUTS
219 : !! fit_data<fit_data_type> = fit_data to be free
220 : !! OUTPUT
221 : !!
222 : !! SOURCE
223 :
224 82 : subroutine fit_data_free(fit_data)
225 :
226 : !Arguments ------------------------------------
227 : !scalars
228 : !arrays
229 : type(fit_data_type),intent(inout) :: fit_data
230 : !Local variables-------------------------------
231 : !scalar
232 : !arrays
233 : ! *************************************************************************
234 :
235 : ! Reset integer values
236 82 : fit_data%ntime = 0
237 82 : fit_data%natom = 0
238 :
239 : ! Deallocate arrays
240 82 : if(allocated(fit_data%energy_diff)) then
241 42 : ABI_FREE(fit_data%energy_diff)
242 : end if
243 82 : if(allocated(fit_data%fcart_diff)) then
244 42 : ABI_FREE(fit_data%fcart_diff)
245 : end if
246 82 : if(allocated(fit_data%strten_diff)) then
247 42 : ABI_FREE(fit_data%strten_diff)
248 : end if
249 82 : call training_set_free(fit_data%training_set)
250 :
251 82 : end subroutine fit_data_free
252 : !!***
253 :
254 : !!****f* m_fit_data/fit_data_compute
255 : !!
256 : !! NAME
257 : !! fit_data_compute
258 : !!
259 : !! FUNCTION
260 : !! Conpute the strain of each configuration.
261 : !! Compute the displacmeent of each configuration.
262 : !! Compute the variation of the displacement due to strain of each configuration.
263 : !! Compute fixed forces and stresse and get the standard deviation.
264 : !! Compute Sheppard and al Factors \Omega^{2} see J.Chem Phys 136, 074103 (2012) [[cite:Sheppard2012]]
265 : !!
266 : !! INPUTS
267 : !! eff_pot<type(effective_potential)> = effective potential
268 : !! hist<type(abihist)> = The history of the MD (or snapshot of DFT
269 : !! comm = MPI communicator
270 : !! verbose = optional, flag for the verbose mode
271 : !!
272 : !! OUTPUT
273 : !! fit_data<fit_data_type> = fit_data is now filled
274 : !!
275 : !! SOURCE
276 :
277 42 : subroutine fit_data_compute(fit_data,eff_pot,hist,comm,verbose)
278 :
279 : use m_strain,only : strain_type,strain_get
280 : use m_effective_potential,only : effective_potential_type,effective_potential_evaluate
281 : use m_effective_potential,only : effective_potential_getDisp
282 : use m_abihist, only : abihist
283 : use m_strain,only : strain_type,strain_get
284 :
285 : !Arguments ------------------------------------
286 : !scalars
287 : integer,intent(in) :: comm
288 : logical,optional,intent(in) :: verbose
289 : !arrays
290 : type(fit_data_type),intent(inout) :: fit_data
291 : type(effective_potential_type),intent(in) :: eff_pot
292 : type(abihist),intent(in) :: hist
293 : !Local variables-------------------------------
294 : !scalar
295 : integer :: ii,itime,natom,ntime
296 : real(dp):: energy
297 : logical :: need_verbose
298 : !arrays
299 : character(len=500) :: message
300 : real(dp) :: gmet(3,3),gprimd(3,3),rmet(3,3)
301 42 : real(dp),allocatable :: energy_diff(:)
302 42 : real(dp),allocatable :: du_delta(:,:,:,:),displacement(:,:,:),strain(:,:)
303 42 : real(dp),allocatable :: fcart_diff(:,:,:),gred_fixed(:,:,:),fcart_fixed(:,:,:)
304 42 : real(dp),allocatable :: strten_diff(:,:),strten_fixed(:,:),sqomega(:),ucvol(:)
305 : type(strain_type) :: strain_t
306 42 : type(training_set_type) :: ts
307 : ! *************************************************************************
308 :
309 : !Initialisation of optional arguments
310 42 : need_verbose = .TRUE.
311 : if(present(verbose)) need_verbose = verbose
312 :
313 : !Checks
314 42 : natom = eff_pot%supercell%natom
315 42 : ntime = hist%mxhist
316 42 : if(natom /= size(hist%xred,2))then
317 : write(message, '(5a)' )&
318 0 : & 'The number of atoms in the hist file does not correspond to the supercell.',ch10,&
319 0 : & 'You should call the routine fit_polynomial_coeff_mapHistToRef before',ch10,&
320 0 : & 'Action: Contact abinit group'
321 0 : ABI_BUG(message)
322 : end if
323 :
324 : !Allocation of temporary arrays
325 168 : ABI_MALLOC(displacement,(3,natom,ntime))
326 168 : ABI_MALLOC(du_delta,(6,3,natom,ntime))
327 126 : ABI_MALLOC(energy_diff,(ntime))
328 126 : ABI_MALLOC(fcart_fixed,(3,natom,ntime))
329 126 : ABI_MALLOC(fcart_diff,(3,natom,ntime))
330 126 : ABI_MALLOC(gred_fixed,(3,natom,ntime))
331 126 : ABI_MALLOC(strain,(6,ntime))
332 84 : ABI_MALLOC(strten_fixed,(6,ntime))
333 84 : ABI_MALLOC(strten_diff,(6,ntime))
334 84 : ABI_MALLOC(sqomega,(ntime))
335 84 : ABI_MALLOC(ucvol,(ntime))
336 :
337 154763 : displacement = zero
338 846683 : du_delta = zero
339 6769 : strain = zero;
340 154763 : fcart_fixed = zero
341 154763 : gred_fixed = zero
342 6769 : strain = zero
343 6769 : strten_fixed = zero
344 6769 : strten_diff = zero
345 1003 : sqomega = zero
346 1003 : ucvol = zero
347 :
348 1003 : do itime=1,ntime
349 : ! Get strain
350 : call strain_get(strain_t,rprim=eff_pot%supercell%rprimd,&
351 961 : & rprim_def=hist%rprimd(:,:,itime),symmetrized=.FALSE.)
352 961 : if (strain_t%name /= "reference") then
353 3844 : do ii=1,3
354 3844 : strain(ii,itime) = strain_t%strain(ii,ii)
355 : end do
356 961 : strain(4,itime) = (strain_t%strain(2,3) + strain_t%strain(3,2))
357 961 : strain(5,itime) = (strain_t%strain(3,1) + strain_t%strain(1,3))
358 961 : strain(6,itime) = (strain_t%strain(2,1) + strain_t%strain(1,2))
359 : else
360 0 : strain(:,itime) = zero
361 : end if
362 :
363 : ! Get displacement and du_delta
364 : call effective_potential_getDisp(displacement(:,:,itime),du_delta(:,:,:,itime),natom,&
365 : & hist%rprimd(:,:,itime),eff_pot%supercell%rprimd,comm,&
366 : & xred_hist=hist%xred(:,:,itime),xcart_ref=eff_pot%supercell%xcart,&
367 961 : & compute_displacement=.TRUE.,compute_duDelta=.TRUE.)
368 :
369 : ! Get forces and stresses from harmonic part (fixed part)
370 : call effective_potential_evaluate(eff_pot,energy,fcart_fixed(:,:,itime),gred_fixed(:,:,itime),&
371 : & strten_fixed(:,itime),natom,hist%rprimd(:,:,itime),&
372 : & displacement=displacement(:,:,itime),&
373 : & du_delta=du_delta(:,:,:,itime),strain=strain(:,itime),&
374 961 : & compute_anharmonic=.true.,verbose=.FALSE.)
375 :
376 : ! Compute \Omega^{2} and ucvol for each time
377 961 : call metric(gmet,gprimd,-1,rmet,hist%rprimd(:,:,itime),ucvol(itime))
378 : ! Formula: sqomega(itime) = (((ucvol(itime)**(-2.))* ((natom)**(0.5)))**(-1.0/3.0))**2
379 : ! Compact form:
380 961 : sqomega(itime) = ((ucvol(itime)**(4.0/3.0)) / ((natom)**(1/3.0)))
381 :
382 : ! Compute the difference between History and model (fixed part)
383 154721 : fcart_diff(:,:,itime) = hist%fcart(:,:,itime) - fcart_fixed(:,:,itime)
384 961 : energy_diff(itime) = hist%etot(itime) - energy
385 6769 : strten_diff(:,itime) = hist%strten(:,itime) - strten_fixed(:,itime)
386 : end do ! End Loop itime
387 :
388 : !Set the training set
389 42 : call training_set_init(ts,displacement,du_delta,natom,ntime,strain,sqomega,ucvol)
390 : !Set the fit_data
391 42 : call fit_data_init(fit_data,energy_diff,fcart_diff,natom,ntime,strten_diff,ts)
392 :
393 : !Free space
394 42 : call training_set_free(ts)
395 42 : ABI_FREE(displacement)
396 42 : ABI_FREE(du_delta)
397 42 : ABI_FREE(energy_diff)
398 42 : ABI_FREE(fcart_fixed)
399 42 : ABI_FREE(fcart_diff)
400 42 : ABI_FREE(gred_fixed)
401 42 : ABI_FREE(strain)
402 42 : ABI_FREE(strten_fixed)
403 42 : ABI_FREE(strten_diff)
404 42 : ABI_FREE(sqomega)
405 42 : ABI_FREE(ucvol)
406 :
407 42 : end subroutine fit_data_compute
408 : !!***
409 :
410 :
411 :
412 :
413 : !!****f* m_fit_data/training_set_init
414 : !!
415 : !! NAME
416 : !! training_set_init
417 : !!
418 : !! FUNCTION
419 : !! Initialize training_set datatype
420 : !!
421 : !! INPUTS
422 : !! du_delta(6,3,natom,ntime) = Variation to displacements wrt to the strain (Bohr)
423 : !! displacement(3,natom,ntime)= Atomic displacement wrt to the reference (Bohr)
424 : !! natom = number of atoms
425 : !! ntime = number of time step
426 : !! strain(6,ntime) = Strain
427 : !! sqomega = Sheppard and al Factors \Omega^{2} see J.Chem Phys 136, 074103 (2012) [[cite:Sheppard2012]]
428 : !! ucvol(ntime) = Volume of the supercell for each time (Bohr^3)
429 : !!
430 : !! OUTPUT
431 : !! ts<training_set_type> = training set to be initialized
432 : !!
433 : !! SOURCE
434 :
435 84 : subroutine training_set_init(ts,displacement,du_delta,natom,ntime,strain,sqomega,ucvol)
436 :
437 : !Arguments ------------------------------------
438 : !scalars
439 : integer,intent(in) :: natom,ntime
440 : !arrays
441 : real(dp),intent(in) :: displacement(3,natom,ntime),du_delta(6,3,natom,ntime)
442 : real(dp),intent(in) :: strain(6,ntime),sqomega(ntime),ucvol(ntime)
443 : type(training_set_type),intent(inout) :: ts
444 : !Local variables-------------------------------
445 : !scalar
446 : !arrays
447 : ! *************************************************************************
448 :
449 : !Free the output
450 84 : call training_set_free(ts)
451 :
452 : !Set integer values
453 84 : ts%ntime = ntime
454 84 : ts%natom = natom
455 :
456 : !allocate arrays
457 336 : ABI_MALLOC(ts%displacement,(3,natom,ntime))
458 309526 : ts%displacement(:,:,:) = displacement(:,:,:)
459 :
460 336 : ABI_MALLOC(ts%du_delta,(6,3,natom,ntime))
461 1693366 : ts%du_delta(:,:,:,:) = du_delta(:,:,:,:)
462 :
463 252 : ABI_MALLOC(ts%strain,(6,ntime))
464 13538 : ts%strain(:,:) = strain(:,:)
465 :
466 252 : ABI_MALLOC(ts%sqomega,(ntime))
467 2006 : ts%sqomega(:) = sqomega(:)
468 :
469 168 : ABI_MALLOC(ts%ucvol,(ntime))
470 2006 : ts%ucvol(:) = ucvol(:)
471 :
472 42 : end subroutine training_set_init
473 : !!***
474 :
475 : !!****f* m_fit_data/training_set_free
476 : !!
477 : !! NAME
478 : !! training_set_free
479 : !!
480 : !! FUNCTION
481 : !! Free the training_set datatype
482 : !!
483 : !! INPUTS
484 : !! training_set<training_set_type> = training_set to be free
485 : !! OUTPUT
486 : !!
487 : !! SOURCE
488 :
489 208 : subroutine training_set_free(ts)
490 :
491 : !Arguments ------------------------------------
492 : !scalars
493 : !arrays
494 : type(training_set_type),intent(inout) :: ts
495 : !Local variables-------------------------------
496 : !scalar
497 : !arrays
498 : ! *************************************************************************
499 :
500 : ! Reset integer values
501 208 : ts%ntime = 0
502 208 : ts%natom = 0
503 :
504 : ! Deallocate arrays
505 208 : if(allocated(ts%displacement)) then
506 84 : ABI_FREE(ts%displacement)
507 : end if
508 208 : if(allocated(ts%du_delta)) then
509 84 : ABI_FREE(ts%du_delta)
510 : end if
511 208 : if(allocated(ts%strain)) then
512 84 : ABI_FREE(ts%strain)
513 : end if
514 208 : if(allocated(ts%sqomega))then
515 84 : ABI_FREE(ts%sqomega)
516 : end if
517 208 : if(allocated(ts%ucvol)) then
518 84 : ABI_FREE(ts%ucvol)
519 : end if
520 :
521 208 : end subroutine training_set_free
522 : !!***
523 :
524 0 : end module m_fit_data
525 : !!***
|
