TABLE OF CONTENTS


ABINIT/m_fit_data [ Modules ]

[ Top ] [ Modules ]

NAME

 m_fit_data

FUNCTION

COPYRIGHT

 Copyright (C) 2010-2018 ABINIT group (AM)
 This file is distributed under the terms of the
 GNU General Public Licence, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

SOURCE

19 #if defined HAVE_CONFIG_H
20 #include "config.h"
21 #endif
22 
23 #include "abi_common.h"
24 
25 module m_fit_data
26 
27  use defs_basis
28  use m_errors
29  use m_abicore
30 
31  use m_geometry,     only : metric
32  
33  implicit none

m_fit_data/fit_data_compute [ Functions ]

[ Top ] [ m_fit_data ] [ Functions ]

NAME

 fit_data_compute

FUNCTION

 Conpute the strain of each configuration.
 Compute the displacmeent of each configuration.
 Compute the variation of the displacement due to strain of each configuration.
 Compute fixed forces and stresse and get the standard deviation.
 Compute Sheppard and al Factors  \Omega^{2} see J.Chem Phys 136, 074103 (2012) [[cite:Sheppard2012]]

INPUTS

 eff_pot<type(effective_potential)> = effective potential
 hist<type(abihist)> = The history of the MD (or snapshot of DFT
 comm = MPI communicator
 verbose  = optional, flag for the verbose mode

OUTPUT

 fit_data<fit_data_type> = fit_data is now filled

PARENTS

      m_fit_polynomial_coeff

CHILDREN

SOURCE

309 subroutine fit_data_compute(fit_data,eff_pot,hist,comm,verbose)
310 
311  use m_strain,only : strain_type,strain_get
312  use m_effective_potential,only : effective_potential_type,effective_potential_evaluate
313  use m_effective_potential,only : effective_potential_getDisp
314  use m_abihist, only : abihist
315  use m_strain,only : strain_type,strain_get
316 
317 !This section has been created automatically by the script Abilint (TD).
318 !Do not modify the following lines by hand.
319 #undef ABI_FUNC
320 #define ABI_FUNC 'fit_data_compute'
321 !End of the abilint section
322 
323  implicit none
324   
325 !Arguments ------------------------------------
326 !scalars
327  integer,intent(in) :: comm
328  logical,optional,intent(in) :: verbose
329  !arrays
330  type(fit_data_type),intent(inout) :: fit_data
331  type(effective_potential_type),intent(in) :: eff_pot
332  type(abihist),intent(in) :: hist
333 !Local variables-------------------------------
334 !scalar
335  integer :: ii,itime,natom,ntime
336  real(dp):: energy
337  logical :: need_verbose
338 !arrays
339  character(len=500) :: message
340  real(dp) :: gmet(3,3),gprimd(3,3),rmet(3,3)
341  real(dp),allocatable :: energy_diff(:)
342  real(dp),allocatable :: du_delta(:,:,:,:),displacement(:,:,:),strain(:,:)
343  real(dp),allocatable :: fcart_diff(:,:,:),fred_fixed(:,:,:),fcart_fixed(:,:,:)
344  real(dp),allocatable :: strten_diff(:,:),strten_fixed(:,:),sqomega(:),ucvol(:)
345  type(strain_type) :: strain_t
346  type(training_set_type) :: ts
347 ! *************************************************************************
348 
349 !Initialisation of optional arguments
350  need_verbose = .TRUE.
351  if(present(verbose)) need_verbose = verbose
352 
353 !Checks
354  natom = eff_pot%supercell%natom
355  ntime = hist%mxhist
356  if(natom /= size(hist%xred,2))then
357    write(message, '(5a)' )&
358 &      'The number of atoms in the hist file does not correspond to the supercell.',ch10,&
359 &      'You should call the routine fit_polynomial_coeff_mapHistToRef before',ch10,&
360 &      'Action: Contact abinit group'
361    MSG_BUG(message)
362  end if
363 
364 !Allocation of temporary arrays
365  ABI_ALLOCATE(displacement,(3,natom,ntime))
366  ABI_ALLOCATE(du_delta,(6,3,natom,ntime))
367  ABI_ALLOCATE(energy_diff,(ntime))
368  ABI_ALLOCATE(fcart_fixed,(3,natom,ntime))
369  ABI_ALLOCATE(fcart_diff,(3,natom,ntime))
370  ABI_ALLOCATE(fred_fixed,(3,natom,ntime))
371  ABI_ALLOCATE(strain,(6,ntime))
372  ABI_ALLOCATE(strten_fixed,(6,ntime))
373  ABI_ALLOCATE(strten_diff,(6,ntime))
374  ABI_ALLOCATE(sqomega,(ntime))
375  ABI_ALLOCATE(ucvol,(ntime))
376 
377  displacement = zero 
378  du_delta = zero 
379  strain = zero; 
380  fcart_fixed  = zero 
381  fred_fixed  = zero
382  strain = zero 
383  strten_fixed = zero 
384  strten_diff = zero
385  sqomega = zero 
386  ucvol = zero
387 
388  do itime=1,ntime
389 !  Get strain
390    call strain_get(strain_t,rprim=eff_pot%supercell%rprimd,&
391 &                  rprim_def=hist%rprimd(:,:,itime),symmetrized=.FALSE.)
392    if (strain_t%name /= "reference")  then
393      do ii=1,3
394        strain(ii,itime) = strain_t%strain(ii,ii)
395      end do
396      strain(4,itime) = (strain_t%strain(2,3) + strain_t%strain(3,2))
397      strain(5,itime) = (strain_t%strain(3,1) + strain_t%strain(1,3))
398      strain(6,itime) = (strain_t%strain(2,1) + strain_t%strain(1,2))
399    else
400      strain(:,itime) = zero
401    end if
402 
403 !  Get displacement and du_delta
404    call effective_potential_getDisp(displacement(:,:,itime),du_delta(:,:,:,itime),natom,&
405 &                                   hist%rprimd(:,:,itime),eff_pot%supercell%rprimd,comm,&
406 &                                   xred_hist=hist%xred(:,:,itime),xcart_ref=eff_pot%supercell%xcart,&
407 &                                   compute_displacement=.TRUE.,compute_duDelta=.TRUE.)
408 
409 !  Get forces and stresses from harmonic part (fixed part)     
410    call effective_potential_evaluate(eff_pot,energy,fcart_fixed(:,:,itime),fred_fixed(:,:,itime),&
411 &                                    strten_fixed(:,itime),natom,hist%rprimd(:,:,itime),&
412 &                                    displacement=displacement(:,:,itime),&
413 &                                    du_delta=du_delta(:,:,:,itime),strain=strain(:,itime),&
414 &                                    compute_anharmonic=.true.,verbose=.FALSE.)
415    
416 !  Compute \Omega^{2} and ucvol for each time
417    call metric(gmet,gprimd,-1,rmet,hist%rprimd(:,:,itime),ucvol(itime))
418 !  Formula: sqomega(itime) = (((ucvol(itime)**(-2.))* ((natom)**(0.5)))**(-1.0/3.0))**2
419 !   Compact form:
420    sqomega(itime) = ((ucvol(itime)**(4.0/3.0)) / ((natom)**(1/3.0)))
421 
422 !  Compute the difference between History and model (fixed part)
423    fcart_diff(:,:,itime) =  hist%fcart(:,:,itime) - fcart_fixed(:,:,itime)
424    energy_diff(itime)    =  hist%etot(itime) - energy
425    strten_diff(:,itime)  =  hist%strten(:,itime) - strten_fixed(:,itime)
426  end do
427    
428 !Set the training set
429  call training_set_init(ts,displacement,du_delta,natom,ntime,strain,sqomega,ucvol)
430 !Set the fit_data
431  call fit_data_init(fit_data,energy_diff,fcart_diff,natom,ntime,strten_diff,ts)
432 
433 !Free space
434  call training_set_free(ts)
435  ABI_DEALLOCATE(displacement)
436  ABI_DEALLOCATE(du_delta)
437  ABI_DEALLOCATE(energy_diff)
438  ABI_DEALLOCATE(fcart_fixed)
439  ABI_DEALLOCATE(fcart_diff)
440  ABI_DEALLOCATE(fred_fixed)
441  ABI_DEALLOCATE(strain)
442  ABI_DEALLOCATE(strten_fixed)
443  ABI_DEALLOCATE(strten_diff)
444  ABI_DEALLOCATE(sqomega)
445  ABI_DEALLOCATE(ucvol)
446 
447 end subroutine fit_data_compute

m_fit_data/fit_data_free [ Functions ]

[ Top ] [ m_fit_data ] [ Functions ]

NAME

 fit_data_free

FUNCTION

 Free the fit_data datatype

INPUTS

 fit_data<fit_data_type> = fit_data to be free

OUTPUT

PARENTS

      m_fit_data,m_fit_polynomial_coeff

CHILDREN

SOURCE

242 subroutine fit_data_free(fit_data)
243 
244 
245 !This section has been created automatically by the script Abilint (TD).
246 !Do not modify the following lines by hand.
247 #undef ABI_FUNC
248 #define ABI_FUNC 'fit_data_free'
249 !End of the abilint section
250 
251  implicit none
252   
253 !Arguments ------------------------------------
254 !scalars
255 !arrays
256  type(fit_data_type),intent(inout) :: fit_data
257 !Local variables-------------------------------
258 !scalar
259 !arrays
260 ! *************************************************************************
261 
262 ! Reset integer values
263   fit_data%ntime = 0
264   fit_data%natom = 0
265 
266 ! Deallocate arrays  
267   if(allocated(fit_data%energy_diff)) then
268     ABI_DEALLOCATE(fit_data%energy_diff)
269   end if
270   if(allocated(fit_data%fcart_diff)) then
271     ABI_DEALLOCATE(fit_data%fcart_diff)
272   end if
273   if(allocated(fit_data%strten_diff)) then
274     ABI_DEALLOCATE(fit_data%strten_diff)
275   end if
276   call training_set_free(fit_data%training_set)
277    
278 end subroutine fit_data_free

m_fit_data/fit_data_init [ Functions ]

[ Top ] [ m_fit_data ] [ Functions ]

NAME

 fit_data_init

FUNCTION

 Initialize fit_data datatype

INPUTS

 energy_diff(3,natom,ntime) = Difference of energy between DFT calculation and 
                             fixed part of the model (more often harmonic part)
 fcart_diff(3,natom,ntime) = Difference of cartesian forces between DFT calculation and 
                             fixed part of the model (more often harmonic part)
 natom = Number of atoms
 ntime = Number of time (number of snapshot, number of md step...)
 strten_diff(6,natom) = Difference of stress tensor between DFT calculation and 
                        fixed part of the model (more often harmonic part)
 sqomega(ntime) =  Sheppard and al Factors \Omega^{2} see J.Chem Phys 136, 074103 (2012) [[cite:Sheppard2012]]
 ucvol(ntime) = Volume of the system for each time
 ts<training_set_type> = datatype with the information about the training set

OUTPUT

 fit_data<fit_data_type> = fit_data datatype to be initialized

PARENTS

      m_fit_data

CHILDREN

SOURCE

165 subroutine fit_data_init(fit_data,energy_diff,fcart_diff,natom,ntime,strten_diff,ts)
166 
167 !Arguments ------------------------------------
168 !scalars
169 
170 !This section has been created automatically by the script Abilint (TD).
171 !Do not modify the following lines by hand.
172 #undef ABI_FUNC
173 #define ABI_FUNC 'fit_data_init'
174 !End of the abilint section
175 
176  integer,intent(in) :: natom,ntime
177 !arrays
178  real(dp),intent(in) :: energy_diff(ntime),fcart_diff(3,natom,ntime)
179  real(dp),intent(in) :: strten_diff(6,ntime)
180  type(training_set_type),intent(in) :: ts
181  type(fit_data_type),intent(inout) :: fit_data
182 !Local variables-------------------------------
183 !scalar
184 !arrays
185  character(len=500) :: message
186 ! *************************************************************************
187 
188  if(natom /= ts%natom)then
189    write(message, '(a)')&
190 &      ' The number of atoms does not correspond to the training set'
191    MSG_BUG(message)
192  end if
193 
194  if(ntime /= ts%ntime)then
195    write(message, '(a)')&
196 &      ' The number of time does not correspond to the training set'
197    MSG_BUG(message)
198  end if
199  
200 !Free the output 
201  call fit_data_free(fit_data)
202 
203 !Set integer values 
204  fit_data%ntime = ntime
205  fit_data%natom = natom
206 
207 !allocate arrays
208  ABI_ALLOCATE(fit_data%fcart_diff,(3,natom,ntime))
209  fit_data%fcart_diff(:,:,:) = fcart_diff(:,:,:)
210  
211  ABI_ALLOCATE(fit_data%strten_diff,(6,ntime))
212  fit_data%strten_diff(:,:) = strten_diff(:,:)
213 
214  ABI_ALLOCATE(fit_data%energy_diff,(ntime))
215  fit_data%energy_diff(:) = energy_diff
216  
217  call training_set_init(fit_data%training_set,ts%displacement,ts%du_delta,&
218 &                       natom,ntime,ts%strain,ts%sqomega,ts%ucvol)
219  
220 end subroutine fit_data_init

m_fit_data/fit_data_type [ Types ]

[ Top ] [ m_fit_data ] [ Types ]

NAME

 fit_data_type

FUNCTION

SOURCE

 98  type, public :: fit_data_type
 99 
100    integer :: ntime
101 !   Number of time in the training set
102    
103    integer :: natom
104 !   Number of atoms in the training set
105 
106    real(dp),allocatable :: energy_diff(:)
107 !   energy(ntime)
108 !   Array with the diffence between energy from training set and energy from initial model.
109 !   The model constains only harmonic part
110 
111    real(dp),allocatable :: fcart_diff(:,:,:)
112 !   fcart_diff(3,natom,ntime)
113 !   Array with the diffence between cartesian forces from training set and forces from initial model.
114 !   The model constains only harmonic part
115 
116    real(dp),allocatable :: strten_diff(:,:)
117 !   strten_diff(6,ntime)
118 !   Array with the diffence between strain from training set and strain from initial model.
119 !   The model constains only harmonic part
120 
121    type(training_set_type) :: training_set
122 !    datatype with the informations of the training set
123    
124  end type fit_data_type
125  
126 !routine for fit_data
127  public :: fit_data_compute
128  public :: fit_data_init
129  public :: fit_data_free 

m_fit_data/training_set_free [ Functions ]

[ Top ] [ m_fit_data ] [ Functions ]

NAME

 training_set_free

FUNCTION

 Free the training_set datatype

INPUTS

 training_set<training_set_type> = training_set to be free

OUTPUT

PARENTS

      m_fit_data

CHILDREN

SOURCE

543 subroutine training_set_free(ts)
544 
545 
546 !This section has been created automatically by the script Abilint (TD).
547 !Do not modify the following lines by hand.
548 #undef ABI_FUNC
549 #define ABI_FUNC 'training_set_free'
550 !End of the abilint section
551 
552  implicit none
553   
554 !Arguments ------------------------------------
555 !scalars
556 !arrays
557  type(training_set_type),intent(inout) :: ts  
558 !Local variables-------------------------------
559 !scalar
560 !arrays
561 ! *************************************************************************
562 
563 ! Reset integer values
564   ts%ntime = 0
565   ts%natom = 0
566 
567 ! Deallocate arrays  
568   if(allocated(ts%displacement)) then
569     ABI_DEALLOCATE(ts%displacement)
570   end if
571   if(allocated(ts%du_delta)) then
572     ABI_DEALLOCATE(ts%du_delta)
573   end if
574   if(allocated(ts%strain)) then
575     ABI_DEALLOCATE(ts%strain)
576   end if
577   if(allocated(ts%sqomega))then
578     ABI_DEALLOCATE(ts%sqomega)
579   end if
580   if(allocated(ts%ucvol)) then
581     ABI_DEALLOCATE(ts%ucvol)
582   end if
583    
584 end subroutine training_set_free

m_fit_data/training_set_init [ Functions ]

[ Top ] [ m_fit_data ] [ Functions ]

NAME

 training_set_init

FUNCTION

 Initialize training_set datatype

INPUTS

 du_delta(6,3,natom,ntime)  = Variation to displacements wrt to the strain (Bohr)
 displacement(3,natom,ntime)= Atomic displacement wrt to the reference (Bohr)
 natom = number of atoms
 ntime = number of time step
 strain(6,ntime) = Strain
 sqomega =  Sheppard and al Factors \Omega^{2} see J.Chem Phys 136, 074103 (2012) [[cite:Sheppard2012]]
 ucvol(ntime) = Volume of the supercell for each time (Bohr^3)

OUTPUT

 ts<training_set_type> = training set to be initialized

PARENTS

      m_fit_data

CHILDREN

SOURCE

477 subroutine training_set_init(ts,displacement,du_delta,natom,ntime,strain,sqomega,ucvol)
478 
479 !Arguments ------------------------------------
480 !scalars
481 
482 !This section has been created automatically by the script Abilint (TD).
483 !Do not modify the following lines by hand.
484 #undef ABI_FUNC
485 #define ABI_FUNC 'training_set_init'
486 !End of the abilint section
487 
488  integer,intent(in) :: natom,ntime
489 !arrays
490  real(dp),intent(in) :: displacement(3,natom,ntime),du_delta(6,3,natom,ntime)
491  real(dp),intent(in) :: strain(6,ntime),sqomega(ntime),ucvol(ntime)
492  type(training_set_type),intent(inout) :: ts  
493 !Local variables-------------------------------
494 !scalar
495 !arrays
496 ! *************************************************************************
497 
498 !Free the output 
499  call training_set_free(ts)
500 
501 !Set integer values 
502  ts%ntime = ntime
503  ts%natom = natom
504 
505 !allocate arrays
506  ABI_ALLOCATE(ts%displacement,(3,natom,ntime))
507  ts%displacement(:,:,:) = displacement(:,:,:)
508  
509  ABI_ALLOCATE(ts%du_delta,(6,3,natom,ntime))
510  ts%du_delta(:,:,:,:) = du_delta(:,:,:,:)
511  
512  ABI_ALLOCATE(ts%strain,(6,ntime))
513  ts%strain(:,:) = strain(:,:)
514  
515  ABI_ALLOCATE(ts%sqomega,(ntime))
516  ts%sqomega(:) = sqomega(:)
517  
518  ABI_ALLOCATE(ts%ucvol,(ntime))
519  ts%ucvol(:) = ucvol(:)
520  
521 end subroutine training_set_init

m_fit_data/training_set_type [ Types ]

[ Top ] [ m_fit_data ] [ Types ]

NAME

 training_set_type

FUNCTION

 datatype for with all the information for the fit process
 This structure contains all the informations of the training set:
   - ntime
   - displacement
   - du_delta
   - strain
   - sqomega
   - ucvol

SOURCE

52  type, public :: training_set_type
53 
54    integer :: ntime
55 !    Number of time in the training set
56    
57    integer :: natom
58 !    Number of atoms in the training set
59    
60    real(dp), allocatable :: displacement(:,:,:)
61 !    displacement(3,natom,ntime)
62 !    displacement array, difference of the position in cartisian coordinates
63 !    between training set and reference
64 
65    real(dp), allocatable :: du_delta(:,:,:,:)
66 !    du_delta(6,3,natom,ntime)
67 !    du_delta array, variation of displacement wrt to the strain
68 
69    real(dp), allocatable :: strain(:,:)
70 !    strain(6,ntime)
71 !    strain array, strain in the training set
72 
73    real(dp), allocatable :: sqomega(:)
74 !    sqomega(ntime)
75 !    sqomega(itime) = (((ucvol(itime)**(-2.))* ((natom_sc)**(0.5)))**(-1.0/3.0))**2
76 
77    real(dp), allocatable :: ucvol(:)
78 !    ucvol(ntime)
79 !    ucvol array, volume of the cell in the training set
80 
81  end type training_set_type
82  
83 !routine for training_set
84  public :: training_set_init
85  public :: training_set_free