TABLE OF CONTENTS


ABINIT/m_abihist [ Modules ]

[ Top ] [ Modules ]

NAME

 m_abihist

FUNCTION

 This module contains definition the type abihist and its related routines

 Datatypes:

 * abihist: Historical record of atomic positions forces and cell parameters

 Subroutines:

 * abihist_init
 * abihist_free
 * abihist_bcast
 * abihist_compare
 * hist2var
 * var2hist
 * vel2hist

COPYRIGHT

 Copyright (C) 2001-2022 ABINIT group (XG, SE)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

SOURCE

30 #if defined HAVE_CONFIG_H
31 #include "config.h"
32 #endif
33 
34 #include "abi_common.h"
35 
36 module m_abihist
37 
38  use defs_basis
39  use m_abicore
40  use m_errors
41  use m_xmpi
42  use m_nctk
43 #if defined HAVE_NETCDF
44  use netcdf
45 #endif
46 
47  use m_geometry,  only : fcart2gred, xred2xcart
48 
49  implicit none
50 
51  private

m_abihist/abihist [ Types ]

[ Top ] [ m_abihist ] [ Types ]

NAME

 abihist

FUNCTION

 This type has several vectors, and index scalars to store
 a proper history of previous evaluations of forces and
 stresses, velocities, positions and energies.

 It contains:
 * mxhist                  : Maximum size of history
 * ihist                   : index of history
 * acell(3,mxhist)         : Acell
 * rprimd(3,3,mxhist)      : Rprimd
 * xred(3,natom,mxhist)    : Xred
 * fcart(3,natom,mxhist)   : Fcart
 * strten(6,mxhist)        : STRten
 * vel(3,natom,mxhist)     : Velocities of atoms
 * vel_cell(3,natom,mxhist): Velocities of cell
 * etot(mxhist)            : Electronic total Energy
 * ekin(mxhist)            : Ionic Kinetic Energy
 * entropy(mxhist)         : Entropy
 * time(mxhist)            : Time (or iteration number for GO)

NOTES

 The vectors are not allocated because in some cases
 not all the vectors are needed, in particular a history
 of stresses is only needed if optcell/=0, and a history
 of velocities is needed for ionmov==1

 Store acell, rprimd and strten even with optcell/=0
 represent a waste of 12x (dp)[Usually 8 Bytes] per
 iteration, the reason to store all the records is
 because some routines (eg bfgs.F90) uses the metric (gmet)
 for initialize the hessian and we need rprimd for that.

SOURCE

 94  type, public :: abihist
 95 
 96 ! scalars
 97 ! Index of the last element on all records
 98     integer :: ihist = 0
 99 ! Maximun size of the historical records
100     integer :: mxhist = 0
101 ! Booleans to know if some arrays are changing
102     logical :: isVused  ! If velocities are changing
103     logical :: isARused ! If Acell and Rprimd are changing
104 
105 ! arrays
106 ! Vector of (x,y,z)x(mxhist) values of cell dimensions
107     real(dp), allocatable :: acell(:,:)
108 ! Vector of (x,y,z)x(x,y,z)x(mxhist) values of primitive vectors
109     real(dp), allocatable :: rprimd(:,:,:)
110 ! Vector of (x,y,z)x(natom)x(mxhist) values of reduced coordinates
111     real(dp), allocatable :: xred(:,:,:)
112 ! Vector of (x,y,z)x(natom)x(mxhist) values of cartesian forces
113     real(dp), allocatable :: fcart(:,:,:)
114 ! Vector of (6)x(mxhist) values of stress tensor
115     real(dp), allocatable :: strten(:,:)
116 ! Vector of (x,y,z)x(natom)x(mxhist) values of atomic velocities
117     real(dp), allocatable :: vel(:,:,:)
118 ! Vector of (x,y,z)x(x,y,z)x(mxhist) values of cell velocities
119     real(dp), allocatable :: vel_cell(:,:,:)
120 ! Vector of (mxhist) values of electronic total energy
121     real(dp), allocatable :: etot(:)
122 ! Vector of (mxhist) values of ionic kinetic energy
123     real(dp), allocatable :: ekin(:)
124 ! Vector of (mxhist) values of Entropy
125     real(dp), allocatable :: entropy(:)
126 ! Vector of (mxhist) values of time (relevant for MD calculations)
127     real(dp), allocatable :: time(:)
128 
129  end type abihist
130 
131  public :: abihist_init             ! Initialize the object
132  public :: abihist_free             ! Destroy the object
133  public :: abihist_bcast            ! Broadcast the object
134  public :: abihist_copy             ! Copy 2 HIST records
135  public :: abihist_compare_and_copy ! Compare 2 HIST records; if similar copy
136  public :: hist2var                 ! Get xred, acell and rprimd from the history.
137  public :: abihist_findIndex        ! Shift history indexes
138  public :: var2hist                 ! Append xred, acell and rprimd
139  public :: vel2hist                 ! Append velocities and Kinetic Energy
140  public :: write_md_hist            ! Write the history into a netcdf file
141  public :: write_md_hist_img        ! Write the history into a netcdf file (with images)
142  public :: read_md_hist             ! Read the history from a netcdf file
143  public :: read_md_hist_img         ! Read the history from a netcdf file (with images)
144  public :: get_dims_hist
145  public :: read_csts_hist
146 
147  interface abihist_init
148    module procedure abihist_init_0D
149    module procedure abihist_init_1D
150  end interface abihist_init
151 
152  interface abihist_free
153    module procedure abihist_free_0D
154    module procedure abihist_free_1D
155  end interface abihist_free
156 
157  interface abihist_bcast
158    module procedure abihist_bcast_0D
159    module procedure abihist_bcast_1D
160  end interface abihist_bcast

m_abihist/abihist_bcast_0D [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 abihist_bcast_0D

FUNCTION

 Broadcast a hist datastructure (from a root process to all others) - Target: scalar

INPUTS

  master=ID of the sending node in comm
  comm=MPI Communicator

OUTPUT

SIDE EFFECTS

  hist <type(abihist)> = The hist to broadcast

SOURCE

351 subroutine abihist_bcast_0D(hist,master,comm)
352 
353 !Arguments ------------------------------------
354 !scalars
355  integer,intent(in) :: master,comm
356  type(abihist),intent(inout) :: hist
357 
358 !Local variables-------------------------------
359 !scalars
360  integer :: bufsize,ierr,indx,nproc,rank
361  integer :: sizeA,sizeA1,sizeA2
362  integer :: sizeEt,sizeEk,sizeEnt,sizeT
363  integer :: sizeR,sizeR1,sizeR2,sizeR3
364  integer :: sizeS,sizeS1,sizeS2
365  integer :: sizeV,sizeV1,sizeV2,sizeV3
366  integer :: sizeVc,sizeVc1,sizeVc2,sizeVc3
367  integer :: sizeX,sizeX1,sizeX2,sizeX3
368  integer :: sizeF,sizeF1,sizeF2,sizeF3
369 !arrays
370  integer,allocatable :: buffer_i(:)
371  real(dp),allocatable :: buffer_r(:)
372 
373 ! ***************************************************************
374 
375  ierr=0
376  nproc=xmpi_comm_size(comm)
377  if (nproc<=1) return
378 
379  rank=xmpi_comm_rank(comm)
380 
381 !=== Broadcast integers and logicals
382  ABI_MALLOC(buffer_i,(4))
383  if (rank==master) then
384    buffer_i(1)=hist%ihist
385    buffer_i(2)=hist%mxhist
386    buffer_i(3)=0;if (hist%isVused)  buffer_i(3)=1
387    buffer_i(4)=0;if (hist%isARused) buffer_i(4)=1
388  end if
389  call xmpi_bcast(buffer_i,master,comm,ierr)
390  if (rank/=master) then
391    hist%ihist=buffer_i(1)
392    hist%mxhist=buffer_i(2)
393    hist%isVused  =(buffer_i(3)==1)
394    hist%isARused =(buffer_i(4)==1)
395  end if
396  ABI_FREE(buffer_i)
397 
398 !If history is empty, return
399  if (hist%mxhist==0.or.hist%ihist==0) return
400 
401 !=== Broadcast sizes of arrays
402  ABI_MALLOC(buffer_i,(23))
403  if (rank==master) then
404    sizeA1=size(hist%acell,1);sizeA2=size(hist%acell,2)
405    sizeEt=size(hist%etot,1);sizeEk=size(hist%ekin,1);
406    sizeEnt=size(hist%entropy,1);sizeT=size(hist%time,1)
407    sizeR1=size(hist%rprimd,1);sizeR2=size(hist%rprimd,2);sizeR3=size(hist%rprimd,3)
408    sizeS1=size(hist%strten,1);sizeS2=size(hist%strten,2)
409    sizeV1=size(hist%vel,1);sizeV2=size(hist%vel,2);sizeV3=size(hist%vel,3)
410    sizeVc1=size(hist%vel_cell,1);sizeVc2=size(hist%vel_cell,2);sizeVc3=size(hist%vel_cell,3)
411    sizeX1=size(hist%xred,1);sizeX2=size(hist%xred,2);sizeX3=size(hist%xred,3)
412    sizeF1=size(hist%fcart,1);sizeF2=size(hist%fcart,2);sizeF3=size(hist%fcart,3)
413    buffer_i(1)=sizeA1  ;buffer_i(2)=sizeA2
414    buffer_i(3)=sizeEt  ;buffer_i(4)=sizeEk
415    buffer_i(5)=sizeEnt ;buffer_i(6)=sizeT
416    buffer_i(7)=sizeR1  ;buffer_i(8)=sizeR2
417    buffer_i(9)=sizeR3  ;buffer_i(10)=sizeS1
418    buffer_i(11)=sizeS2 ;buffer_i(12)=sizeV1
419    buffer_i(13)=sizeV2 ;buffer_i(14)=sizeV3
420    buffer_i(15)=sizeVc1;buffer_i(16)=sizeVc2
421    buffer_i(17)=sizeVc3;buffer_i(18)=sizeX1
422    buffer_i(19)=sizeX2 ;buffer_i(20)=sizeX3
423    buffer_i(21)=sizeF1 ;buffer_i(22)=sizeF2
424    buffer_i(23)=sizeF3
425  end if
426  call xmpi_bcast(buffer_i,master,comm,ierr)
427 
428  if (rank/=master) then
429    sizeA1 =buffer_i(1) ;sizeA2 =buffer_i(2)
430    sizeEt =buffer_i(3) ;sizeEk =buffer_i(4)
431    sizeEnt=buffer_i(5);sizeT   =buffer_i(6)
432    sizeR1 =buffer_i(7) ;sizeR2 =buffer_i(8)
433    sizeR3 =buffer_i(9) ;sizeS1 =buffer_i(10)
434    sizeS2 =buffer_i(11);sizeV1 =buffer_i(12)
435    sizeV2 =buffer_i(13);sizeV3 =buffer_i(14)
436    sizeVc1=buffer_i(15);sizeVc2=buffer_i(16)
437    sizeVc3=buffer_i(17);sizeX1 =buffer_i(18)
438    sizeX2 =buffer_i(19);sizeX3 =buffer_i(20)
439    sizeF1 =buffer_i(21);sizeF2 =buffer_i(22)
440    sizeF3 =buffer_i(23)
441  end if
442  ABI_FREE(buffer_i)
443 
444 !=== Broadcast reals
445  sizeA=sizeA1*sizeA2;sizeR=sizeR1*sizeR2*sizeR3;sizeS=sizeS1*sizeS2
446  sizeV=sizeV1*sizeV2*sizeV3;sizeVc=sizeVc1*sizeVc2*sizeVc3;
447  sizeX=sizeX1*sizeX2*sizeX3;sizeF=sizeF1*sizeF2*sizeF3
448  bufsize=sizeA+sizeEt+sizeEk+sizeEnt+sizeT+sizeR+sizeS+sizeV+sizeVc+sizeX+sizeF
449  ABI_MALLOC(buffer_r,(bufsize))
450  if (rank==master) then
451    indx=0
452    buffer_r(indx+1:indx+sizeA)=reshape(hist%acell(1:sizeA1,1:sizeA2),(/sizeA/))
453    indx=indx+sizeA
454    buffer_r(indx+1:indx+sizeEt)=hist%etot(1:sizeEt)
455    indx=indx+sizeEt
456    buffer_r(indx+1:indx+sizeEk)=hist%ekin(1:sizeEk)
457    indx=indx+sizeEk
458    buffer_r(indx+1:indx+sizeEnt)=hist%entropy(1:sizeEnt)
459    indx=indx+sizeEnt
460    buffer_r(indx+1:indx+sizeT)=hist%time(1:sizeT)
461    indx=indx+sizeT
462    buffer_r(indx+1:indx+sizeR)=reshape(hist%rprimd(1:sizeR1,1:sizeR2,1:sizeR3),(/sizeR/))
463    indx=indx+sizeR
464    buffer_r(indx+1:indx+sizeS)=reshape(hist%strten(1:sizeS1,1:sizeS2),(/sizeS/))
465    indx=indx+sizeS
466    buffer_r(indx+1:indx+sizeV)=reshape(hist%vel(1:sizeV1,1:sizeV2,1:sizeV3),(/sizeV/))
467    indx=indx+sizeV
468    buffer_r(indx+1:indx+sizeVc)=reshape(hist%vel_cell(1:sizeVc1,1:sizeVc2,1:sizeVc3),(/sizeVc/))
469    indx=indx+sizeVc
470    buffer_r(indx+1:indx+sizeX)=reshape(hist%xred(1:sizeX1,1:sizeX2,1:sizeX3),(/sizeX/))
471    indx=indx+sizeX
472    buffer_r(indx+1:indx+sizeF)=reshape(hist%fcart(1:sizeF1,1:sizeF2,1:sizeF3),(/sizeF/))
473  else
474    call abihist_free(hist)
475    ABI_MALLOC(hist%acell,(sizeA1,sizeA2))
476    ABI_MALLOC(hist%etot,(sizeEt))
477    ABI_MALLOC(hist%ekin,(sizeEk))
478    ABI_MALLOC(hist%entropy,(sizeEnt))
479    ABI_MALLOC(hist%time,(sizeT))
480    ABI_MALLOC(hist%rprimd,(sizeR1,sizeR2,sizeR3))
481    ABI_MALLOC(hist%strten,(sizeS1,sizeS2))
482    ABI_MALLOC(hist%vel,(sizeV1,sizeV2,sizeV3))
483    ABI_MALLOC(hist%vel_cell,(sizeVc1,sizeVc2,sizeVc3))
484    ABI_MALLOC(hist%xred,(sizeX1,sizeX2,sizeX3))
485    ABI_MALLOC(hist%fcart,(sizeF1,sizeF2,sizeF3))
486  end if
487  call xmpi_bcast(buffer_r,master,comm,ierr)
488 
489  if (rank/=master) then
490    indx=0
491    hist%acell(1:sizeA1,1:sizeA2)=reshape(buffer_r(indx+1:indx+sizeA), &
492 &                                       (/sizeA1,sizeA2/))
493    indx=indx+sizeA
494    hist%etot(1:sizeEt)=buffer_r(indx+1:indx+sizeEt)
495    indx=indx+sizeEt
496    hist%ekin(1:sizeEk)=buffer_r(indx+1:indx+sizeEk)
497    indx=indx+sizeEk
498    hist%entropy(1:sizeEnt)=buffer_r(indx+1:indx+sizeEnt)
499    indx=indx+sizeEnt
500    hist%time(1:sizeT)=buffer_r(indx+1:indx+sizeT)
501    indx=indx+sizeT
502    hist%rprimd(1:sizeR1,1:sizeR2,1:sizeR3)=reshape(buffer_r(indx+1:indx+sizeR), &
503 &                                                 (/sizeR1,sizeR2,sizeR3/))
504    indx=indx+sizeR
505    hist%strten(1:sizeS1,1:sizeS2)=reshape(buffer_r(indx+1:indx+sizeS), &
506 &                                        (/sizeS1,sizeS2/))
507    indx=indx+sizeS
508    hist%vel(1:sizeV1,1:sizeV2,1:sizeV3)=reshape(buffer_r(indx+1:indx+sizeV), &
509 &                                              (/sizeV1,sizeV2,sizeV3/))
510    indx=indx+sizeV
511    hist%vel_cell(1:sizeVc1,1:sizeVc2,1:sizeVc3)=reshape(buffer_r(indx+1:indx+sizeVc), &
512 &                                                      (/sizeVc1,sizeVc2,sizeVc3/))
513    indx=indx+sizeVc
514    hist%xred(1:sizeX1,1:sizeX2,1:sizeX3)=reshape(buffer_r(indx+1:indx+sizeX), &
515 &                                               (/sizeX1,sizeX2,sizeX3/))
516    indx=indx+sizeX
517    hist%fcart(1:sizeF1,1:sizeF2,1:sizeF3)=reshape(buffer_r(indx+1:indx+sizeF), &
518 &                                                (/sizeF1,sizeF2,sizeF3/))
519  end if
520  ABI_FREE(buffer_r)
521 
522 end subroutine abihist_bcast_0D

m_abihist/abihist_bcast_1D [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 abihist_bcast_1D

FUNCTION

 Broadcast a hist datastructure (from a root process to all others) - Target: 1D array

INPUTS

  master=ID of the sending node in comm
  comm=MPI Communicator

OUTPUT

SIDE EFFECTS

  hist(:) <type(abihist)> = The hist to broadcast

SOURCE

545 subroutine abihist_bcast_1D(hist,master,comm)
546 
547 !Arguments ------------------------------------
548 !scalars
549  integer,intent(in) :: master,comm
550  type(abihist),intent(inout) :: hist(:)
551 
552 !Local variables-------------------------------
553  integer :: ii
554 
555 ! ***************************************************************
556 
557  do ii=1,size(hist)
558    call abihist_bcast_0D(hist(ii),master,comm)
559  end do
560 
561 end subroutine abihist_bcast_1D

m_abihist/abihist_compare_and_copy [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 abihist_compare

FUNCTION

 Compare 2 HIST records

INPUTS

  hist_in <type(abihist)>
  tolerance
  store_all = flag to know if we need to increment ihist (store all the history)
              or just call shift (store just the last step)

OUTPUT

  similar= 1 the records are consistent
           0 the records are not consistent

SIDE EFFECTS

  hist_out <type(abihist)>

SOURCE

891 subroutine abihist_compare_and_copy(hist_in,hist_out,natom,similar,tolerance,store_all)
892 
893 !Arguments ------------------------------------
894 !scalars
895 integer,intent(in) :: natom
896 integer,intent(out) :: similar
897 real(dp),intent(in) :: tolerance
898 type(abihist),intent(in) :: hist_in
899 type(abihist),intent(inout) :: hist_out
900 logical,intent(in) :: store_all
901 !Local variables-------------------------------
902 !scalars
903 integer :: kk,jj
904 real(dp) :: maxdiff,diff
905 real(dp) :: x,y
906 !array
907 character(len= 500) :: msg
908 ! ***************************************************************
909 
910  ABI_UNUSED(store_all)
911 
912  similar=1
913 
914  write(msg,'(a,I0,4a)')  'Using values from history, iteration:',hist_in%ihist,ch10,&
915 &                     'Differences between present history and values stored',ch10,&
916 &                     'on the previous history.(Relative difference)'
917  call wrtout(std_out,msg,'COLL')
918 
919  x=hist_out%xred(1,1,hist_out%ihist)
920  y=hist_in%xred(1,1,hist_in%ihist)
921  maxdiff=2*abs(x-y)/(abs(x)+abs(y))
922  do kk=1,natom
923    do jj=1,3
924      x=hist_out%xred(jj,kk,hist_out%ihist)
925      y=hist_in%xred(jj,kk,hist_in%ihist)
926      diff=2*abs(x-y)/(abs(x)+abs(y))
927      if (diff>maxdiff) maxdiff=diff
928    end do
929  end do
930  write(msg,'(a,e12.5)') 'xred:     ',maxdiff
931  call wrtout(std_out,msg,'COLL')
932 
933  if (maxdiff>tolerance) similar=0
934 
935 
936  x=hist_out%rprimd(1,1,hist_out%ihist)
937  y=hist_in%rprimd(1,1,hist_in%ihist)
938  maxdiff=2*abs(x-y)/(abs(x)+abs(y))
939  do kk=1,3
940    do jj=1,3
941      x=hist_out%rprimd(jj,kk,hist_out%ihist)
942      y=hist_in%rprimd(jj,kk,hist_in%ihist)
943      diff=2*abs(x-y)/(abs(x)+abs(y))
944      if (diff>maxdiff) maxdiff=diff
945    end do
946  end do
947  write(msg,'(a,e12.5)') 'rprimd:   ',maxdiff
948  call wrtout(std_out,msg,'COLL')
949  if (maxdiff>tolerance) similar=0
950 
951 
952  x=hist_out%acell(1,hist_out%ihist)
953  y=hist_in%acell(1,hist_in%ihist)
954  maxdiff=2*abs(x-y)/(abs(x)+abs(y))
955  do kk=1,3
956    x=hist_out%acell(kk,hist_out%ihist)
957    y=hist_in%acell(kk,hist_in%ihist)
958    diff=2*abs(x-y)/(abs(x)+abs(y))
959    if (diff>maxdiff) maxdiff=diff
960  end do
961  write(msg,'(a,e12.5)') 'acell:    ',maxdiff
962  call wrtout(std_out,msg,'COLL')
963  if (maxdiff>tolerance) similar=0
964 
965  if (similar==1) then
966    hist_out%acell(:,hist_out%ihist)     =hist_in%acell(:,hist_in%ihist)
967    hist_out%rprimd(:,:,hist_out%ihist)  =hist_in%rprimd(:,:,hist_in%ihist)
968    hist_out%xred(:,:,hist_out%ihist)    =hist_in%xred(:,:,hist_in%ihist)
969    hist_out%fcart(:,:,hist_out%ihist)   =hist_in%fcart(:,:,hist_in%ihist)
970    hist_out%strten(:,hist_out%ihist)    =hist_in%strten(:,hist_in%ihist)
971    hist_out%vel(:,:,hist_out%ihist)     =hist_in%vel(:,:,hist_in%ihist)
972    hist_out%vel_cell(:,:,hist_out%ihist)=hist_in%vel_cell(:,:,hist_in%ihist)
973    hist_out%etot(hist_out%ihist)        =hist_in%etot(hist_in%ihist)
974    hist_out%ekin(hist_out%ihist)        =hist_in%ekin(hist_in%ihist)
975    hist_out%entropy(hist_out%ihist)     =hist_in%entropy(hist_in%ihist)
976    hist_out%time(hist_out%ihist)        =hist_in%time(hist_in%ihist)
977  end if
978 
979 end subroutine abihist_compare_and_copy

m_abihist/abihist_copy [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 abihist_copy

FUNCTION

 Copy one HIST record in another

INPUTS

  hist_in <type(abihist)>

OUTPUT

SIDE EFFECTS

  hist_out <type(abihist)>

SOURCE

829 subroutine abihist_copy(hist_in,hist_out)
830 
831 !Arguments ------------------------------------
832 !scalars
833 type(abihist),intent(in) :: hist_in
834 type(abihist),intent(inout) :: hist_out
835 
836 !Local variables-------------------------------
837 !scalars
838 ! character(len=500) :: msg
839 
840 ! ***************************************************************
841 
842 !Check
843  if (size(hist_in%xred,2)/=size(hist_out%xred,2)) then
844    ABI_BUG('Incompatible sizes for hist_in and hist_out!')
845  end if
846 
847 !Copy scalars (except ihist and mxhist)
848  hist_out%isVused  =hist_in%isVused
849  hist_out%isARused =hist_in%isARused
850 
851 !Copy arrays
852  hist_out%acell(:,hist_out%ihist)     = hist_in%acell(:,hist_in%ihist)
853  hist_out%rprimd(:,:,hist_out%ihist)  = hist_in%rprimd(:,:,hist_in%ihist)
854  hist_out%xred(:,:,hist_out%ihist)    = hist_in%xred(:,:,hist_in%ihist)
855  hist_out%fcart(:,:,hist_out%ihist)   = hist_in%fcart(:,:,hist_in%ihist)
856  hist_out%strten(:,hist_out%ihist)    = hist_in%strten(:,hist_in%ihist)
857  hist_out%vel(:,:,hist_out%ihist)     = hist_in%vel(:,:,hist_in%ihist)
858  hist_out%vel_cell(:,:,hist_out%ihist)= hist_in%vel_cell(:,:,hist_in%ihist)
859  hist_out%etot(hist_out%ihist)        = hist_in%etot(hist_in%ihist)
860  hist_out%ekin(hist_out%ihist)        = hist_in%ekin(hist_in%ihist)
861  hist_out%entropy(hist_out%ihist)     = hist_in%entropy(hist_in%ihist)
862  hist_out%time(hist_out%ihist)        = hist_in%time(hist_in%ihist)
863 
864 end subroutine abihist_copy

m_abihist/abihist_findIndex [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 abihist_findIndex

FUNCTION

INPUTS

 hist<type abihist>=Historical record of positions, forces, acell, stresses, and energies
 step = value of the needed step

OUTPUT

 index = index of the step in the hist file

SOURCE

643 function abihist_findIndex(hist,step) result(index)
644 
645 !Arguments ------------------------------------
646 !scalars
647  integer,intent(in) :: step
648  integer :: index
649 !arrays
650  type(abihist),intent(in) :: hist
651 !Local variables-------------------------------
652 !scalars
653  integer :: ii,mxhist
654 !arrays
655  character(len=500) :: msg
656 ! *************************************************************
657 
658  mxhist = hist%mxhist
659 
660  if ((mxhist ==1.and.step/=+1) .or. (mxhist /=1.and.abs(step) >=mxhist)) then
661    write(msg,'(a,I0,2a)')' The requested step must be less than ',mxhist,ch10,&
662                          'Action: increase the number of history stored in the history'
663    ABI_BUG(msg)
664  end if
665 
666  ii = hist%ihist + step
667 
668  do while (ii > mxhist)
669    ii = ii - mxhist
670  end do
671  do while (ii <= 0)
672    ii = ii + mxhist
673  end do
674 
675  index = ii
676 
677 end function abihist_findIndex

m_abihist/abihist_free_0D [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 abihist_free_0D

FUNCTION

 Deallocate dynamic memory in a hist structure - Target: scalar

SOURCE

281 subroutine abihist_free_0D(hist)
282 
283 !Arguments ------------------------------------
284  type(abihist),intent(inout) :: hist
285 
286 ! ***************************************************************
287 
288  ABI_SFREE(hist%acell)
289  ABI_SFREE(hist%rprimd)
290  ABI_SFREE(hist%xred)
291  ABI_SFREE(hist%fcart)
292  ABI_SFREE(hist%strten)
293  ABI_SFREE(hist%vel)
294  ABI_SFREE(hist%vel_cell)
295  ABI_SFREE(hist%etot)
296  ABI_SFREE(hist%ekin)
297  ABI_SFREE(hist%entropy)
298  ABI_SFREE(hist%time)
299 
300 end subroutine abihist_free_0D

m_abihist/abihist_free_1D [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 abihist_free_1D

FUNCTION

 Deallocate dynamic memory in a hist structure - Target: 1D array

SOURCE

314 subroutine abihist_free_1D(hist)
315 
316 !Arguments ------------------------------------
317  type(abihist),intent(inout) :: hist(:)
318 
319 !Local variables-------------------------------
320  integer :: ii
321 
322 ! ***************************************************************
323 
324  do ii=1,size(hist)
325    call abihist_free_0D(hist(ii))
326  end do
327 
328 end subroutine abihist_free_1D

m_abihist/abihist_init_0D [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 abihist_init_0D

FUNCTION

 Initialize a hist structure - Target: scalar

INPUTS

  natom = Number of atoms per unit cell
  mxhist = Maximal number of records to store
  isVUsed,isARUsed=flags used to initialize hsit structure

OUTPUT

  hist <type(abihist)> = The hist to initialize

SOURCE

185 subroutine abihist_init_0D(hist,natom,mxhist,isVused,isARused)
186 
187 !Arguments ------------------------------------
188  integer,intent(in) :: natom,mxhist
189  logical,intent(in) :: isVUsed,isARused
190  type(abihist),intent(inout) :: hist
191 
192 ! ***************************************************************
193 
194 !Initialize indexes
195  hist%ihist=1
196  hist%mxhist=mxhist
197 
198 !Initialize flags
199  hist%isVused=isVUsed
200  hist%isARused=isARUsed
201 
202 
203 !Allocate all the histories
204  ABI_MALLOC(hist%acell,(3,mxhist))
205  ABI_MALLOC(hist%rprimd,(3,3,mxhist))
206  ABI_MALLOC(hist%xred,(3,natom,mxhist))
207  ABI_MALLOC(hist%fcart,(3,natom,mxhist))
208  ABI_MALLOC(hist%strten,(6,mxhist))
209  ABI_MALLOC(hist%vel,(3,natom,mxhist))
210  ABI_MALLOC(hist%vel_cell,(3,3,mxhist))
211  ABI_MALLOC(hist%etot,(mxhist))
212  ABI_MALLOC(hist%ekin,(mxhist))
213  ABI_MALLOC(hist%entropy,(mxhist))
214  ABI_MALLOC(hist%time,(mxhist))
215 
216  hist%etot(1)=zero
217  hist%ekin(1)=zero
218  hist%entropy(1)=zero
219  hist%time(1)=zero
220 
221  hist%acell(:,1)=zero
222  hist%rprimd(:,:,1)=zero
223  hist%xred(:,:,1)=zero
224  hist%fcart(:,:,1)=zero
225  hist%strten(:,1)=zero
226  hist%vel(:,:,1)=zero
227  hist%vel_cell(:,:,1)=zero
228 
229 end subroutine abihist_init_0D

m_abihist/abihist_init_1D [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 abihist_init_1D

FUNCTION

 Initialize a hist structure - Target: 1D array

INPUTS

  natom = Number of atoms per unitary cell
  mxhist = Maximal number of records to store
  isVUsed,isARUsed=flags used to initialize hsit structure

OUTPUT

  hist(:) <type(abihist)> = The hist to initialize

SOURCE

251 subroutine abihist_init_1D(hist,natom,mxhist,isVUsed,isARUsed)
252 
253 !Arguments ------------------------------------
254  integer,intent(in) :: natom,mxhist
255  logical,intent(in) :: isVUsed,isARUsed
256  type(abihist),intent(inout) :: hist(:)
257 
258 !Local variables-------------------------------
259  integer :: ii
260 
261 ! ***************************************************************
262 
263  do ii=1,size(hist)
264    call abihist_init_0D(hist(ii),natom,mxhist,isVUsed,isARUsed)
265  end do
266 
267 end subroutine abihist_init_1D

m_abihist/def_file_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 def_file_hist

FUNCTION

INPUTS

OUTPUT

SOURCE

1475 subroutine def_file_hist(ncid,natom,nimage,ntypat,npsp,has_nimage)
1476 
1477 !Arguments ------------------------------------
1478 !scalars
1479  integer,intent(in) :: ncid
1480  integer,intent(in) :: natom,nimage,ntypat,npsp
1481  logical,intent(in) :: has_nimage
1482 
1483 !Local variables-------------------------------
1484 #if defined HAVE_NETCDF
1485 !scalars
1486  integer :: ncerr
1487  integer :: natom_id,nimage_id,ntypat_id,npsp_id,time_id,xyz_id,six_id
1488  integer :: xcart_id,xred_id,fcart_id,gred_id,vel_id,vel_cell_id
1489  integer :: rprimd_id,acell_id,strten_id
1490  integer :: etotal_id,ekin_id,entropy_id,mdtime_id
1491  integer :: typat_id,znucl_id,amu_id,dtion_id,imgmov_id, two_id,mdtemp_id
1492  !character(len=500) :: msg
1493 !arrays
1494  integer :: dim0(0),dim1(1),dim2(2),dim3(3),dim4(4)
1495 #endif
1496 
1497 ! *************************************************************************
1498 
1499 #if defined HAVE_NETCDF
1500 
1501 !1.Define the dimensions
1502 
1503  if (npsp/=ntypat) then
1504    ABI_WARNING('HIST file does not support alchemical mixing!')
1505  end if
1506 
1507  ncerr = nf90_def_dim(ncid,"natom",natom,natom_id)
1508  NCF_CHECK_MSG(ncerr," define dimension natom")
1509 
1510  ncerr = nf90_def_dim(ncid,"ntypat",ntypat,ntypat_id)
1511  NCF_CHECK_MSG(ncerr," define dimension ntypat")
1512 
1513  if (has_nimage) then
1514    ncerr = nf90_def_dim(ncid,"nimage",nimage,nimage_id)
1515    NCF_CHECK_MSG(ncerr," define dimension nimage")
1516  end if
1517 
1518  ncerr = nf90_def_dim(ncid,"npsp",npsp,npsp_id)
1519  NCF_CHECK_MSG(ncerr," define dimension npsp")
1520 
1521  ncerr = nf90_def_dim(ncid,"xyz",3,xyz_id)
1522  NCF_CHECK_MSG(ncerr," define dimension xyz")
1523 
1524  ncerr = nf90_def_dim(ncid,"six",6,six_id)
1525  NCF_CHECK_MSG(ncerr," define dimension six")
1526 
1527  ncerr = nf90_def_dim(ncid,"time",NF90_UNLIMITED,time_id)
1528  NCF_CHECK_MSG(ncerr," define dimension time")
1529 
1530  ncerr = nf90_def_dim(ncid,"two",2,two_id)
1531  NCF_CHECK_MSG(ncerr," define dimension two")
1532 
1533 !2.Define the constant variables
1534 
1535  dim1=(/natom_id/)
1536  call ab_define_var(ncid,dim1,typat_id,NF90_DOUBLE,&
1537 &  "typat","types of atoms","dimensionless" )
1538 
1539  dim1=(/npsp_id/)
1540  call ab_define_var(ncid,dim1,znucl_id,NF90_DOUBLE,&
1541 &  "znucl","atomic charges","atomic units" )
1542 
1543  dim1=(/ntypat_id/)
1544  call ab_define_var(ncid,dim1,amu_id,NF90_DOUBLE,&
1545 &  "amu","atomic masses","atomic units" )
1546 
1547  call ab_define_var(ncid,dim0,dtion_id,NF90_DOUBLE,&
1548 &  "dtion","time step","atomic units" )
1549 
1550 !mdtemp
1551  dim1=(/two_id/)
1552  call ab_define_var(ncid,dim1,mdtemp_id,NF90_DOUBLE,&
1553 &  "mdtemp","Molecular Dynamics Thermostat Temperatures","Kelvin" )
1554 
1555 !mdtime
1556  dim1=(/time_id/)
1557  call ab_define_var(ncid,dim1,mdtime_id,NF90_DOUBLE,&
1558 & "mdtime","Molecular Dynamics or Relaxation TIME","hbar/Ha" )
1559 
1560 !3.Define the evolving variables
1561 
1562 !xcart,xred,fcart,gred,vel
1563  if (has_nimage) then
1564    call ab_define_var(ncid,dim0,imgmov_id,NF90_INT,&
1565   &  "imgmov","Image mover","Not relevant" )
1566    dim4=(/xyz_id,natom_id,nimage_id,time_id/)
1567    call ab_define_var(ncid,dim4,xcart_id,NF90_DOUBLE,&
1568 &   "xcart","vectors (X) of atom positions in CARTesian coordinates","bohr" )
1569    call ab_define_var(ncid,dim4,xred_id,NF90_DOUBLE,&
1570 &   "xred","vectors (X) of atom positions in REDuced coordinates","dimensionless" )
1571    call ab_define_var(ncid,dim4,fcart_id,NF90_DOUBLE,&
1572 &   "fcart","atom Forces in CARTesian coordinates","Ha/bohr" )
1573    call ab_define_var(ncid,dim4,gred_id,NF90_DOUBLE,&
1574 &   "fred","atom Forces in REDuced coordinates","dimensionless" ) ! XG210315 : should be "Gradients in REDuced ..."
1575    call ab_define_var(ncid,dim4,vel_id,NF90_DOUBLE,&
1576 &   "vel","VELocities of atoms","bohr*Ha/hbar" )
1577  else
1578    dim3=(/xyz_id,natom_id,time_id/)
1579    call ab_define_var(ncid,dim3,xcart_id,NF90_DOUBLE,&
1580 &   "xcart","vectors (X) of atom positions in CARTesian coordinates","bohr" )
1581    call ab_define_var(ncid,dim3,xred_id,NF90_DOUBLE,&
1582 &   "xred","vectors (X) of atom positions in REDuced coordinates","dimensionless" )
1583    call ab_define_var(ncid,dim3,fcart_id,NF90_DOUBLE,&
1584 &   "fcart","atom Forces in CARTesian coordinates","Ha/bohr" )
1585    call ab_define_var(ncid,dim3,gred_id,NF90_DOUBLE,&
1586 &   "fred","atom Forces in REDuced coordinates","dimensionless" ) ! XG210315 : should be "Gradients in REDuced ..."
1587    call ab_define_var(ncid,dim3,vel_id,NF90_DOUBLE,&
1588 &   "vel","VELocities of atoms","bohr*Ha/hbar" )
1589  end if
1590 
1591 !rprimd,vel_cell
1592  if (has_nimage) then
1593    dim4=(/xyz_id,xyz_id,nimage_id,time_id/)
1594    call ab_define_var(ncid,dim4,rprimd_id,NF90_DOUBLE,&
1595 &   "rprimd","Real space PRIMitive translations, Dimensional","bohr" )
1596    call ab_define_var(ncid,dim4,vel_cell_id,NF90_DOUBLE,&
1597 &   "vel_cell","VELocities of CELl","bohr*Ha/hbar" )
1598  else
1599    dim3=(/xyz_id,xyz_id,time_id/)
1600    call ab_define_var(ncid,dim3,rprimd_id,NF90_DOUBLE,&
1601 &   "rprimd","Real space PRIMitive translations, Dimensional","bohr" )
1602    call ab_define_var(ncid,dim3,vel_cell_id,NF90_DOUBLE,&
1603 &   "vel_cell","VELocities of cell","bohr*Ha/hbar" )
1604  end if
1605 
1606 !acell
1607  if (has_nimage) then
1608    dim3=(/xyz_id,nimage_id,time_id/)
1609    call ab_define_var(ncid,dim3,acell_id,NF90_DOUBLE,&
1610 &   "acell","CELL lattice vector scaling","bohr" )
1611  else
1612    dim2=(/xyz_id,time_id/)
1613    call ab_define_var(ncid,dim2,acell_id,NF90_DOUBLE,&
1614 &   "acell","CELL lattice vector scaling","bohr" )
1615  end if
1616 
1617 !strten
1618  if (has_nimage) then
1619    dim3=(/six_id,nimage_id,time_id/)
1620    call ab_define_var(ncid,dim3,strten_id,NF90_DOUBLE,&
1621 &   "strten","STRess tensor","Ha/bohr^3" )
1622  else
1623    dim2=(/six_id,time_id/)
1624    call ab_define_var(ncid,dim2,strten_id,NF90_DOUBLE,&
1625 &   "strten","STRess tensor","Ha/bohr^3" )
1626  end if
1627 
1628 !etotal,ekin,entropy
1629  if (has_nimage) then
1630    dim2=(/nimage_id,time_id/)
1631    call ab_define_var(ncid,dim2,etotal_id,NF90_DOUBLE,&
1632 &   "etotal","TOTAL Energy","Ha" )
1633    call ab_define_var(ncid,dim2,ekin_id,NF90_DOUBLE,&
1634 &   "ekin","Energy KINetic ionic","Ha" )
1635    call ab_define_var(ncid,dim2,entropy_id,NF90_DOUBLE,&
1636 &   "entropy","Entropy","" )
1637  else
1638    dim1=(/time_id/)
1639    call ab_define_var(ncid,dim1,etotal_id,NF90_DOUBLE,&
1640 &   "etotal","TOTAL Energy","Ha" )
1641    call ab_define_var(ncid,dim1,ekin_id,NF90_DOUBLE,&
1642 &   "ekin","Energy KINetic ionic","Ha" )
1643    call ab_define_var(ncid,dim1,entropy_id,NF90_DOUBLE,&
1644 &   "entropy","Entropy","" )
1645  end if
1646 
1647 !4.End define mode
1648 
1649  ncerr = nf90_enddef(ncid)
1650  NCF_CHECK_MSG(ncerr," end define mode")
1651 
1652 #endif
1653 
1654 end subroutine def_file_hist

m_abihist/get_dims_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 get_dims_hist

FUNCTION

INPUTS

OUTPUT

SOURCE

1672 subroutine get_dims_hist(ncid,natom,ntypat,nimage,time,&
1673 &          natom_id,ntypat_id,nimage_id,time_id,xyz_id,six_id,has_nimage)
1674 
1675 !Arguments ------------------------------------
1676 !scalars
1677  integer,intent(in) :: ncid
1678  integer,intent(out) :: natom,nimage,time,ntypat
1679  integer,intent(out) :: natom_id,nimage_id,time_id,xyz_id,six_id, ntypat_id
1680  logical,intent(out) :: has_nimage
1681 
1682 !Local variables-------------------------------
1683 #if defined HAVE_NETCDF
1684 !scalars
1685  integer :: ncerr
1686  character(len=5) :: char_tmp
1687 #endif
1688 
1689 ! *************************************************************************
1690 
1691 #if defined HAVE_NETCDF
1692 !Inquire dimensions IDs
1693 
1694  ncerr = nf90_inq_dimid(ncid,"natom",natom_id)
1695  NCF_CHECK_MSG(ncerr," inquire dimension ID for natom")
1696 
1697  ncerr = nf90_inq_dimid(ncid,"npsp",ntypat_id)
1698  NCF_CHECK_MSG(ncerr," inquire dimension ID for npsp")
1699 
1700  ncerr = nf90_inq_dimid(ncid,"xyz",xyz_id)
1701  NCF_CHECK_MSG(ncerr," inquire dimension ID for xyz")
1702 
1703  ncerr = nf90_inq_dimid(ncid,"time",time_id)
1704  NCF_CHECK_MSG(ncerr," inquire dimension ID for time")
1705 
1706  ncerr = nf90_inq_dimid(ncid,"six",six_id)
1707  NCF_CHECK_MSG(ncerr," inquire dimension ID for six")
1708 
1709  ncerr = nf90_inq_dimid(ncid,"nimage",nimage_id)
1710  has_nimage=(ncerr==nf90_noerr)
1711 
1712 !Inquire dimensions lengths
1713 
1714  if (has_nimage) then
1715    ncerr = nf90_inquire_dimension(ncid,nimage_id,char_tmp,nimage)
1716    has_nimage=(ncerr==nf90_noerr)
1717  end if
1718  if (.not.has_nimage) nimage=1
1719 
1720  ncerr = nf90_inquire_dimension(ncid,natom_id,char_tmp,natom)
1721  NCF_CHECK_MSG(ncerr," inquire dimension natom")
1722 
1723  ncerr = nf90_inquire_dimension(ncid,ntypat_id,char_tmp,ntypat)
1724  NCF_CHECK_MSG(ncerr," inquire dimension ntypat")
1725 
1726  ncerr = nf90_inquire_dimension(ncid,time_id,char_tmp,time)
1727  NCF_CHECK_MSG(ncerr," inquire dimension time")
1728 
1729 #endif
1730 
1731 end subroutine get_dims_hist

m_abihist/get_varid_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 get_varid_hist

FUNCTION

INPUTS

OUTPUT

SOURCE

1749 subroutine get_varid_hist(ncid,xcart_id,xred_id,fcart_id,gred_id,vel_id,vel_cell_id,&
1750 &          rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id,has_nimage)
1751 
1752 !Arguments ------------------------------------
1753 !scalars
1754  integer,intent(in) :: ncid
1755  integer,intent(out) :: xcart_id,xred_id,fcart_id,gred_id,vel_id
1756  integer,intent(out) :: vel_cell_id,rprimd_id,acell_id,strten_id
1757  integer,intent(out) :: etotal_id,ekin_id,entropy_id,mdtime_id
1758  logical,intent(in)  :: has_nimage
1759 !Local variables-------------------------------
1760 #if defined HAVE_NETCDF
1761 !scalars
1762  integer :: ncerr
1763 #endif
1764 
1765 ! *************************************************************************
1766 
1767 #if defined HAVE_NETCDF
1768 
1769  ncerr = nf90_inq_varid(ncid, "mdtime", mdtime_id)
1770  NCF_CHECK_MSG(ncerr," get the id for mdtime")
1771 
1772  ncerr = nf90_inq_varid(ncid, "xcart", xcart_id)
1773  NCF_CHECK_MSG(ncerr," get the id for xcart")
1774 
1775  ncerr = nf90_inq_varid(ncid, "xred", xred_id)
1776  NCF_CHECK_MSG(ncerr," get the id for xred")
1777 
1778  ncerr = nf90_inq_varid(ncid, "fcart", fcart_id)
1779  NCF_CHECK_MSG(ncerr," get the id for fcart")
1780 
1781  ncerr = nf90_inq_varid(ncid, "fred", gred_id)
1782  NCF_CHECK_MSG(ncerr," get the id for fred") ! XG210315 : should be "gred"
1783 
1784  ncerr = nf90_inq_varid(ncid, "vel", vel_id)
1785  NCF_CHECK_MSG(ncerr," get the id for vel")
1786 
1787  ncerr = nf90_inq_varid(ncid, "vel_cell", vel_cell_id)
1788  if(has_nimage) then
1789    NCF_CHECK_MSG(ncerr," get the id for vel_cell")
1790  end if
1791 
1792  ncerr = nf90_inq_varid(ncid, "rprimd", rprimd_id)
1793  NCF_CHECK_MSG(ncerr," get the id for rprimd")
1794 
1795  ncerr = nf90_inq_varid(ncid, "acell", acell_id)
1796  NCF_CHECK_MSG(ncerr," get the id for acell")
1797 
1798  ncerr = nf90_inq_varid(ncid, "strten", strten_id)
1799  NCF_CHECK_MSG(ncerr," get the id for strten")
1800 
1801  ncerr = nf90_inq_varid(ncid, "etotal", etotal_id)
1802  NCF_CHECK_MSG(ncerr," get the id for etotal")
1803 
1804  ncerr = nf90_inq_varid(ncid, "ekin", ekin_id)
1805  NCF_CHECK_MSG(ncerr," get the id for ekin")
1806 
1807  ncerr = nf90_inq_varid(ncid, "entropy", entropy_id)
1808  NCF_CHECK_MSG(ncerr," get the id for entropy")
1809 
1810 #endif
1811 
1812 end subroutine get_varid_hist

m_abihist/hist2var [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 hist2var

FUNCTION

 Return the last values of xred, acell and rprimd stored in the history.

INPUTS

 natom = Number of atoms
 hist<type abihist>=Historical record of positions, forces, acell, stresses, and energies,
 zDebug = If true some output will be printed

OUTPUT

  xred(3,natom) = reduced dimensionless atomic coordinates
  acell(3)    = length scales of primitive translations (bohr)
  rprimd(3,3) = dimensional real space primitive translations (bohr)

SOURCE

702 subroutine hist2var(acell,hist,natom,rprimd,xred,zDEBUG)
703 
704 !Arguments ------------------------------------
705 !scalars
706 integer,intent(in) :: natom
707 type(abihist),intent(in) :: hist
708 logical,intent(in) :: zDEBUG
709 !arrays
710 real(dp),intent(out) :: acell(3)
711 real(dp),intent(out) :: rprimd(3,3)
712 real(dp),intent(out) :: xred(3,natom)
713 
714 !Local variables-------------------------------
715 !scalars
716 integer :: kk
717 
718 ! *************************************************************
719 
720  xred  (:,:)=hist%xred(:,:,hist%ihist)
721  acell (:  )=hist%acell(:,hist%ihist)
722  rprimd(:,:)=hist%rprimd(:,:,hist%ihist)
723 
724  if(zDEBUG)then
725    write (std_out,*) 'Atom positions and cell parameters '
726    write (std_out,*) 'ihist: ',hist%ihist
727    write (std_out,*) 'xred:'
728    do kk=1,natom
729      write (std_out,*) xred(:,kk)
730    end do
731    write(std_out,*) 'rprimd:'
732    do kk=1,3
733      write(std_out,*) rprimd(:,kk)
734    end do
735    write(std_out,*) 'acell:'
736    write(std_out,*) acell(:)
737  end if
738 
739 end subroutine hist2var

m_abihist/read_csts_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 read_csts_hist

FUNCTION

INPUTS

OUTPUT

SOURCE

1830 subroutine read_csts_hist(ncid,dtion,typat,znucl,amu)
1831 
1832 !Arguments ------------------------------------
1833 !scalars
1834  integer,intent(in) :: ncid
1835  real(dp),intent(out) :: dtion
1836 !arrays
1837  integer,intent(out) :: typat(:)
1838  real(dp),intent(out) :: amu(:),znucl(:)
1839 
1840 !Local variables-------------------------------
1841 #if defined HAVE_NETCDF
1842 !scalars
1843  integer :: ncerr
1844  integer :: typat_id,znucl_id,amu_id,dtion_id
1845 #endif
1846 
1847 ! *************************************************************************
1848 
1849 #if defined HAVE_NETCDF
1850 
1851 !1.Get the IDs
1852  ncerr = nf90_inq_varid(ncid, "typat", typat_id)
1853  NCF_CHECK_MSG(ncerr," get the id for typat")
1854 
1855  ncerr = nf90_inq_varid(ncid, "znucl", znucl_id)
1856  NCF_CHECK_MSG(ncerr," get the id for znucl")
1857 
1858  ncerr = nf90_inq_varid(ncid, "amu", amu_id)
1859  NCF_CHECK_MSG(ncerr," get the id for amu")
1860 
1861  ncerr = nf90_inq_varid(ncid, "dtion", dtion_id)
1862  NCF_CHECK_MSG(ncerr," get the id for dtion")
1863 
1864 !2.Write the constants
1865  ncerr = nf90_get_var(ncid, typat_id, typat)
1866  NCF_CHECK_MSG(ncerr," get variable typat")
1867 
1868  ncerr = nf90_get_var(ncid, znucl_id, znucl)
1869  NCF_CHECK_MSG(ncerr," get variable znucl")
1870 
1871  ncerr = nf90_get_var(ncid, amu_id, amu)
1872  NCF_CHECK_MSG(ncerr," get variable amu")
1873 
1874  ncerr = nf90_get_var(ncid, dtion_id, dtion)
1875  NCF_CHECK_MSG(ncerr," get variable dtion")
1876 
1877 #endif
1878 
1879 end subroutine read_csts_hist

m_abihist/read_md_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 read_md_hist

FUNCTION

 Read the history file from a netcdf file and store it into a hist dataset structure
 This version is not compatible with multiple images of the simulation cell.

INPUTS

  filename = Filename of the NetCDF to read
  isVUsed,isARUsed=flags used to initialize hist structure

OUTPUT

  hist<type abihist>=Historical record of positions, forces, stresses, cell dims and energies,

SOURCE

1276 subroutine read_md_hist(filename,hist,isVUsed,isARUsed,readOnlyLast)
1277 
1278 !Arguments ------------------------------------
1279 !scalars
1280  logical,intent(in) :: isVUsed,isARUsed,readOnlyLast
1281  character(len=*),intent(in) :: filename
1282 !arrays
1283  type(abihist),intent(inout),target :: hist
1284 
1285 !Local variables-------------------------------
1286 #if defined HAVE_NETCDF
1287 !scalars
1288  integer :: ncerr,ncid,nimage,natom,time,start_time, ntypat
1289  integer :: nimage_id,natom_id,xyz_id,time_id,six_id, ntypat_id
1290  integer :: xcart_id,xred_id,fcart_id,gred_id,ekin_id,entropy_id
1291  integer :: mdtime_id,vel_id,vel_cell_id,etotal_id
1292  integer :: acell_id,rprimd_id,strten_id
1293  logical :: has_nimage
1294 #endif
1295 
1296 ! *************************************************************************
1297 
1298 #if defined HAVE_NETCDF
1299 
1300  hist%ihist=0 ; hist%mxhist=0
1301 
1302 !Open netCDF file
1303  ncerr=nf90_open(path=trim(filename),mode=NF90_NOWRITE,ncid=ncid)
1304  if(ncerr /= NF90_NOERR) then
1305    write(std_out,*) 'Could no open ',trim(filename),', starting from scratch'
1306    return
1307  else
1308    write(std_out,*) 'Succesfully open ',trim(filename),' for reading'
1309    write(std_out,*) 'Extracting information from NetCDF file...'
1310  end if
1311 
1312 !Inquire dimensions IDs and lengths
1313  call get_dims_hist(ncid,natom,ntypat,nimage,time,&
1314 &     natom_id,ntypat_id,nimage_id,time_id,xyz_id,six_id,has_nimage)
1315 
1316 !If only the last step is needing (restarxf==-3 for example)
1317  if(readOnlyLast)then
1318    start_time = time
1319    time = 1
1320  else
1321    start_time = 1
1322  end if
1323 
1324 !Allocate hist structure
1325  call abihist_init(hist,natom,time,isVused,isARused)
1326 
1327 !Get the ID of a variables from their name
1328  call get_varid_hist(ncid,xcart_id,xred_id,fcart_id,gred_id,vel_id,vel_cell_id,&
1329 &     rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id,has_nimage)
1330 
1331 !Read variables from the dataset and write them into hist
1332  call read_vars_hist(ncid,hist,natom,time,has_nimage,1,start_time,&
1333 &     xred_id,fcart_id,vel_id,vel_cell_id,rprimd_id,acell_id,&
1334 &     strten_id,etotal_id,ekin_id,entropy_id,mdtime_id)
1335 
1336 !Close NetCDF file
1337  ncerr = nf90_close(ncid)
1338  NCF_CHECK_MSG(ncerr," close netcdf history file")
1339 
1340 #endif
1341 
1342 end subroutine read_md_hist

m_abihist/read_md_hist_img [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 read_md_hist_img

FUNCTION

 Read the history file from a netcdf file and store it into a hist dataset structure
 This version is compatible with multiple images of the simulation cell.

INPUTS

  filename = Filename of the NetCDF to read
  isVUsed,isARUsed=flags used to initialize hist structure
  [imgtab(:)]= In case of multiple images,indexes of images to be read
               Default is 1,2,3,...
               Size must be equal to size(hist)

OUTPUT

  hist(:)<type abihist>=Historical record of positions, forces, stresses, cell dims and energies,
    Size(hist) is equal to a number of images to be read

SOURCE

1369 subroutine read_md_hist_img(filename,hist,isVUsed,isARused,imgtab)
1370 
1371 !Arguments ------------------------------------
1372 !scalars
1373  logical,intent(in) :: isVUsed,isARused
1374  character(len=*),intent(in) :: filename
1375 !arrays
1376  integer,intent(in),optional :: imgtab(:)
1377  type(abihist),intent(inout),target :: hist(:)
1378 
1379 !Local variables-------------------------------
1380 #if defined HAVE_NETCDF
1381 !scalars
1382  integer :: iimage,iimg,my_nimage,ncerr,ncid,nimage,natom,time
1383  integer :: nimage_id,natom_id,xyz_id,time_id,six_id, ntypat, ntypat_id
1384  integer :: xcart_id,xred_id,fcart_id,gred_id,ekin_id,entropy_id
1385  integer :: mdtime_id,vel_id,vel_cell_id,etotal_id
1386  integer :: acell_id,rprimd_id,strten_id
1387  logical :: has_nimage
1388  !character(len=500) :: msg
1389  type(abihist),pointer :: hist_
1390  integer,allocatable :: my_imgtab(:)
1391 #endif
1392 
1393 ! *************************************************************************
1394 
1395 #if defined HAVE_NETCDF
1396 
1397  hist%ihist=0 ; hist%mxhist=0
1398 
1399 !Open netCDF file
1400  ncerr=nf90_open(path=trim(filename),mode=NF90_NOWRITE,ncid=ncid)
1401  if(ncerr /= NF90_NOERR) then
1402    write(std_out,*) 'Could no open ',trim(filename),', starting from scratch'
1403    return
1404  else
1405    write(std_out,*) 'Succesfully open ',trim(filename),' for reading'
1406    write(std_out,*) 'Extracting information from NetCDF file...'
1407  end if
1408 
1409  !Manage multiple images of the cell
1410  my_nimage=size(hist)
1411  if (my_nimage==0) return
1412  ABI_MALLOC(my_imgtab,(my_nimage))
1413  if (present(imgtab)) then
1414   if (size(my_imgtab)/=my_nimage) then
1415     ABI_BUG('Inconsistency between hist and imgtab!')
1416   end if
1417   my_imgtab(:)=imgtab(:)
1418  else
1419    my_imgtab(:)=(/(iimage,iimage=1,my_nimage)/)
1420  end if
1421 
1422 !Inquire dimensions IDs and lengths
1423  call get_dims_hist(ncid,natom,ntypat,nimage,time,&
1424 &     natom_id,ntypat_id,nimage_id,time_id,xyz_id,six_id,has_nimage)
1425 
1426  if (nimage<maxval(my_imgtab)) then
1427    ABI_ERROR('Not enough images in the HIST file!')
1428  end if
1429 
1430 !Loop over images
1431  do iimage=1,my_nimage
1432    iimg=my_imgtab(iimage)
1433    hist_ => hist(iimage)
1434 
1435 !  Allocate hist structure
1436    call abihist_init(hist_,natom,time,isVused,isARused)
1437 
1438 !  Get the ID of a variables from their name
1439    call get_varid_hist(ncid,xcart_id,xred_id,fcart_id,gred_id,vel_id,vel_cell_id,&
1440 &       rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id,has_nimage)
1441 
1442 !  Read variables from the dataset and write them into hist
1443    call read_vars_hist(ncid,hist_,natom,time,has_nimage,iimg,1,&
1444 &       xred_id,fcart_id,vel_id,vel_cell_id,rprimd_id,acell_id,&
1445 &       strten_id,etotal_id,ekin_id,entropy_id,mdtime_id)
1446 
1447  end do
1448 
1449 !Close NetCDF file
1450  ncerr = nf90_close(ncid)
1451  NCF_CHECK_MSG(ncerr," close netcdf history file")
1452 
1453  ABI_FREE(my_imgtab)
1454 
1455 #endif
1456 
1457 end subroutine read_md_hist_img

m_abihist/read_vars_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 read_vars_hist

FUNCTION

INPUTS

OUTPUT

SOURCE

2136 subroutine read_vars_hist(ncid,hist,natom,time,has_nimage,iimg,start_time,&
2137 &          xred_id,fcart_id,vel_id,vel_cell_id,rprimd_id,acell_id,&
2138 &          strten_id,etotal_id,ekin_id,entropy_id,mdtime_id)
2139 
2140 !Arguments ------------------------------------
2141 !scalars
2142  integer,intent(in) :: ncid,natom,time,iimg
2143  integer,intent(in) :: xred_id,fcart_id,vel_id,vel_cell_id,rprimd_id
2144  integer,intent(in) :: acell_id,strten_id
2145  integer,intent(in) :: etotal_id,ekin_id,entropy_id,mdtime_id
2146  integer,intent(in) :: start_time
2147  logical,intent(in) :: has_nimage
2148  type(abihist),intent(inout),target :: hist
2149 
2150 !Local variables-------------------------------
2151 #if defined HAVE_NETCDF
2152 !scalars
2153  integer :: ncerr
2154 !arrays
2155  integer :: count1(1),count2(2),count3(3),count4(4)
2156  integer :: start1(1),start2(2),start3(3),start4(4)
2157 #endif
2158 
2159 ! *************************************************************************
2160 
2161 #if defined HAVE_NETCDF
2162 
2163 !Variables not depending on imes
2164 
2165 !mdtime
2166  start1=(/start_time/);count1=(/time/)
2167  ncerr = nf90_get_var(ncid,mdtime_id,hist%time(:),count=count1,start=start1)
2168  NCF_CHECK_MSG(ncerr," read variable mdtime")
2169 
2170 !Variables depending on images
2171 
2172 !xred,fcart,vel
2173  if (has_nimage) then
2174    start4=(/1,1,iimg,start_time/);count4=(/3,natom,1,time/)
2175    ncerr = nf90_get_var(ncid,xred_id  ,hist%xred(:,:,:),count=count4,start=start4)
2176    NCF_CHECK_MSG(ncerr," read variable xred")
2177    ncerr = nf90_get_var(ncid,fcart_id ,hist%fcart(:,:,:),count=count4,start=start4)
2178    NCF_CHECK_MSG(ncerr," read variable fcart")
2179    ncerr = nf90_get_var(ncid,vel_id,hist%vel(:,:,:),count=count4,start=start4)
2180    NCF_CHECK_MSG(ncerr," read variable vel")
2181  else
2182    start3=(/1,1,start_time/);count3=(/3,natom,time/)
2183    ncerr = nf90_get_var(ncid,xred_id  ,hist%xred(:,:,:),count=count3,start=start3)
2184    NCF_CHECK_MSG(ncerr," read variable xred")
2185    ncerr = nf90_get_var(ncid,fcart_id ,hist%fcart(:,:,:),count=count3,start=start3)
2186    NCF_CHECK_MSG(ncerr," read variable fcart")
2187    ncerr = nf90_get_var(ncid,vel_id,hist%vel(:,:,:),count=count3,start=start3)
2188    NCF_CHECK_MSG(ncerr," read variable vel")
2189  end if
2190 
2191 !rprimd,vel_cell
2192  if (has_nimage) then
2193    start4=(/1,1,iimg,start_time/);count4=(/3,3,start_time,time/)
2194    ncerr = nf90_get_var(ncid,rprimd_id,hist%rprimd(:,:,:),count=count4,start=start4)
2195    NCF_CHECK_MSG(ncerr," read variable rprimd")
2196    ncerr = nf90_get_var(ncid,vel_cell_id,hist%vel_cell(:,:,:),count=count4,start=start4)
2197    NCF_CHECK_MSG(ncerr," read variable vel_cell")
2198  else
2199    start3=(/1,1,start_time/);count3=(/3,3,time/)
2200    ncerr = nf90_get_var(ncid,rprimd_id,hist%rprimd(:,:,:),count=count3,start=start3)
2201    NCF_CHECK_MSG(ncerr," read variable rprimd")
2202  end if
2203 
2204 !acell
2205  if (has_nimage) then
2206    start3=(/1,iimg,start_time/);count3=(/3,1,time/)
2207    ncerr = nf90_get_var(ncid,acell_id,hist%acell(:,:),count=count3,start=start3)
2208    NCF_CHECK_MSG(ncerr," read variable acell")
2209  else
2210    start2=(/1,start_time/);count2=(/3,time/)
2211    ncerr = nf90_get_var(ncid,acell_id,hist%acell(:,:),count=count2,start=start2)
2212    NCF_CHECK_MSG(ncerr," read variable acell")
2213  end if
2214 
2215 !strten
2216  if (has_nimage) then
2217    start3=(/1,iimg,start_time/);count3=(/6,1,time/)
2218    ncerr = nf90_get_var(ncid, strten_id,hist%strten(:,:),count=count3,start=start3)
2219    NCF_CHECK_MSG(ncerr," read variable strten")
2220  else
2221    start2=(/1,start_time/);count2=(/6,time/)
2222    ncerr = nf90_get_var(ncid, strten_id,hist%strten(:,:),count=count2,start=start2)
2223    NCF_CHECK_MSG(ncerr," read variable strten")
2224  end if
2225 
2226 !etotal,ekin,entropy
2227  if (has_nimage) then
2228    start2=(/1,start_time/);count2=(/1,time/)
2229    ncerr = nf90_get_var(ncid,etotal_id,hist%etot(:),count=count2,start=start2)
2230    NCF_CHECK_MSG(ncerr," read variable etotal")
2231    ncerr = nf90_get_var(ncid,ekin_id ,hist%ekin(:),count=count2,start=start2)
2232    NCF_CHECK_MSG(ncerr," read variable ekin")
2233    ncerr = nf90_get_var(ncid,entropy_id,hist%entropy(:),count=count2,start=start2)
2234    NCF_CHECK_MSG(ncerr," read variable entropy")
2235  else
2236    start1=(/start_time/);count1=(/time/)
2237    ncerr = nf90_get_var(ncid,etotal_id,hist%etot(:),count=count1,start=start1)
2238    NCF_CHECK_MSG(ncerr," read variable etotal")
2239    ncerr = nf90_get_var(ncid,ekin_id,hist%ekin(:),count=count1,start=start1)
2240    NCF_CHECK_MSG(ncerr," read variable ekin")
2241    ncerr = nf90_get_var(ncid,entropy_id,hist%entropy(:),count=count1,start=start1)
2242    NCF_CHECK_MSG(ncerr," read variable entropy")
2243  end if
2244 
2245 #endif
2246 
2247 end subroutine read_vars_hist

m_abihist/var2hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 var2hist

FUNCTION

 Set the values of the history "hist" with the values of xred, acell and rprimd

INPUTS

 natom = number of atoms
 xred(3,natom) = reduced dimensionless atomic coordinates
 acell(3)    = length scales of primitive translations (bohr)
 rprimd(3,3) = dimensionlal real space primitive translations (bohr)

OUTPUT

SIDE EFFECTS

 hist<type abihist>=Historical record of positions, forces, acell, stresses, and energies,

SOURCE

587 subroutine var2hist(acell,hist,natom,rprimd,xred,zDEBUG)
588 
589 !Arguments ------------------------------------
590 !scalars
591  integer,intent(in) :: natom
592  type(abihist),intent(inout) :: hist
593  logical,intent(in) :: zDEBUG
594 !arrays
595  real(dp),intent(in) :: acell(3)
596  real(dp),intent(in) :: rprimd(3,3)
597  real(dp),intent(in) :: xred(3,natom)
598 
599 !Local variables-------------------------------
600 !scalars
601  integer :: kk
602 
603 ! *************************************************************
604 
605  hist%xred(:,:,hist%ihist)=xred(:,:)
606  hist%rprimd(:,:,hist%ihist)=rprimd(:,:)
607  hist%acell(:,hist%ihist)=acell(:)
608 
609  if(zDEBUG)then
610    write (std_out,*) 'Atom positions and cell parameters '
611    write (std_out,*) 'ihist: ',hist%ihist
612    write (std_out,*) 'xred:'
613    do kk=1,natom
614      write (std_out,*) xred(:,kk)
615    end do
616    write(std_out,*) 'rprimd:'
617    do kk=1,3
618      write(std_out,*) rprimd(:,kk)
619    end do
620    write(std_out,*) 'acell:'
621    write(std_out,*) acell(:)
622  end if
623 
624 end subroutine var2hist

m_abihist/vel2hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 vel2hist

FUNCTION

 Set the values of the history "hist" related with velocities, ie
 The array of velocities and Kinetic Energy

INPUTS

 amass(natom) = mass of the atoms
 vel(3,natom)= Velocities of the atoms
 vel_cell(3,3)= Velocities of the cell

OUTPUT

SIDE EFFECTS

 hist<type abihist>=Historical record of positions, forces, stresses, cell and energies,

SOURCE

765 subroutine vel2hist(amass,hist,vel,vel_cell)
766 
767 !Arguments ------------------------------------
768 !scalars
769 type(abihist),intent(inout) :: hist
770 !arrays
771 real(dp),intent(in) :: amass(:)
772 real(dp),intent(in) :: vel(:,:)
773 real(dp),intent(in) :: vel_cell(:,:)
774 
775 !Local variables-------------------------------
776 !scalars
777 integer :: ii,jj,natom
778 real(dp) :: ekin
779 
780 ! *************************************************************
781 
782  natom=size(vel,2)
783 
784  if (hist%isVused) then
785 
786 !  Store the velocities
787    hist%vel(:,:,hist%ihist)=vel(:,:)
788    hist%vel_cell(:,:,hist%ihist)=vel_cell(:,:)
789 
790 !  Compute the Ionic Kinetic energy
791    ekin=zero
792    do ii=1,natom
793      do jj=1,3
794        ekin=ekin+half*amass(ii)*vel(jj,ii)**2
795      end do
796    end do
797 
798  else
799    hist%vel(:,:,hist%ihist)=zero
800    hist%vel_cell(:,:,hist%ihist)=zero
801    ekin=zero
802  end if
803 
804 !Store the Ionic Kinetic Energy
805  hist%ekin(hist%ihist)=ekin
806 
807 end subroutine vel2hist

m_abihist/write_csts_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 write_csts_hist

FUNCTION

INPUTS

OUTPUT

SOURCE

1897 subroutine write_csts_hist(ncid,dtion,imgmov,typat,znucl,amu,mdtemp)
1898 
1899 !Arguments ------------------------------------
1900 !scalars
1901  integer,intent(in) :: ncid
1902  real(dp),intent(in) :: dtion
1903  integer,intent(in) :: imgmov
1904 !arrays
1905  integer,intent(in) :: typat(:)
1906  real(dp),intent(in) :: amu(:),znucl(:), mdtemp(2)
1907 
1908 !Local variables-------------------------------
1909 #if defined HAVE_NETCDF
1910 !scalars
1911  integer :: ncerr
1912  integer :: typat_id,znucl_id,amu_id,dtion_id, imgmov_id, mdtemp_id
1913 #endif
1914 
1915 ! *************************************************************************
1916 
1917 #if defined HAVE_NETCDF
1918 
1919 !1.Get the IDs
1920 
1921  ncerr = nf90_inq_varid(ncid, "typat", typat_id)
1922  NCF_CHECK_MSG(ncerr," get the id for typat")
1923 
1924  ncerr = nf90_inq_varid(ncid, "znucl", znucl_id)
1925  NCF_CHECK_MSG(ncerr," get the id for znucl")
1926 
1927  ncerr = nf90_inq_varid(ncid, "amu", amu_id)
1928  NCF_CHECK_MSG(ncerr," get the id for amu")
1929 
1930  ncerr = nf90_inq_varid(ncid, "dtion", dtion_id)
1931  NCF_CHECK_MSG(ncerr," get the id for dtion")
1932 
1933  if ( nf90_noerr == nf90_inq_varid(ncid, "imgmov", imgmov_id) ) then
1934    ncerr = nf90_put_var(ncid, imgmov_id, imgmov)
1935    NCF_CHECK_MSG(ncerr," write variable imgmov")
1936  end if
1937 
1938  if ( nf90_noerr == nf90_inq_varid(ncid, "mdtemp", mdtemp_id) ) then
1939    ncerr = nf90_put_var(ncid, mdtemp_id, mdtemp)
1940    NCF_CHECK_MSG(ncerr," write variable mdtemp")
1941  end if
1942 
1943 !2.Write the constants
1944 
1945  ncerr = nf90_put_var(ncid, typat_id, typat)
1946  NCF_CHECK_MSG(ncerr," write variable typat")
1947 
1948  ncerr = nf90_put_var(ncid, znucl_id, znucl)
1949  NCF_CHECK_MSG(ncerr," write variable znucl")
1950 
1951  ncerr = nf90_put_var(ncid, amu_id, amu)
1952  NCF_CHECK_MSG(ncerr," write variable amu")
1953 
1954  ncerr = nf90_put_var(ncid, dtion_id, dtion)
1955  NCF_CHECK_MSG(ncerr," write variable dtion")
1956 
1957 #endif
1958 
1959 end subroutine write_csts_hist

m_abihist/write_md_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 write_md_hist

FUNCTION

 Write the history file into a netcdf file
 This version is not compatible with multiple images.

INPUTS

  filname=filename of the file where the history will be stored
  hist<type abihist>=Historical record of positions, forces, stresses, cell dims and energies,
  ifirst=1 if first access to the file
  itime = index of the step in the hist file
  natom=Number of atoms.
  nctime=NetCdf TIME between output of molecular dynamics information
  ntypat=Number of type of atoms.
  typat(natom)=Type of each natom
   amu(ntypat)=mass of the atoms (atomic mass unit)
  znucl(:)=Nuclear charge for each type of pseudopotential.
           WARNING: alchemical mixing is not supported. We assume npsp == ntypat
  dtion=time step for Molecular Dynamics

TODO

  Have you ever heard about ETSF-IO specifications?

OUTPUT

  (only writing)

SOURCE

1015 subroutine write_md_hist(hist,filename,ifirst,itime,natom,nctime,ntypat,&
1016 &                        typat,amu,znucl,dtion,mdtemp)
1017 
1018 !Arguments ------------------------------------
1019 !scalars
1020  integer,intent(in) :: ifirst,itime,natom,nctime,ntypat
1021  real(dp),intent(in) :: dtion
1022  character(len=*),intent(in) :: filename
1023 !arrays
1024  integer,intent(in) :: typat(natom)
1025  real(dp),intent(in) :: amu(ntypat),znucl(:),mdtemp(2)
1026  type(abihist),intent(inout),target :: hist
1027 
1028 !Local variables-------------------------------
1029 #if defined HAVE_NETCDF
1030 !scalars
1031  integer :: itime_file,ncerr,ncid,npsp
1032  integer :: xcart_id,xred_id,fcart_id,gred_id
1033  integer :: vel_id,vel_cell_id,etotal_id,acell_id,rprimd_id,strten_id
1034  integer :: ekin_id,entropy_id,mdtime_id
1035  logical :: has_nimage=.false.,need_to_write
1036  integer, parameter :: imgmov=0
1037 !arrays
1038 #endif
1039 
1040 ! *************************************************************************
1041 
1042 #if defined HAVE_NETCDF
1043 
1044  need_to_write = .FALSE.
1045  if(nctime==0 .or. ifirst==1) need_to_write = .TRUE.
1046  if (itime > nctime .and. nctime /= 0) then
1047      if (mod(itime,nctime) == 0) need_to_write = .TRUE.
1048  end if
1049 !Return if we don't need to write the HIST file at this step
1050  if (.not. need_to_write) return
1051 
1052  if (ifirst==1) then
1053 !##### First access: Create NetCDF file and write defs
1054 
1055    write(std_out,*) 'Write iteration in HIST netCDF file (also create it)'
1056    npsp=size(znucl)
1057 
1058 !  Create netCDF file
1059    ncerr = nf90_create(path=trim(filename),cmode=NF90_CLOBBER,ncid=ncid)
1060    NCF_CHECK_MSG(ncerr," create netcdf history file")
1061 
1062 !  Define all dims and vars
1063    call def_file_hist(ncid,natom,1,ntypat,npsp,has_nimage)
1064 
1065 !  Write variables that do not change
1066 !  (they are not read in a hist structure).
1067    call write_csts_hist(ncid,dtion,imgmov,typat,znucl,amu,mdtemp)
1068 
1069 !  Compute the itime for the hist file
1070    itime_file = 1
1071  else
1072 !##### itime>2 access: just open NetCDF file
1073 
1074    if(need_to_write) then
1075 
1076      write(std_out,*) 'Write iteration in HIST netCDF file'
1077 
1078 !    Open netCDF file
1079      ncerr = nf90_open(path=trim(filename),mode=NF90_WRITE, ncid=ncid)
1080      NCF_CHECK_MSG(ncerr," open netcdf history file")
1081 
1082 !    Compute the itime for the hist file
1083      itime_file = itime
1084      if(nctime > 0) itime_file = int(anint(real(itime / nctime,sp)))
1085 
1086    end if
1087  endif
1088 
1089  if(need_to_write) then
1090   !##### Write variables into the dataset
1091   !Get the IDs
1092    call get_varid_hist(ncid,xcart_id,xred_id,fcart_id,gred_id,vel_id,vel_cell_id,&
1093 &       rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id,has_nimage)
1094 !Write
1095    call write_vars_hist(ncid,hist,natom,has_nimage,1,itime_file,&
1096 &       xcart_id,xred_id,fcart_id,gred_id,vel_id,vel_cell_id,&
1097 &       rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id)
1098 
1099 !##### Close the file
1100    ncerr = nf90_close(ncid)
1101    NCF_CHECK_MSG(ncerr," close netcdf history file")
1102  end if
1103 #endif
1104 
1105 end subroutine write_md_hist

m_abihist/write_md_hist_img [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 write_md_hist_img

FUNCTION

 Write the history file into a netcdf file
 This version is compatible with multiple images of the

INPUTS

  filname= filename of the file where the history will be stored
  hist(:)<type abihist>= Historical record of positions, forces, stresses, cell dims and energies,
    Size(hist) is equal to a number of images to be written
  ifirst= 1 if first access to the file
  itime = index of the step in the hist file
  natom= Number of atoms.
  ntypat= Number of type of atoms.
  typat(natom)= Type of each natom
  amu(ntypat)= Mass of the atoms (atomic mass unit)
  znucl(:)= Nuclear charge for each type of pseudopotential.
           WARNING: alchemical mixing is not supported. We assume npsp == ntypat
  dtion= Time step for Molecular Dynamics
  [nimage]= Total number of images of the cell
  [comm_img]= MPI communicator over images of the cell
  [imgtab(:)]= In case of multiple images, indexes of images to be read
               Default is 1,2,3,...
               Size must be equal to size(hist)

TODO

  Have you ever heard about ETSF-IO specifications?

OUTPUT

  (only writing)

SOURCE

1146 subroutine write_md_hist_img(hist,filename,ifirst,itime,natom,ntypat,&
1147 &                            typat,amu,znucl,dtion,&
1148 &                            nimage,imgmov,mdtemp,comm_img,imgtab) ! optional arguments
1149 
1150 !Arguments ------------------------------------
1151 !scalars
1152  integer,intent(in) :: ifirst,itime,natom,ntypat
1153  integer,intent(in),optional :: nimage,imgmov,comm_img
1154  real(dp),intent(in) :: dtion
1155  character(len=*),intent(in) :: filename
1156 !arrays
1157  integer,intent(in) :: typat(natom)
1158  integer,intent(in),optional :: imgtab(:)
1159  real(dp),intent(in) :: amu(ntypat),znucl(:),mdtemp(2)
1160  type(abihist),intent(inout),target :: hist(:)
1161 
1162 !Local variables-------------------------------
1163 #if defined HAVE_NETCDF
1164 !scalars
1165  integer :: ii,iimage,iimg,me_img,my_comm_img,my_nimage,ncerr
1166  integer :: ncid,nimage_,nproc_img,npsp,imgmov_
1167  integer :: xcart_id,xred_id,fcart_id,gred_id
1168  integer :: vel_id,vel_cell_id,etotal_id
1169  integer :: acell_id,rprimd_id,strten_id
1170  integer :: ekin_id,entropy_id,mdtime_id
1171  logical :: has_nimage, has_imgmov
1172  !character(len=500) :: msg
1173  type(abihist),pointer :: hist_
1174 !arrays
1175  integer,allocatable :: my_imgtab(:)
1176 #endif
1177 
1178 ! *************************************************************************
1179 
1180 #if defined HAVE_NETCDF
1181 
1182 !Manage multiple images of the cell
1183  has_nimage=present(nimage)
1184  has_imgmov=present(imgmov)
1185  nimage_=merge(nimage,1,has_nimage)
1186  imgmov_=merge(imgmov,0,has_imgmov)
1187  my_nimage=size(hist) ; if (my_nimage==0) return
1188  my_comm_img=xmpi_comm_self;if(present(comm_img)) my_comm_img=comm_img
1189  nproc_img=xmpi_comm_size(my_comm_img)
1190  me_img=xmpi_comm_rank(my_comm_img)
1191  ABI_MALLOC(my_imgtab,(my_nimage))
1192  if (present(imgtab)) then
1193   if (size(my_imgtab)/=my_nimage) then
1194     ABI_BUG('Inconsistency between hist and imgtab!')
1195   end if
1196   my_imgtab(:)=imgtab(:)
1197  else
1198    my_imgtab(:)=(/(iimage,iimage=1,my_nimage)/)
1199  end if
1200 
1201 !Has to access the HIST file sequentially, proc by proc
1202  do ii=0,nproc_img-1
1203    call xmpi_barrier(my_comm_img)
1204    if (me_img==ii) then
1205 
1206 !    ##### First access: Create NetCDF file and write defs
1207      if (ifirst==1.and.me_img==0) then
1208        npsp=size(znucl)
1209 !      Create netCDF file
1210        ncerr = nf90_create(path=trim(filename),cmode=NF90_CLOBBER,ncid=ncid)
1211        NCF_CHECK_MSG(ncerr," create netcdf history file")
1212 !      Define all dims and vars
1213        call def_file_hist(ncid,natom,nimage_,ntypat,npsp,has_nimage)
1214 !      Write variables that do not change
1215 !      (they are not read in a hist structure).
1216        call write_csts_hist(ncid,dtion,imgmov_,typat,znucl,amu,mdtemp)
1217      end if
1218 
1219 !    ##### itime>2 access: just open NetCDF file
1220      if (ifirst/=1.or.me_img/=0) then
1221 !      Open netCDF file
1222        ncerr = nf90_open(path=trim(filename),mode=NF90_WRITE, ncid=ncid)
1223        NCF_CHECK_MSG(ncerr," open netcdf history file")
1224      end if
1225 
1226      write(std_out,*) 'Write iteration in HIST netCDF file'
1227 
1228 !    ##### Write variables into the dataset (loop over images)
1229 !    Get the IDs
1230      call get_varid_hist(ncid,xcart_id,xred_id,fcart_id,gred_id,vel_id,vel_cell_id,&
1231 &         rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id,has_nimage)
1232 
1233 !    Write
1234      do iimage=1,my_nimage
1235        iimg=my_imgtab(iimage)
1236        hist_ => hist(iimage)
1237        call write_vars_hist(ncid,hist_,natom,has_nimage,iimg,itime,&
1238 &           xcart_id,xred_id,fcart_id,gred_id,vel_id,vel_cell_id,&
1239 &           rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id)
1240      end do
1241 
1242 !    ##### Close the file
1243      ncerr = nf90_close(ncid)
1244      NCF_CHECK_MSG(ncerr," close netcdf history file")
1245      ABI_FREE(my_imgtab)
1246 
1247 !  End loop on MPI processes
1248    end if
1249  end do
1250 
1251 #endif
1252 
1253 end subroutine write_md_hist_img

m_abihist/write_vars_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 write_vars_hist

FUNCTION

INPUTS

OUTPUT

SOURCE

1977 subroutine write_vars_hist(ncid,hist,natom,has_nimage,iimg,itime,&
1978 &          xcart_id,xred_id,fcart_id,gred_id,vel_id,vel_cell_id,rprimd_id,&
1979 &          acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id)
1980 
1981 !Arguments ------------------------------------
1982 !scalars
1983  integer,intent(in) :: ncid,natom,iimg,itime
1984  integer,intent(in) :: xcart_id,xred_id,fcart_id,gred_id,vel_id
1985  integer,intent(in) :: vel_cell_id,rprimd_id,acell_id,strten_id
1986  integer,intent(in) :: etotal_id,ekin_id,entropy_id,mdtime_id
1987  logical,intent(in) :: has_nimage
1988  type(abihist),intent(inout),target :: hist
1989 
1990 !Local variables-------------------------------
1991 #if defined HAVE_NETCDF
1992 !scalars
1993  integer :: ncerr
1994 !arrays
1995  integer :: count2(2),count3(3),count4(4)
1996  integer :: start1(1),start2(2),start3(3),start4(4)
1997  real(dp),allocatable :: conv(:,:)
1998  real(dp),pointer :: xred(:,:),fcart(:,:),rprimd(:,:),vel(:,:),vel_cell(:,:)
1999 #endif
2000 
2001 ! *************************************************************************
2002 
2003 #if defined HAVE_NETCDF
2004 
2005  xred     => hist%xred(:,:,hist%ihist)
2006  fcart    => hist%fcart(:,:,hist%ihist)
2007  vel      => hist%vel(:,:,hist%ihist)
2008  vel_cell => hist%vel_cell(:,:,hist%ihist)
2009  rprimd   => hist%rprimd(:,:,hist%ihist)
2010 
2011 !Variables not depending on images
2012 
2013 !mdtime
2014  start1=(/itime/)
2015  ncerr = nf90_put_var(ncid,mdtime_id,hist%time(hist%ihist),start=start1)
2016  NCF_CHECK_MSG(ncerr," write variable mdtime")
2017 
2018 !Variables depending on images
2019 
2020  ABI_MALLOC(conv,(3,natom))
2021 
2022 !xcart,xred,fcart,gred,vel
2023  if (has_nimage) then
2024    start4=(/1,1,iimg,itime/);count4=(/3,natom,1,1/)
2025    call xred2xcart(natom,rprimd,conv,xred)
2026    ncerr = nf90_put_var(ncid,xcart_id,conv, start = start4,count = count4)
2027    NCF_CHECK_MSG(ncerr," write variable xcart")
2028    ncerr = nf90_put_var(ncid,xred_id,xred, start = start4,count = count4)
2029    NCF_CHECK_MSG(ncerr," write variable xred")
2030    ncerr = nf90_put_var(ncid,fcart_id,fcart,start = start4,count = count4)
2031    NCF_CHECK_MSG(ncerr," write variable fcart")
2032    call fcart2gred(fcart,conv,rprimd,natom)
2033    ncerr = nf90_put_var(ncid,gred_id,conv, start = start4,count = count4)
2034    NCF_CHECK_MSG(ncerr," write variable fred") ! XG210315 : should be "gred"
2035    ncerr = nf90_put_var(ncid,vel_id,vel, start = start4,count = count4)
2036    NCF_CHECK_MSG(ncerr," write variable vel")
2037  else
2038    start3=(/1,1,itime/);count3=(/3,natom,1/)
2039    call xred2xcart(natom,rprimd,conv,xred)
2040    ncerr = nf90_put_var(ncid,xcart_id,conv, start = start3,count = count3)
2041    NCF_CHECK_MSG(ncerr," write variable xcart")
2042    ncerr = nf90_put_var(ncid,xred_id,xred, start = start3,count = count3)
2043    NCF_CHECK_MSG(ncerr," write variable xred")
2044    ncerr = nf90_put_var(ncid,fcart_id,fcart,start = start3,count = count3)
2045    NCF_CHECK_MSG(ncerr," write variable fcart")
2046    call fcart2gred(fcart,conv,rprimd,natom)
2047    ncerr = nf90_put_var(ncid,gred_id,conv, start = start3,count = count3)
2048    NCF_CHECK_MSG(ncerr," write variable fred") ! XG210315 : should be "gred"
2049    ncerr = nf90_put_var(ncid,vel_id,vel, start = start3,count = count3)
2050    NCF_CHECK_MSG(ncerr," write variable vel")
2051  end if
2052 
2053  ABI_FREE(conv)
2054 
2055 !rprimd,vel_cell
2056  if (has_nimage) then
2057    start4=(/1,1,iimg,itime/);count4=(/3,3,1,1/)
2058    ncerr = nf90_put_var(ncid,rprimd_id,hist%rprimd(:,:,hist%ihist),&
2059 &                       start = start4,count = count4)
2060    NCF_CHECK_MSG(ncerr," write variable rprimd")
2061    ncerr = nf90_put_var(ncid,vel_cell_id,hist%vel_cell(:,:,hist%ihist),&
2062 &                       start = start4,count = count4)
2063    NCF_CHECK_MSG(ncerr," write variable vel_cell")
2064  else
2065    start3=(/1,1,itime/);count3=(/3,3,1/)
2066    ncerr = nf90_put_var(ncid,rprimd_id,hist%rprimd(:,:,hist%ihist),&
2067 &                       start = start3,count = count3)
2068    NCF_CHECK_MSG(ncerr," write variable rprimd")
2069  end if
2070 
2071 !acell
2072  if (has_nimage) then
2073    start3=(/1,iimg,itime/);count3=(/3,1,1/)
2074    ncerr = nf90_put_var(ncid,acell_id,hist%acell(:,hist%ihist),&
2075 &                       start = start3,count = count3)
2076    NCF_CHECK_MSG(ncerr," write variable acell")
2077  else
2078    start2=(/1,itime/);count2=(/3,1/)
2079    ncerr = nf90_put_var(ncid,acell_id,hist%acell(:,hist%ihist),&
2080 &                       start = start2,count = count2)
2081    NCF_CHECK_MSG(ncerr," write variable acell")
2082  end if
2083 
2084 !strten
2085  if (has_nimage) then
2086    start3=(/1,iimg,itime/);count3=(/6,1,1/)
2087    ncerr = nf90_put_var(ncid,strten_id,hist%strten(:,hist%ihist),&
2088 &                       start = start3,count = count3)
2089    NCF_CHECK_MSG(ncerr," write variable strten")
2090  else
2091    start2=(/1,itime/);count2=(/6,1/)
2092    ncerr = nf90_put_var(ncid,strten_id,hist%strten(:,hist%ihist),&
2093 &                       start = start2,count = count2)
2094    NCF_CHECK_MSG(ncerr," write variable strten")
2095  end if
2096 
2097 !etotal,ekin,entropy
2098  if (has_nimage) then
2099    start2=(/iimg,itime/)
2100    ncerr = nf90_put_var(ncid,etotal_id,hist%etot(hist%ihist),start=start2)
2101    NCF_CHECK_MSG(ncerr," write variable etotal")
2102    ncerr = nf90_put_var(ncid,ekin_id,hist%ekin(hist%ihist),start=start2)
2103    NCF_CHECK_MSG(ncerr," write variable ekin")
2104    ncerr = nf90_put_var(ncid,entropy_id,hist%entropy(hist%ihist),start=start2)
2105    NCF_CHECK_MSG(ncerr," write variable entropy")
2106  else
2107    start1=(/itime/)
2108    ncerr = nf90_put_var(ncid,etotal_id,hist%etot(hist%ihist),start=start1)
2109    NCF_CHECK_MSG(ncerr," write variable etotal")
2110    ncerr = nf90_put_var(ncid,ekin_id,hist%ekin(hist%ihist),start=start1)
2111    NCF_CHECK_MSG(ncerr," write variable ekin")
2112    ncerr = nf90_put_var(ncid,entropy_id,hist%entropy(hist%ihist),start=start1)
2113    NCF_CHECK_MSG(ncerr," write variable entropy")
2114  end if
2115 
2116 #endif
2117 
2118 end subroutine write_vars_hist