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-2018 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

32 #if defined HAVE_CONFIG_H
33 #include "config.h"
34 #endif
35 
36 #include "abi_common.h"
37 
38 module m_abihist
39 
40  use defs_basis
41  use m_abicore
42  use m_errors
43  use m_xmpi
44  use m_nctk
45 #if defined HAVE_NETCDF
46  use netcdf
47 #endif
48 
49  use m_geometry,  only : fcart2fred, xred2xcart
50 
51  implicit none
52 
53  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

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

PARENTS

      m_abihist

CHILDREN

NOTES

SOURCE

472 subroutine abihist_bcast_0D(hist,master,comm)
473 
474 
475 !This section has been created automatically by the script Abilint (TD).
476 !Do not modify the following lines by hand.
477 #undef ABI_FUNC
478 #define ABI_FUNC 'abihist_bcast_0D'
479 !End of the abilint section
480 
481  implicit none
482 
483 !Arguments ------------------------------------
484 !scalars
485  integer,intent(in) :: master,comm
486  type(abihist),intent(inout) :: hist
487 
488 !Local variables-------------------------------
489 !scalars
490  integer :: bufsize,ierr,indx,nproc,rank
491  integer :: sizeA,sizeA1,sizeA2
492  integer :: sizeEt,sizeEk,sizeEnt,sizeT
493  integer :: sizeR,sizeR1,sizeR2,sizeR3
494  integer :: sizeS,sizeS1,sizeS2
495  integer :: sizeV,sizeV1,sizeV2,sizeV3
496  integer :: sizeVc,sizeVc1,sizeVc2,sizeVc3
497  integer :: sizeX,sizeX1,sizeX2,sizeX3
498  integer :: sizeF,sizeF1,sizeF2,sizeF3
499 !arrays
500  integer,allocatable :: buffer_i(:)
501  real(dp),allocatable :: buffer_r(:)
502 
503 ! ***************************************************************
504 
505  ierr=0
506  nproc=xmpi_comm_size(comm)
507  if (nproc<=1) return
508 
509  rank=xmpi_comm_rank(comm)
510 
511 !=== Broadcast integers and logicals
512  ABI_ALLOCATE(buffer_i,(4))
513  if (rank==master) then
514    buffer_i(1)=hist%ihist
515    buffer_i(2)=hist%mxhist
516    buffer_i(3)=0;if (hist%isVused)  buffer_i(3)=1
517    buffer_i(4)=0;if (hist%isARused) buffer_i(4)=1
518  end if
519  call xmpi_bcast(buffer_i,master,comm,ierr)
520  if (rank/=master) then
521    hist%ihist=buffer_i(1)
522    hist%mxhist=buffer_i(2)
523    hist%isVused  =(buffer_i(3)==1)
524    hist%isARused =(buffer_i(4)==1)
525  end if
526  ABI_DEALLOCATE(buffer_i)
527 
528 !If history is empty, return
529  if (hist%mxhist==0.or.hist%ihist==0) return
530 
531 !=== Broadcast sizes of arrays
532  ABI_ALLOCATE(buffer_i,(23))
533  if (rank==master) then
534    sizeA1=size(hist%acell,1);sizeA2=size(hist%acell,2)
535    sizeEt=size(hist%etot,1);sizeEk=size(hist%ekin,1);
536    sizeEnt=size(hist%entropy,1);sizeT=size(hist%time,1)
537    sizeR1=size(hist%rprimd,1);sizeR2=size(hist%rprimd,2);sizeR3=size(hist%rprimd,3)
538    sizeS1=size(hist%strten,1);sizeS2=size(hist%strten,2)
539    sizeV1=size(hist%vel,1);sizeV2=size(hist%vel,2);sizeV3=size(hist%vel,3)
540    sizeVc1=size(hist%vel_cell,1);sizeVc2=size(hist%vel_cell,2);sizeVc3=size(hist%vel_cell,3)
541    sizeX1=size(hist%xred,1);sizeX2=size(hist%xred,2);sizeX3=size(hist%xred,3)
542    sizeF1=size(hist%fcart,1);sizeF2=size(hist%fcart,2);sizeF3=size(hist%fcart,3)
543    buffer_i(1)=sizeA1  ;buffer_i(2)=sizeA2
544    buffer_i(3)=sizeEt  ;buffer_i(4)=sizeEk
545    buffer_i(5)=sizeEnt ;buffer_i(6)=sizeT
546    buffer_i(7)=sizeR1  ;buffer_i(8)=sizeR2
547    buffer_i(9)=sizeR3  ;buffer_i(10)=sizeS1
548    buffer_i(11)=sizeS2 ;buffer_i(12)=sizeV1
549    buffer_i(13)=sizeV2 ;buffer_i(14)=sizeV3
550    buffer_i(15)=sizeVc1;buffer_i(16)=sizeVc2
551    buffer_i(17)=sizeVc3;buffer_i(18)=sizeX1
552    buffer_i(19)=sizeX2 ;buffer_i(20)=sizeX3
553    buffer_i(21)=sizeF1 ;buffer_i(22)=sizeF2
554    buffer_i(23)=sizeF3
555  end if
556  call xmpi_bcast(buffer_i,master,comm,ierr)
557 
558  if (rank/=master) then
559    sizeA1 =buffer_i(1) ;sizeA2 =buffer_i(2)
560    sizeEt =buffer_i(3) ;sizeEk =buffer_i(4)
561    sizeEnt=buffer_i(5);sizeT   =buffer_i(6)
562    sizeR1 =buffer_i(7) ;sizeR2 =buffer_i(8)
563    sizeR3 =buffer_i(9) ;sizeS1 =buffer_i(10)
564    sizeS2 =buffer_i(11);sizeV1 =buffer_i(12)
565    sizeV2 =buffer_i(13);sizeV3 =buffer_i(14)
566    sizeVc1=buffer_i(15);sizeVc2=buffer_i(16)
567    sizeVc3=buffer_i(17);sizeX1 =buffer_i(18)
568    sizeX2 =buffer_i(19);sizeX3 =buffer_i(20)
569    sizeF1 =buffer_i(21);sizeF2 =buffer_i(22)
570    sizeF3 =buffer_i(23)
571  end if
572  ABI_DEALLOCATE(buffer_i)
573 
574 !=== Broadcast reals
575  sizeA=sizeA1*sizeA2;sizeR=sizeR1*sizeR2*sizeR3;sizeS=sizeS1*sizeS2
576  sizeV=sizeV1*sizeV2*sizeV3;sizeVc=sizeVc1*sizeVc2*sizeVc3;
577  sizeX=sizeX1*sizeX2*sizeX3;sizeF=sizeF1*sizeF2*sizeF3
578  bufsize=sizeA+sizeEt+sizeEk+sizeEnt+sizeT+sizeR+sizeS+sizeV+sizeVc+sizeX+sizeF
579  ABI_ALLOCATE(buffer_r,(bufsize))
580  if (rank==master) then
581    indx=0
582    buffer_r(indx+1:indx+sizeA)=reshape(hist%acell(1:sizeA1,1:sizeA2),(/sizeA/))
583    indx=indx+sizeA
584    buffer_r(indx+1:indx+sizeEt)=hist%etot(1:sizeEt)
585    indx=indx+sizeEt
586    buffer_r(indx+1:indx+sizeEk)=hist%ekin(1:sizeEk)
587    indx=indx+sizeEk
588    buffer_r(indx+1:indx+sizeEnt)=hist%entropy(1:sizeEnt)
589    indx=indx+sizeEnt
590    buffer_r(indx+1:indx+sizeT)=hist%time(1:sizeT)
591    indx=indx+sizeT
592    buffer_r(indx+1:indx+sizeR)=reshape(hist%rprimd(1:sizeR1,1:sizeR2,1:sizeR3),(/sizeR/))
593    indx=indx+sizeR
594    buffer_r(indx+1:indx+sizeS)=reshape(hist%strten(1:sizeS1,1:sizeS2),(/sizeS/))
595    indx=indx+sizeS
596    buffer_r(indx+1:indx+sizeV)=reshape(hist%vel(1:sizeV1,1:sizeV2,1:sizeV3),(/sizeV/))
597    indx=indx+sizeV
598    buffer_r(indx+1:indx+sizeVc)=reshape(hist%vel_cell(1:sizeVc1,1:sizeVc2,1:sizeVc3),(/sizeVc/))
599    indx=indx+sizeVc
600    buffer_r(indx+1:indx+sizeX)=reshape(hist%xred(1:sizeX1,1:sizeX2,1:sizeX3),(/sizeX/))
601    indx=indx+sizeX
602    buffer_r(indx+1:indx+sizeF)=reshape(hist%fcart(1:sizeF1,1:sizeF2,1:sizeF3),(/sizeF/))
603  else
604    call abihist_free(hist)
605    ABI_ALLOCATE(hist%acell,(sizeA1,sizeA2))
606    ABI_ALLOCATE(hist%etot,(sizeEt))
607    ABI_ALLOCATE(hist%ekin,(sizeEk))
608    ABI_ALLOCATE(hist%entropy,(sizeEnt))
609    ABI_ALLOCATE(hist%time,(sizeT))
610    ABI_ALLOCATE(hist%rprimd,(sizeR1,sizeR2,sizeR3))
611    ABI_ALLOCATE(hist%strten,(sizeS1,sizeS2))
612    ABI_ALLOCATE(hist%vel,(sizeV1,sizeV2,sizeV3))
613    ABI_ALLOCATE(hist%vel_cell,(sizeVc1,sizeVc2,sizeVc3))
614    ABI_ALLOCATE(hist%xred,(sizeX1,sizeX2,sizeX3))
615    ABI_ALLOCATE(hist%fcart,(sizeF1,sizeF2,sizeF3))
616  end if
617  call xmpi_bcast(buffer_r,master,comm,ierr)
618 
619  if (rank/=master) then
620    indx=0
621    hist%acell(1:sizeA1,1:sizeA2)=reshape(buffer_r(indx+1:indx+sizeA), &
622 &                                       (/sizeA1,sizeA2/))
623    indx=indx+sizeA
624    hist%etot(1:sizeEt)=buffer_r(indx+1:indx+sizeEt)
625    indx=indx+sizeEt
626    hist%ekin(1:sizeEk)=buffer_r(indx+1:indx+sizeEk)
627    indx=indx+sizeEk
628    hist%entropy(1:sizeEnt)=buffer_r(indx+1:indx+sizeEnt)
629    indx=indx+sizeEnt
630    hist%time(1:sizeT)=buffer_r(indx+1:indx+sizeT)
631    indx=indx+sizeT
632    hist%rprimd(1:sizeR1,1:sizeR2,1:sizeR3)=reshape(buffer_r(indx+1:indx+sizeR), &
633 &                                                 (/sizeR1,sizeR2,sizeR3/))
634    indx=indx+sizeR
635    hist%strten(1:sizeS1,1:sizeS2)=reshape(buffer_r(indx+1:indx+sizeS), &
636 &                                        (/sizeS1,sizeS2/))
637    indx=indx+sizeS
638    hist%vel(1:sizeV1,1:sizeV2,1:sizeV3)=reshape(buffer_r(indx+1:indx+sizeV), &
639 &                                              (/sizeV1,sizeV2,sizeV3/))
640    indx=indx+sizeV
641    hist%vel_cell(1:sizeVc1,1:sizeVc2,1:sizeVc3)=reshape(buffer_r(indx+1:indx+sizeVc), &
642 &                                                      (/sizeVc1,sizeVc2,sizeVc3/))
643    indx=indx+sizeVc
644    hist%xred(1:sizeX1,1:sizeX2,1:sizeX3)=reshape(buffer_r(indx+1:indx+sizeX), &
645 &                                               (/sizeX1,sizeX2,sizeX3/))
646    indx=indx+sizeX
647    hist%fcart(1:sizeF1,1:sizeF2,1:sizeF3)=reshape(buffer_r(indx+1:indx+sizeF), &
648 &                                                (/sizeF1,sizeF2,sizeF3/))
649  end if
650  ABI_DEALLOCATE(buffer_r)
651 
652 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

PARENTS

CHILDREN

NOTES

SOURCE

681 subroutine abihist_bcast_1D(hist,master,comm)
682 
683 
684 !This section has been created automatically by the script Abilint (TD).
685 !Do not modify the following lines by hand.
686 #undef ABI_FUNC
687 #define ABI_FUNC 'abihist_bcast_1D'
688 !End of the abilint section
689 
690  implicit none
691 
692 !Arguments ------------------------------------
693 !scalars
694  integer,intent(in) :: master,comm
695  type(abihist),intent(inout) :: hist(:)
696 
697 !Local variables-------------------------------
698  integer :: ii
699 
700 ! ***************************************************************
701 
702  do ii=1,size(hist)
703    call abihist_bcast_0D(hist(ii),master,comm)
704  end do
705 
706 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)>

PARENTS

      mover

CHILDREN

SOURCE

1128 subroutine abihist_compare_and_copy(hist_in,hist_out,natom,similar,tolerance,store_all)
1129 
1130 
1131 !This section has been created automatically by the script Abilint (TD).
1132 !Do not modify the following lines by hand.
1133 #undef ABI_FUNC
1134 #define ABI_FUNC 'abihist_compare_and_copy'
1135 !End of the abilint section
1136 
1137  implicit none
1138 
1139 !Arguments ------------------------------------
1140 !scalars
1141 integer,intent(in) :: natom
1142 integer,intent(out) :: similar
1143 real(dp),intent(in) :: tolerance
1144 type(abihist),intent(in) :: hist_in
1145 type(abihist),intent(inout) :: hist_out
1146 logical,intent(in) :: store_all
1147 !Local variables-------------------------------
1148 !scalars
1149 integer :: kk,jj
1150 real(dp) :: maxdiff,diff
1151 real(dp) :: x,y
1152 !array
1153 character(len= 500) :: msg
1154 ! ***************************************************************
1155 
1156  ABI_UNUSED(store_all)
1157 
1158  similar=1
1159 
1160  write(msg,'(a,I0,4a)')  'Using values from history, iteration:',hist_in%ihist,ch10,&
1161 &                     'Differences between present history and values stored',ch10,&
1162 &                     'on the previous history.(Relative difference)'
1163  call wrtout(std_out,msg,'COLL')
1164 
1165  x=hist_out%xred(1,1,hist_out%ihist)
1166  y=hist_in%xred(1,1,hist_in%ihist)
1167  maxdiff=2*abs(x-y)/(abs(x)+abs(y))
1168  do kk=1,natom
1169    do jj=1,3
1170      x=hist_out%xred(jj,kk,hist_out%ihist)
1171      y=hist_in%xred(jj,kk,hist_in%ihist)
1172      diff=2*abs(x-y)/(abs(x)+abs(y))
1173      if (diff>maxdiff) maxdiff=diff
1174    end do
1175  end do
1176  write(msg,'(a,e12.5)') 'xred:     ',maxdiff
1177  call wrtout(std_out,msg,'COLL')
1178 
1179  if (maxdiff>tolerance) similar=0
1180 
1181 
1182  x=hist_out%rprimd(1,1,hist_out%ihist)
1183  y=hist_in%rprimd(1,1,hist_in%ihist)
1184  maxdiff=2*abs(x-y)/(abs(x)+abs(y))
1185  do kk=1,3
1186    do jj=1,3
1187      x=hist_out%rprimd(jj,kk,hist_out%ihist)
1188      y=hist_in%rprimd(jj,kk,hist_in%ihist)
1189      diff=2*abs(x-y)/(abs(x)+abs(y))
1190      if (diff>maxdiff) maxdiff=diff
1191    end do
1192  end do
1193  write(msg,'(a,e12.5)') 'rprimd:   ',maxdiff
1194  call wrtout(std_out,msg,'COLL')
1195  if (maxdiff>tolerance) similar=0
1196 
1197 
1198  x=hist_out%acell(1,hist_out%ihist)
1199  y=hist_in%acell(1,hist_in%ihist)
1200  maxdiff=2*abs(x-y)/(abs(x)+abs(y))
1201  do kk=1,3
1202    x=hist_out%acell(kk,hist_out%ihist)
1203    y=hist_in%acell(kk,hist_in%ihist)
1204    diff=2*abs(x-y)/(abs(x)+abs(y))
1205    if (diff>maxdiff) maxdiff=diff
1206  end do
1207  write(msg,'(a,e12.5)') 'acell:    ',maxdiff
1208  call wrtout(std_out,msg,'COLL')
1209  if (maxdiff>tolerance) similar=0
1210 
1211  if (similar==1) then
1212    hist_out%acell(:,hist_out%ihist)     =hist_in%acell(:,hist_in%ihist)
1213    hist_out%rprimd(:,:,hist_out%ihist)  =hist_in%rprimd(:,:,hist_in%ihist)
1214    hist_out%xred(:,:,hist_out%ihist)    =hist_in%xred(:,:,hist_in%ihist)
1215    hist_out%fcart(:,:,hist_out%ihist)   =hist_in%fcart(:,:,hist_in%ihist)
1216    hist_out%strten(:,hist_out%ihist)    =hist_in%strten(:,hist_in%ihist)
1217    hist_out%vel(:,:,hist_out%ihist)     =hist_in%vel(:,:,hist_in%ihist)
1218    hist_out%vel_cell(:,:,hist_out%ihist)=hist_in%vel_cell(:,:,hist_in%ihist)
1219    hist_out%etot(hist_out%ihist)        =hist_in%etot(hist_in%ihist)
1220    hist_out%ekin(hist_out%ihist)        =hist_in%ekin(hist_in%ihist)
1221    hist_out%entropy(hist_out%ihist)     =hist_in%entropy(hist_in%ihist)
1222    hist_out%time(hist_out%ihist)        =hist_in%time(hist_in%ihist)
1223  end if
1224 
1225 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)>

PARENTS

      gstateimg,m_effective_potential_file

CHILDREN

SOURCE

1051 subroutine abihist_copy(hist_in,hist_out)
1052 
1053 
1054 !This section has been created automatically by the script Abilint (TD).
1055 !Do not modify the following lines by hand.
1056 #undef ABI_FUNC
1057 #define ABI_FUNC 'abihist_copy'
1058 !End of the abilint section
1059 
1060  implicit none
1061 
1062 !Arguments ------------------------------------
1063 !scalars
1064 type(abihist),intent(in) :: hist_in
1065 type(abihist),intent(inout) :: hist_out
1066 
1067 !Local variables-------------------------------
1068 !scalars
1069  character(len=500) :: msg
1070 
1071 ! ***************************************************************
1072 
1073 !Check
1074  if (size(hist_in%xred,2)/=size(hist_out%xred,2)) then
1075    msg='Incompatible sizes for hist_in and hist_out!'
1076    MSG_BUG(msg)
1077  end if
1078 
1079 !Copy scalars (except ihist and mxhist)
1080  hist_out%isVused  =hist_in%isVused
1081  hist_out%isARused =hist_in%isARused
1082 
1083 !Copy arrays
1084  hist_out%acell(:,hist_out%ihist)     = hist_in%acell(:,hist_in%ihist)
1085  hist_out%rprimd(:,:,hist_out%ihist)  = hist_in%rprimd(:,:,hist_in%ihist)
1086  hist_out%xred(:,:,hist_out%ihist)    = hist_in%xred(:,:,hist_in%ihist)
1087  hist_out%fcart(:,:,hist_out%ihist)   = hist_in%fcart(:,:,hist_in%ihist)
1088  hist_out%strten(:,hist_out%ihist)    = hist_in%strten(:,hist_in%ihist)
1089  hist_out%vel(:,:,hist_out%ihist)     = hist_in%vel(:,:,hist_in%ihist)
1090  hist_out%vel_cell(:,:,hist_out%ihist)= hist_in%vel_cell(:,:,hist_in%ihist)
1091  hist_out%etot(hist_out%ihist)        = hist_in%etot(hist_in%ihist)
1092  hist_out%ekin(hist_out%ihist)        = hist_in%ekin(hist_in%ihist)
1093  hist_out%entropy(hist_out%ihist)     = hist_in%entropy(hist_in%ihist)
1094  hist_out%time(hist_out%ihist)        = hist_in%time(hist_in%ihist)
1095 
1096 end subroutine abihist_copy

m_abihist/abihist_fi [ 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

SIDE EFFECTS

PARENTS

CHILDREN

SOURCE

816 function abihist_findIndex(hist,step) result(index)
817 
818 
819 !This section has been created automatically by the script Abilint (TD).
820 !Do not modify the following lines by hand.
821 #undef ABI_FUNC
822 #define ABI_FUNC 'abihist_findIndex'
823 !End of the abilint section
824 
825  implicit none
826 
827 !Arguments ------------------------------------
828 !scalars
829  integer,intent(in) :: step
830  integer :: index
831 !arrays
832  type(abihist),intent(in) :: hist
833 !Local variables-------------------------------
834 !scalars
835  integer :: ii,mxhist
836 !arrays
837  character(len=500) :: msg
838 ! *************************************************************
839 
840  mxhist = hist%mxhist
841 
842  if ((mxhist ==1.and.step/=+1).or.&
843 &    (mxhist /=1.and.abs(step) >=mxhist)) then
844    write(msg,'(a,I0,2a)')' The requested step must be less than ',mxhist,ch10,&
845 &                     'Action: increase the number of history stored in the hist'
846    MSG_BUG(msg)
847  end if
848 
849  ii = hist%ihist + step
850 
851  do while (ii > mxhist)
852    ii = ii - mxhist
853  end do
854  do while (ii <= 0)
855    ii = ii + mxhist
856  end do
857 
858  index = ii
859 
860 end function abihist_findIndex

m_abihist/abihist_free_0D [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 abihist_free_0D

FUNCTION

 Deallocate all the pointers in a hist structure - Target: scalar

INPUTS

OUTPUT

SIDE EFFECTS

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

PARENTS

      m_abihist

CHILDREN

NOTES

SOURCE

331 subroutine abihist_free_0D(hist)
332 
333 
334 !This section has been created automatically by the script Abilint (TD).
335 !Do not modify the following lines by hand.
336 #undef ABI_FUNC
337 #define ABI_FUNC 'abihist_free_0D'
338 !End of the abilint section
339 
340  implicit none
341 
342 !Arguments ------------------------------------
343  type(abihist),intent(inout) :: hist
344 
345 ! ***************************************************************
346 
347 !Vector of (mxhist) values of cell dimensions
348  if (allocated(hist%acell))  then
349    ABI_DEALLOCATE(hist%acell)
350  end if
351 !Vector of (mxhist) values of cell primitive translations
352  if (allocated(hist%rprimd))  then
353    ABI_DEALLOCATE(hist%rprimd)
354  end if
355 !Vector of (mxhist) values of reduced coordinates
356  if (allocated(hist%xred))  then
357    ABI_DEALLOCATE(hist%xred)
358  end if
359 !Vector of (mxhist) values of cartesian forces
360  if (allocated(hist%fcart))  then
361    ABI_DEALLOCATE(hist%fcart)
362  end if
363 !Vector of (mxhist) values of stress tensor
364  if (allocated(hist%strten))  then
365    ABI_DEALLOCATE(hist%strten)
366  end if
367 ! Vector of (mxhist) values of atomic velocities
368  if (allocated(hist%vel))  then
369    ABI_DEALLOCATE(hist%vel)
370  end if
371 ! Vector of (mxhist) values of cell velocities
372  if (allocated(hist%vel_cell))  then
373    ABI_DEALLOCATE(hist%vel_cell)
374  end if
375 ! Vector of (mxhist) values of electronic total energy
376  if (allocated(hist%etot))  then
377    ABI_DEALLOCATE(hist%etot)
378  end if
379 ! Vector of (mxhist) values of ionic kinetic energy
380  if (allocated(hist%ekin))  then
381    ABI_DEALLOCATE(hist%ekin)
382  end if
383 ! Vector of (mxhist) values of Entropy
384  if (allocated(hist%entropy))  then
385    ABI_DEALLOCATE(hist%entropy)
386  end if
387 ! Vector of (mxhist) values of time (relevant for MD calculations)
388  if (allocated(hist%time))  then
389    ABI_DEALLOCATE(hist%time)
390  end if
391 
392 end subroutine abihist_free_0D

m_abihist/abihist_free_1D [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 abihist_free_1D

FUNCTION

 Deallocate all the pointers in a hist structure - Target: 1D array

INPUTS

OUTPUT

SIDE EFFECTS

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

PARENTS

CHILDREN

NOTES

SOURCE

419 subroutine abihist_free_1D(hist)
420 
421 
422 !This section has been created automatically by the script Abilint (TD).
423 !Do not modify the following lines by hand.
424 #undef ABI_FUNC
425 #define ABI_FUNC 'abihist_free_1D'
426 !End of the abilint section
427 
428  implicit none
429 
430 !Arguments ------------------------------------
431  type(abihist),intent(inout) :: hist(:)
432 
433 !Local variables-------------------------------
434  integer :: ii
435 
436 ! ***************************************************************
437 
438  do ii=1,size(hist)
439    call abihist_free_0D(hist(ii))
440  end do
441 
442 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 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

SIDE EFFECTS

PARENTS

      m_abihist

CHILDREN

NOTES

SOURCE

196 subroutine abihist_init_0D(hist,natom,mxhist,isVused,isARused)
197 
198 
199 !This section has been created automatically by the script Abilint (TD).
200 !Do not modify the following lines by hand.
201 #undef ABI_FUNC
202 #define ABI_FUNC 'abihist_init_0D'
203 !End of the abilint section
204 
205  implicit none
206 
207 !Arguments ------------------------------------
208  integer,intent(in) :: natom,mxhist
209  logical,intent(in) :: isVUsed,isARused
210  type(abihist),intent(inout) :: hist
211 
212 ! ***************************************************************
213 
214 !Initialize indexes
215  hist%ihist=1
216  hist%mxhist=mxhist
217 
218 !Initialize flags
219  hist%isVused=isVUsed
220  hist%isARused=isARUsed
221 
222 !Allocate all the histories
223  ABI_ALLOCATE(hist%acell,(3,mxhist))
224  ABI_ALLOCATE(hist%rprimd,(3,3,mxhist))
225  ABI_ALLOCATE(hist%xred,(3,natom,mxhist))
226  ABI_ALLOCATE(hist%fcart,(3,natom,mxhist))
227  ABI_ALLOCATE(hist%strten,(6,mxhist))
228  ABI_ALLOCATE(hist%vel,(3,natom,mxhist))
229  ABI_ALLOCATE(hist%vel_cell,(3,3,mxhist))
230  ABI_ALLOCATE(hist%etot,(mxhist))
231  ABI_ALLOCATE(hist%ekin,(mxhist))
232  ABI_ALLOCATE(hist%entropy,(mxhist))
233  ABI_ALLOCATE(hist%time,(mxhist))
234 
235  hist%etot(1)=zero
236  hist%ekin(1)=zero
237  hist%entropy(1)=zero
238  hist%time(1)=zero
239 
240  hist%acell(:,1)=zero
241  hist%rprimd(:,:,1)=zero
242  hist%xred(:,:,1)=zero
243  hist%fcart(:,:,1)=zero
244  hist%strten(:,1)=zero
245  hist%vel(:,:,1)=zero
246  hist%vel_cell(:,:,1)=zero
247 
248 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

SIDE EFFECTS

PARENTS

CHILDREN

NOTES

SOURCE

278 subroutine abihist_init_1D(hist,natom,mxhist,isVUsed,isARUsed)
279 
280 
281 !This section has been created automatically by the script Abilint (TD).
282 !Do not modify the following lines by hand.
283 #undef ABI_FUNC
284 #define ABI_FUNC 'abihist_init_1D'
285 !End of the abilint section
286 
287  implicit none
288 
289 !Arguments ------------------------------------
290  integer,intent(in) :: natom,mxhist
291  logical,intent(in) :: isVUsed,isARUsed
292  type(abihist),intent(inout) :: hist(:)
293 
294 !Local variables-------------------------------
295  integer :: ii
296 
297 ! ***************************************************************
298 
299  do ii=1,size(hist)
300    call abihist_init_0D(hist(ii),natom,mxhist,isVUsed,isARUsed)
301  end do
302 
303 end subroutine abihist_init_1D

m_abihist/def_file_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 def_file_hist

FUNCTION

INPUTS

OUTPUT

PARENTS

      m_abihist

CHILDREN

SOURCE

1786 subroutine def_file_hist(ncid,natom,nimage,ntypat,npsp,has_nimage)
1787 
1788 
1789 !This section has been created automatically by the script Abilint (TD).
1790 !Do not modify the following lines by hand.
1791 #undef ABI_FUNC
1792 #define ABI_FUNC 'def_file_hist'
1793 !End of the abilint section
1794 
1795 implicit none
1796 
1797 !Arguments ------------------------------------
1798 !scalars
1799  integer,intent(in) :: ncid
1800  integer,intent(in) :: natom,nimage,ntypat,npsp
1801  logical,intent(in) :: has_nimage
1802 
1803 !Local variables-------------------------------
1804 #if defined HAVE_NETCDF
1805 !scalars
1806  integer :: ncerr
1807  integer :: natom_id,nimage_id,ntypat_id,npsp_id,time_id,xyz_id,six_id
1808  integer :: xcart_id,xred_id,fcart_id,fred_id,vel_id,vel_cell_id
1809  integer :: rprimd_id,acell_id,strten_id
1810  integer :: etotal_id,ekin_id,entropy_id,mdtime_id
1811  integer :: typat_id,znucl_id,amu_id,dtion_id,imgmov_id, two_id,mdtemp_id
1812  character(len=500) :: msg
1813 !arrays
1814  integer :: dim0(0),dim1(1),dim2(2),dim3(3),dim4(4)
1815 #endif
1816 
1817 ! *************************************************************************
1818 
1819 #if defined HAVE_NETCDF
1820 
1821 !1.Define the dimensions
1822 
1823  if (npsp/=ntypat) then
1824    msg='HIST file does not support alchemical mixing!'
1825    MSG_WARNING(msg)
1826  end if
1827 
1828  ncerr = nf90_def_dim(ncid,"natom",natom,natom_id)
1829  NCF_CHECK_MSG(ncerr," define dimension natom")
1830 
1831  ncerr = nf90_def_dim(ncid,"ntypat",ntypat,ntypat_id)
1832  NCF_CHECK_MSG(ncerr," define dimension ntypat")
1833 
1834  if (has_nimage) then
1835    ncerr = nf90_def_dim(ncid,"nimage",nimage,nimage_id)
1836    NCF_CHECK_MSG(ncerr," define dimension nimage")
1837  end if
1838 
1839  ncerr = nf90_def_dim(ncid,"npsp",npsp,npsp_id)
1840  NCF_CHECK_MSG(ncerr," define dimension npsp")
1841 
1842  ncerr = nf90_def_dim(ncid,"xyz",3,xyz_id)
1843  NCF_CHECK_MSG(ncerr," define dimension xyz")
1844 
1845  ncerr = nf90_def_dim(ncid,"six",6,six_id)
1846  NCF_CHECK_MSG(ncerr," define dimension six")
1847 
1848  ncerr = nf90_def_dim(ncid,"time",NF90_UNLIMITED,time_id)
1849  NCF_CHECK_MSG(ncerr," define dimension time")
1850 
1851  ncerr = nf90_def_dim(ncid,"two",2,two_id)
1852  NCF_CHECK_MSG(ncerr," define dimension two")
1853 
1854 !2.Define the constant variables
1855 
1856  dim1=(/natom_id/)
1857  call ab_define_var(ncid,dim1,typat_id,NF90_DOUBLE,&
1858 &  "typat","types of atoms","dimensionless" )
1859 
1860  dim1=(/npsp_id/)
1861  call ab_define_var(ncid,dim1,znucl_id,NF90_DOUBLE,&
1862 &  "znucl","atomic charges","atomic units" )
1863 
1864  dim1=(/ntypat_id/)
1865  call ab_define_var(ncid,dim1,amu_id,NF90_DOUBLE,&
1866 &  "amu","atomic masses","atomic units" )
1867 
1868  call ab_define_var(ncid,dim0,dtion_id,NF90_DOUBLE,&
1869 &  "dtion","time step","atomic units" )
1870 
1871 !mdtemp
1872  dim1=(/two_id/)
1873  call ab_define_var(ncid,dim1,mdtemp_id,NF90_DOUBLE,&
1874 &  "mdtemp","Molecular Dynamics Thermostat Temperatures","Kelvin" )
1875 
1876 !mdtime
1877  dim1=(/time_id/)
1878  call ab_define_var(ncid,dim1,mdtime_id,NF90_DOUBLE,&
1879 & "mdtime","Molecular Dynamics or Relaxation TIME","hbar/Ha" )
1880 
1881 !3.Define the evolving variables
1882 
1883 !xcart,xred,fcart,fred,vel
1884  if (has_nimage) then
1885    call ab_define_var(ncid,dim0,imgmov_id,NF90_INT,&
1886   &  "imgmov","Image mover","Not relevant" )
1887    dim4=(/xyz_id,natom_id,nimage_id,time_id/)
1888    call ab_define_var(ncid,dim4,xcart_id,NF90_DOUBLE,&
1889 &   "xcart","vectors (X) of atom positions in CARTesian coordinates","bohr" )
1890    call ab_define_var(ncid,dim4,xred_id,NF90_DOUBLE,&
1891 &   "xred","vectors (X) of atom positions in REDuced coordinates","dimensionless" )
1892    call ab_define_var(ncid,dim4,fcart_id,NF90_DOUBLE,&
1893 &   "fcart","atom Forces in CARTesian coordinates","Ha/bohr" )
1894    call ab_define_var(ncid,dim4,fred_id,NF90_DOUBLE,&
1895 &   "fred","atom Forces in REDuced coordinates","dimensionless" )
1896    call ab_define_var(ncid,dim4,vel_id,NF90_DOUBLE,&
1897 &   "vel","VELocities of atoms","bohr*Ha/hbar" )
1898  else
1899    dim3=(/xyz_id,natom_id,time_id/)
1900    call ab_define_var(ncid,dim3,xcart_id,NF90_DOUBLE,&
1901 &   "xcart","vectors (X) of atom positions in CARTesian coordinates","bohr" )
1902    call ab_define_var(ncid,dim3,xred_id,NF90_DOUBLE,&
1903 &   "xred","vectors (X) of atom positions in REDuced coordinates","dimensionless" )
1904    call ab_define_var(ncid,dim3,fcart_id,NF90_DOUBLE,&
1905 &   "fcart","atom Forces in CARTesian coordinates","Ha/bohr" )
1906    call ab_define_var(ncid,dim3,fred_id,NF90_DOUBLE,&
1907 &   "fred","atom Forces in REDuced coordinates","dimensionless" )
1908    call ab_define_var(ncid,dim3,vel_id,NF90_DOUBLE,&
1909 &   "vel","VELocities of atoms","bohr*Ha/hbar" )
1910  end if
1911 
1912 !rprimd,vel_cell
1913  if (has_nimage) then
1914    dim4=(/xyz_id,xyz_id,nimage_id,time_id/)
1915    call ab_define_var(ncid,dim4,rprimd_id,NF90_DOUBLE,&
1916 &   "rprimd","Real space PRIMitive translations, Dimensional","bohr" )
1917    call ab_define_var(ncid,dim4,vel_cell_id,NF90_DOUBLE,&
1918 &   "vel_cell","VELocities of CELl","bohr*Ha/hbar" )
1919  else
1920    dim3=(/xyz_id,xyz_id,time_id/)
1921    call ab_define_var(ncid,dim3,rprimd_id,NF90_DOUBLE,&
1922 &   "rprimd","Real space PRIMitive translations, Dimensional","bohr" )
1923    call ab_define_var(ncid,dim3,vel_cell_id,NF90_DOUBLE,&
1924 &   "vel_cell","VELocities of cell","bohr*Ha/hbar" )
1925  end if
1926 
1927 !acell
1928  if (has_nimage) then
1929    dim3=(/xyz_id,nimage_id,time_id/)
1930    call ab_define_var(ncid,dim3,acell_id,NF90_DOUBLE,&
1931 &   "acell","CELL lattice vector scaling","bohr" )
1932  else
1933    dim2=(/xyz_id,time_id/)
1934    call ab_define_var(ncid,dim2,acell_id,NF90_DOUBLE,&
1935 &   "acell","CELL lattice vector scaling","bohr" )
1936  end if
1937 
1938 !strten
1939  if (has_nimage) then
1940    dim3=(/six_id,nimage_id,time_id/)
1941    call ab_define_var(ncid,dim3,strten_id,NF90_DOUBLE,&
1942 &   "strten","STRess tensor","Ha/bohr^3" )
1943  else
1944    dim2=(/six_id,time_id/)
1945    call ab_define_var(ncid,dim2,strten_id,NF90_DOUBLE,&
1946 &   "strten","STRess tensor","Ha/bohr^3" )
1947  end if
1948 
1949 !etotal,ekin,entropy
1950  if (has_nimage) then
1951    dim2=(/nimage_id,time_id/)
1952    call ab_define_var(ncid,dim2,etotal_id,NF90_DOUBLE,&
1953 &   "etotal","TOTAL Energy","Ha" )
1954    call ab_define_var(ncid,dim2,ekin_id,NF90_DOUBLE,&
1955 &   "ekin","Energy KINetic ionic","Ha" )
1956    call ab_define_var(ncid,dim2,entropy_id,NF90_DOUBLE,&
1957 &   "entropy","Entropy","" )
1958  else
1959    dim1=(/time_id/)
1960    call ab_define_var(ncid,dim1,etotal_id,NF90_DOUBLE,&
1961 &   "etotal","TOTAL Energy","Ha" )
1962    call ab_define_var(ncid,dim1,ekin_id,NF90_DOUBLE,&
1963 &   "ekin","Energy KINetic ionic","Ha" )
1964    call ab_define_var(ncid,dim1,entropy_id,NF90_DOUBLE,&
1965 &   "entropy","Entropy","" )
1966  end if
1967 
1968 !4.End define mode
1969 
1970  ncerr = nf90_enddef(ncid)
1971  NCF_CHECK_MSG(ncerr," end define mode")
1972 
1973 #endif
1974 
1975 end subroutine def_file_hist

m_abihist/get_dims_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 get_dims_hist

FUNCTION

INPUTS

OUTPUT

PARENTS

      m_abihist,m_tdep_readwrite

CHILDREN

SOURCE

1998 subroutine get_dims_hist(ncid,natom,ntypat,nimage,time,&
1999 &          natom_id,ntypat_id,nimage_id,time_id,xyz_id,six_id,has_nimage)
2000 
2001 
2002 !This section has been created automatically by the script Abilint (TD).
2003 !Do not modify the following lines by hand.
2004 #undef ABI_FUNC
2005 #define ABI_FUNC 'get_dims_hist'
2006 !End of the abilint section
2007 
2008 implicit none
2009 
2010 !Arguments ------------------------------------
2011 !scalars
2012  integer,intent(in) :: ncid
2013  integer,intent(out) :: natom,nimage,time,ntypat
2014  integer,intent(out) :: natom_id,nimage_id,time_id,xyz_id,six_id, ntypat_id
2015  logical,intent(out) :: has_nimage
2016 
2017 !Local variables-------------------------------
2018 #if defined HAVE_NETCDF
2019 !scalars
2020  integer :: ncerr
2021  character(len=5) :: char_tmp
2022 #endif
2023 
2024 ! *************************************************************************
2025 
2026 #if defined HAVE_NETCDF
2027 !Inquire dimensions IDs
2028 
2029  ncerr = nf90_inq_dimid(ncid,"natom",natom_id)
2030  NCF_CHECK_MSG(ncerr," inquire dimension ID for natom")
2031 
2032  ncerr = nf90_inq_dimid(ncid,"npsp",ntypat_id)
2033  NCF_CHECK_MSG(ncerr," inquire dimension ID for npsp")
2034 
2035  ncerr = nf90_inq_dimid(ncid,"xyz",xyz_id)
2036  NCF_CHECK_MSG(ncerr," inquire dimension ID for xyz")
2037 
2038  ncerr = nf90_inq_dimid(ncid,"time",time_id)
2039  NCF_CHECK_MSG(ncerr," inquire dimension ID for time")
2040 
2041  ncerr = nf90_inq_dimid(ncid,"six",six_id)
2042  NCF_CHECK_MSG(ncerr," inquire dimension ID for six")
2043 
2044  ncerr = nf90_inq_dimid(ncid,"nimage",nimage_id)
2045  has_nimage=(ncerr==nf90_noerr)
2046 
2047 !Inquire dimensions lengths
2048 
2049  if (has_nimage) then
2050    ncerr = nf90_inquire_dimension(ncid,nimage_id,char_tmp,nimage)
2051    has_nimage=(ncerr==nf90_noerr)
2052  end if
2053  if (.not.has_nimage) nimage=1
2054 
2055  ncerr = nf90_inquire_dimension(ncid,natom_id,char_tmp,natom)
2056  NCF_CHECK_MSG(ncerr," inquire dimension natom")
2057 
2058  ncerr = nf90_inquire_dimension(ncid,ntypat_id,char_tmp,ntypat)
2059  NCF_CHECK_MSG(ncerr," inquire dimension ntypat")
2060 
2061  ncerr = nf90_inquire_dimension(ncid,time_id,char_tmp,time)
2062  NCF_CHECK_MSG(ncerr," inquire dimension time")
2063 
2064 #endif
2065 
2066 end subroutine get_dims_hist

m_abihist/get_varid_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 get_varid_hist

FUNCTION

INPUTS

OUTPUT

PARENTS

      m_abihist

CHILDREN

SOURCE

2089 subroutine get_varid_hist(ncid,xcart_id,xred_id,fcart_id,fred_id,vel_id,vel_cell_id,&
2090 &          rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id,has_nimage)
2091 
2092 
2093 !This section has been created automatically by the script Abilint (TD).
2094 !Do not modify the following lines by hand.
2095 #undef ABI_FUNC
2096 #define ABI_FUNC 'get_varid_hist'
2097 !End of the abilint section
2098 
2099 implicit none
2100 
2101 !Arguments ------------------------------------
2102 !scalars
2103  integer,intent(in) :: ncid
2104  integer,intent(out) :: xcart_id,xred_id,fcart_id,fred_id,vel_id
2105  integer,intent(out) :: vel_cell_id,rprimd_id,acell_id,strten_id
2106  integer,intent(out) :: etotal_id,ekin_id,entropy_id,mdtime_id
2107  logical,intent(in)  :: has_nimage
2108 !Local variables-------------------------------
2109 #if defined HAVE_NETCDF
2110 !scalars
2111  integer :: ncerr
2112 #endif
2113 
2114 ! *************************************************************************
2115 
2116 #if defined HAVE_NETCDF
2117 
2118  ncerr = nf90_inq_varid(ncid, "mdtime", mdtime_id)
2119  NCF_CHECK_MSG(ncerr," get the id for mdtime")
2120 
2121  ncerr = nf90_inq_varid(ncid, "xcart", xcart_id)
2122  NCF_CHECK_MSG(ncerr," get the id for xcart")
2123 
2124  ncerr = nf90_inq_varid(ncid, "xred", xred_id)
2125  NCF_CHECK_MSG(ncerr," get the id for xred")
2126 
2127  ncerr = nf90_inq_varid(ncid, "fcart", fcart_id)
2128  NCF_CHECK_MSG(ncerr," get the id for fcart")
2129 
2130  ncerr = nf90_inq_varid(ncid, "fred", fred_id)
2131  NCF_CHECK_MSG(ncerr," get the id for fred")
2132 
2133  ncerr = nf90_inq_varid(ncid, "vel", vel_id)
2134  NCF_CHECK_MSG(ncerr," get the id for vel")
2135 
2136  ncerr = nf90_inq_varid(ncid, "vel_cell", vel_cell_id)
2137  if(has_nimage) then
2138    NCF_CHECK_MSG(ncerr," get the id for vel_cell")
2139  end if
2140 
2141  ncerr = nf90_inq_varid(ncid, "rprimd", rprimd_id)
2142  NCF_CHECK_MSG(ncerr," get the id for rprimd")
2143 
2144  ncerr = nf90_inq_varid(ncid, "acell", acell_id)
2145  NCF_CHECK_MSG(ncerr," get the id for acell")
2146 
2147  ncerr = nf90_inq_varid(ncid, "strten", strten_id)
2148  NCF_CHECK_MSG(ncerr," get the id for strten")
2149 
2150  ncerr = nf90_inq_varid(ncid, "etotal", etotal_id)
2151  NCF_CHECK_MSG(ncerr," get the id for etotal")
2152 
2153  ncerr = nf90_inq_varid(ncid, "ekin", ekin_id)
2154  NCF_CHECK_MSG(ncerr," get the id for ekin")
2155 
2156  ncerr = nf90_inq_varid(ncid, "entropy", entropy_id)
2157  NCF_CHECK_MSG(ncerr," get the id for entropy")
2158 
2159 #endif
2160 
2161 end subroutine get_varid_hist

m_abihist/hist2var [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 hist2var

FUNCTION

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

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)

SIDE EFFECTS

PARENTS

      m_pred_lotf,mover,pred_bfgs,pred_delocint,pred_diisrelax,pred_hmc
      pred_isokinetic,pred_isothermal,pred_langevin,pred_lbfgs,pred_moldyn
      pred_nose,pred_srkna14,pred_steepdesc,pred_velverlet,pred_verlet

CHILDREN

SOURCE

896 subroutine hist2var(acell,hist,natom,rprimd,xred,zDEBUG)
897 
898 
899 !This section has been created automatically by the script Abilint (TD).
900 !Do not modify the following lines by hand.
901 #undef ABI_FUNC
902 #define ABI_FUNC 'hist2var'
903 !End of the abilint section
904 
905  implicit none
906 
907 !Arguments ------------------------------------
908 !scalars
909 integer,intent(in) :: natom
910 type(abihist),intent(in) :: hist
911 logical,intent(in) :: zDEBUG
912 !arrays
913 real(dp),intent(out) :: acell(3)
914 real(dp),intent(out) :: rprimd(3,3)
915 real(dp),intent(out) :: xred(3,natom)
916 
917 !Local variables-------------------------------
918 !scalars
919 integer :: kk
920 
921 ! *************************************************************
922 
923  xred  (:,:)=hist%xred(:,:,hist%ihist)
924  acell (:  )=hist%acell(:,hist%ihist)
925  rprimd(:,:)=hist%rprimd(:,:,hist%ihist)
926 
927  if(zDEBUG)then
928    write (std_out,*) 'Atom positions and cell parameters '
929    write (std_out,*) 'ihist: ',hist%ihist
930    write (std_out,*) 'xred:'
931    do kk=1,natom
932      write (std_out,*) xred(:,kk)
933    end do
934    write(std_out,*) 'rprimd:'
935    do kk=1,3
936      write(std_out,*) rprimd(:,kk)
937    end do
938    write(std_out,*) 'acell:'
939    write(std_out,*) acell(:)
940  end if
941 
942 end subroutine hist2var

m_abihist/read_csts_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 read_csts_hist

FUNCTION

INPUTS

OUTPUT

PARENTS

      m_tdep_readwrite

CHILDREN

SOURCE

2184 subroutine read_csts_hist(ncid,dtion,typat,znucl,amu)
2185 
2186 
2187 !This section has been created automatically by the script Abilint (TD).
2188 !Do not modify the following lines by hand.
2189 #undef ABI_FUNC
2190 #define ABI_FUNC 'read_csts_hist'
2191 !End of the abilint section
2192 
2193 implicit none
2194 
2195 !Arguments ------------------------------------
2196 !scalars
2197  integer,intent(in) :: ncid
2198  real(dp),intent(out) :: dtion
2199 !arrays
2200  integer,intent(out) :: typat(:)
2201  real(dp),intent(out) :: amu(:),znucl(:)
2202 
2203 !Local variables-------------------------------
2204 #if defined HAVE_NETCDF
2205 !scalars
2206  integer :: ncerr
2207  integer :: typat_id,znucl_id,amu_id,dtion_id
2208 #endif
2209 
2210 ! *************************************************************************
2211 
2212 #if defined HAVE_NETCDF
2213 
2214 !1.Get the IDs
2215  ncerr = nf90_inq_varid(ncid, "typat", typat_id)
2216  NCF_CHECK_MSG(ncerr," get the id for typat")
2217 
2218  ncerr = nf90_inq_varid(ncid, "znucl", znucl_id)
2219  NCF_CHECK_MSG(ncerr," get the id for znucl")
2220 
2221  ncerr = nf90_inq_varid(ncid, "amu", amu_id)
2222  NCF_CHECK_MSG(ncerr," get the id for amu")
2223 
2224  ncerr = nf90_inq_varid(ncid, "dtion", dtion_id)
2225  NCF_CHECK_MSG(ncerr," get the id for dtion")
2226 
2227 !2.Write the constants
2228  ncerr = nf90_get_var(ncid, typat_id, typat)
2229  NCF_CHECK_MSG(ncerr," get variable typat")
2230 
2231  ncerr = nf90_get_var(ncid, znucl_id, znucl)
2232  NCF_CHECK_MSG(ncerr," get variable znucl")
2233 
2234  ncerr = nf90_get_var(ncid, amu_id, amu)
2235  NCF_CHECK_MSG(ncerr," get variable amu")
2236 
2237  ncerr = nf90_get_var(ncid, dtion_id, dtion)
2238  NCF_CHECK_MSG(ncerr," get variable dtion")
2239 
2240 #endif
2241 
2242 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,

PARENTS

      m_effective_potential_file,m_tdep_readwrite,mover

CHILDREN

SOURCE

1556 subroutine read_md_hist(filename,hist,isVUsed,isARUsed,readOnlyLast)
1557 
1558 
1559 !This section has been created automatically by the script Abilint (TD).
1560 !Do not modify the following lines by hand.
1561 #undef ABI_FUNC
1562 #define ABI_FUNC 'read_md_hist'
1563 !End of the abilint section
1564 
1565 implicit none
1566 
1567 !Arguments ------------------------------------
1568 !scalars
1569  logical,intent(in) :: isVUsed,isARUsed,readOnlyLast
1570  character(len=*),intent(in) :: filename
1571 !arrays
1572  type(abihist),intent(inout),target :: hist
1573 
1574 !Local variables-------------------------------
1575 #if defined HAVE_NETCDF
1576 !scalars
1577  integer :: ncerr,ncid,nimage,natom,time,start_time, ntypat
1578  integer :: nimage_id,natom_id,xyz_id,time_id,six_id, ntypat_id
1579  integer :: xcart_id,xred_id,fcart_id,fred_id,ekin_id,entropy_id
1580  integer :: mdtime_id,vel_id,vel_cell_id,etotal_id
1581  integer :: acell_id,rprimd_id,strten_id
1582  logical :: has_nimage
1583 #endif
1584 
1585 ! *************************************************************************
1586 
1587 #if defined HAVE_NETCDF
1588 
1589  hist%ihist=0 ; hist%mxhist=0
1590 
1591 !Open netCDF file
1592  ncerr=nf90_open(path=trim(filename),mode=NF90_NOWRITE,ncid=ncid)
1593  if(ncerr /= NF90_NOERR) then
1594    write(std_out,*) 'Could no open ',trim(filename),', starting from scratch'
1595    return
1596  else
1597    write(std_out,*) 'Succesfully open ',trim(filename),' for reading'
1598    write(std_out,*) 'Extracting information from NetCDF file...'
1599  end if
1600 
1601 !Inquire dimensions IDs and lengths
1602  call get_dims_hist(ncid,natom,ntypat,nimage,time,&
1603 &     natom_id,ntypat_id,nimage_id,time_id,xyz_id,six_id,has_nimage)
1604 
1605 !If only the last step is needing (restarxf==-3 for example)
1606  if(readOnlyLast)then
1607    start_time = time
1608    time = 1
1609  else
1610    start_time = 1
1611  end if
1612 
1613 !Allocate hist structure
1614  call abihist_init(hist,natom,time,isVused,isARused)
1615 
1616 !Get the ID of a variables from their name
1617  call get_varid_hist(ncid,xcart_id,xred_id,fcart_id,fred_id,vel_id,vel_cell_id,&
1618 &     rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id,has_nimage)
1619 
1620 !Read variables from the dataset and write them into hist
1621  call read_vars_hist(ncid,hist,natom,time,has_nimage,1,start_time,&
1622 &     xred_id,fcart_id,vel_id,vel_cell_id,rprimd_id,acell_id,&
1623 &     strten_id,etotal_id,ekin_id,entropy_id,mdtime_id)
1624 
1625 !Close NetCDF file
1626  ncerr = nf90_close(ncid)
1627  NCF_CHECK_MSG(ncerr," close netcdf history file")
1628 
1629 #endif
1630 
1631 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

PARENTS

      gstateimg

CHILDREN

SOURCE

1664 subroutine read_md_hist_img(filename,hist,isVUsed,isARused,imgtab)
1665 
1666 
1667 !This section has been created automatically by the script Abilint (TD).
1668 !Do not modify the following lines by hand.
1669 #undef ABI_FUNC
1670 #define ABI_FUNC 'read_md_hist_img'
1671 !End of the abilint section
1672 
1673 implicit none
1674 
1675 !Arguments ------------------------------------
1676 !scalars
1677  logical,intent(in) :: isVUsed,isARused
1678  character(len=*),intent(in) :: filename
1679 !arrays
1680  integer,intent(in),optional :: imgtab(:)
1681  type(abihist),intent(inout),target :: hist(:)
1682 
1683 !Local variables-------------------------------
1684 #if defined HAVE_NETCDF
1685 !scalars
1686  integer :: iimage,iimg,my_nimage,ncerr,ncid,nimage,natom,time
1687  integer :: nimage_id,natom_id,xyz_id,time_id,six_id, ntypat, ntypat_id
1688  integer :: xcart_id,xred_id,fcart_id,fred_id,ekin_id,entropy_id
1689  integer :: mdtime_id,vel_id,vel_cell_id,etotal_id
1690  integer :: acell_id,rprimd_id,strten_id
1691  logical :: has_nimage
1692  character(len=500) :: msg
1693  type(abihist),pointer :: hist_
1694  integer,allocatable :: my_imgtab(:)
1695 #endif
1696 
1697 ! *************************************************************************
1698 
1699 #if defined HAVE_NETCDF
1700 
1701  hist%ihist=0 ; hist%mxhist=0
1702 
1703 !Open netCDF file
1704  ncerr=nf90_open(path=trim(filename),mode=NF90_NOWRITE,ncid=ncid)
1705  if(ncerr /= NF90_NOERR) then
1706    write(std_out,*) 'Could no open ',trim(filename),', starting from scratch'
1707    return
1708  else
1709    write(std_out,*) 'Succesfully open ',trim(filename),' for reading'
1710    write(std_out,*) 'Extracting information from NetCDF file...'
1711  end if
1712 
1713  !Manage multiple images of the cell
1714  my_nimage=size(hist)
1715  if (my_nimage==0) return
1716  ABI_ALLOCATE(my_imgtab,(my_nimage))
1717  if (present(imgtab)) then
1718   if (size(my_imgtab)/=my_nimage) then
1719     msg='Inconsistency between hist and imgtab!'
1720     MSG_BUG(msg)
1721   end if
1722   my_imgtab(:)=imgtab(:)
1723  else
1724    my_imgtab(:)=(/(iimage,iimage=1,my_nimage)/)
1725  end if
1726 
1727 !Inquire dimensions IDs and lengths
1728  call get_dims_hist(ncid,natom,ntypat,nimage,time,&
1729 &     natom_id,ntypat_id,nimage_id,time_id,xyz_id,six_id,has_nimage)
1730 
1731  if (nimage<maxval(my_imgtab)) then
1732    msg='Not enough images in the HIST file!'
1733    MSG_ERROR(msg)
1734  end if
1735 
1736 !Loop over images
1737  do iimage=1,my_nimage
1738    iimg=my_imgtab(iimage)
1739    hist_ => hist(iimage)
1740 
1741 !  Allocate hist structure
1742    call abihist_init(hist_,natom,time,isVused,isARused)
1743 
1744 !  Get the ID of a variables from their name
1745    call get_varid_hist(ncid,xcart_id,xred_id,fcart_id,fred_id,vel_id,vel_cell_id,&
1746 &       rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id,has_nimage)
1747 
1748 !  Read variables from the dataset and write them into hist
1749    call read_vars_hist(ncid,hist_,natom,time,has_nimage,iimg,1,&
1750 &       xred_id,fcart_id,vel_id,vel_cell_id,rprimd_id,acell_id,&
1751 &       strten_id,etotal_id,ekin_id,entropy_id,mdtime_id)
1752 
1753  end do
1754 
1755 !Close NetCDF file
1756  ncerr = nf90_close(ncid)
1757  NCF_CHECK_MSG(ncerr," close netcdf history file")
1758 
1759  ABI_DEALLOCATE(my_imgtab)
1760 
1761 #endif
1762 
1763 end subroutine read_md_hist_img

m_abihist/read_vars_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 read_vars_hist

FUNCTION

INPUTS

OUTPUT

PARENTS

      m_abihist

CHILDREN

SOURCE

2532 subroutine read_vars_hist(ncid,hist,natom,time,has_nimage,iimg,start_time,&
2533 &          xred_id,fcart_id,vel_id,vel_cell_id,rprimd_id,acell_id,&
2534 &          strten_id,etotal_id,ekin_id,entropy_id,mdtime_id)
2535 
2536 
2537 !This section has been created automatically by the script Abilint (TD).
2538 !Do not modify the following lines by hand.
2539 #undef ABI_FUNC
2540 #define ABI_FUNC 'read_vars_hist'
2541 !End of the abilint section
2542 
2543 implicit none
2544 
2545 !Arguments ------------------------------------
2546 !scalars
2547  integer,intent(in) :: ncid,natom,time,iimg
2548  integer,intent(in) :: xred_id,fcart_id,vel_id,vel_cell_id,rprimd_id
2549  integer,intent(in) :: acell_id,strten_id
2550  integer,intent(in) :: etotal_id,ekin_id,entropy_id,mdtime_id
2551  integer,intent(in) :: start_time
2552  logical,intent(in) :: has_nimage
2553  type(abihist),intent(inout),target :: hist
2554 
2555 !Local variables-------------------------------
2556 #if defined HAVE_NETCDF
2557 !scalars
2558  integer :: ncerr
2559 !arrays
2560  integer :: count1(1),count2(2),count3(3),count4(4)
2561  integer :: start1(1),start2(2),start3(3),start4(4)
2562 #endif
2563 
2564 ! *************************************************************************
2565 
2566 #if defined HAVE_NETCDF
2567 
2568 !Variables not depending on imes
2569 
2570 !mdtime
2571  start1=(/start_time/);count1=(/time/)
2572  ncerr = nf90_get_var(ncid,mdtime_id,hist%time(:),count=count1,start=start1)
2573  NCF_CHECK_MSG(ncerr," read variable mdtime")
2574 
2575 !Variables depending on images
2576 
2577 !xred,fcart,vel
2578  if (has_nimage) then
2579    start4=(/1,1,iimg,start_time/);count4=(/3,natom,1,time/)
2580    ncerr = nf90_get_var(ncid,xred_id  ,hist%xred(:,:,:),count=count4,start=start4)
2581    NCF_CHECK_MSG(ncerr," read variable xred")
2582    ncerr = nf90_get_var(ncid,fcart_id ,hist%fcart(:,:,:),count=count4,start=start4)
2583    NCF_CHECK_MSG(ncerr," read variable fcart")
2584    ncerr = nf90_get_var(ncid,vel_id,hist%vel(:,:,:),count=count4,start=start4)
2585    NCF_CHECK_MSG(ncerr," read variable vel")
2586  else
2587    start3=(/1,1,start_time/);count3=(/3,natom,time/)
2588    ncerr = nf90_get_var(ncid,xred_id  ,hist%xred(:,:,:),count=count3,start=start3)
2589    NCF_CHECK_MSG(ncerr," read variable xred")
2590    ncerr = nf90_get_var(ncid,fcart_id ,hist%fcart(:,:,:),count=count3,start=start3)
2591    NCF_CHECK_MSG(ncerr," read variable fcart")
2592    ncerr = nf90_get_var(ncid,vel_id,hist%vel(:,:,:),count=count3,start=start3)
2593    NCF_CHECK_MSG(ncerr," read variable vel")
2594  end if
2595 
2596 !rprimd,vel_cell
2597  if (has_nimage) then
2598    start4=(/1,1,iimg,start_time/);count4=(/3,3,start_time,time/)
2599    ncerr = nf90_get_var(ncid,rprimd_id,hist%rprimd(:,:,:),count=count4,start=start4)
2600    NCF_CHECK_MSG(ncerr," read variable rprimd")
2601    ncerr = nf90_get_var(ncid,vel_cell_id,hist%vel_cell(:,:,:),count=count4,start=start4)
2602    NCF_CHECK_MSG(ncerr," read variable vel_cell")
2603  else
2604    start3=(/1,1,start_time/);count3=(/3,3,time/)
2605    ncerr = nf90_get_var(ncid,rprimd_id,hist%rprimd(:,:,:),count=count3,start=start3)
2606    NCF_CHECK_MSG(ncerr," read variable rprimd")
2607  end if
2608 
2609 !acell
2610  if (has_nimage) then
2611    start3=(/1,iimg,start_time/);count3=(/3,1,time/)
2612    ncerr = nf90_get_var(ncid,acell_id,hist%acell(:,:),count=count3,start=start3)
2613    NCF_CHECK_MSG(ncerr," read variable acell")
2614  else
2615    start2=(/1,start_time/);count2=(/3,time/)
2616    ncerr = nf90_get_var(ncid,acell_id,hist%acell(:,:),count=count2,start=start2)
2617    NCF_CHECK_MSG(ncerr," read variable acell")
2618  end if
2619 
2620 !strten
2621  if (has_nimage) then
2622    start3=(/1,iimg,start_time/);count3=(/6,1,time/)
2623    ncerr = nf90_get_var(ncid, strten_id,hist%strten(:,:),count=count3,start=start3)
2624    NCF_CHECK_MSG(ncerr," read variable strten")
2625  else
2626    start2=(/1,start_time/);count2=(/6,time/)
2627    ncerr = nf90_get_var(ncid, strten_id,hist%strten(:,:),count=count2,start=start2)
2628    NCF_CHECK_MSG(ncerr," read variable strten")
2629  end if
2630 
2631 !etotal,ekin,entropy
2632  if (has_nimage) then
2633    start2=(/1,start_time/);count2=(/1,time/)
2634    ncerr = nf90_get_var(ncid,etotal_id,hist%etot(:),count=count2,start=start2)
2635    NCF_CHECK_MSG(ncerr," read variable etotal")
2636    ncerr = nf90_get_var(ncid,ekin_id ,hist%ekin(:),count=count2,start=start2)
2637    NCF_CHECK_MSG(ncerr," read variable ekin")
2638    ncerr = nf90_get_var(ncid,entropy_id,hist%entropy(:),count=count2,start=start2)
2639    NCF_CHECK_MSG(ncerr," read variable entropy")
2640  else
2641    start1=(/start_time/);count1=(/time/)
2642    ncerr = nf90_get_var(ncid,etotal_id,hist%etot(:),count=count1,start=start1)
2643    NCF_CHECK_MSG(ncerr," read variable etotal")
2644    ncerr = nf90_get_var(ncid,ekin_id,hist%ekin(:),count=count1,start=start1)
2645    NCF_CHECK_MSG(ncerr," read variable ekin")
2646    ncerr = nf90_get_var(ncid,entropy_id,hist%entropy(:),count=count1,start=start1)
2647    NCF_CHECK_MSG(ncerr," read variable entropy")
2648  end if
2649 
2650 #endif
2651 
2652 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,

PARENTS

      gstateimg,m_pred_lotf,mover,pred_bfgs,pred_delocint,pred_diisrelax
      pred_hmc,pred_isokinetic,pred_isothermal,pred_langevin,pred_lbfgs
      pred_moldyn,pred_nose,pred_srkna14,pred_steepdesc,pred_velverlet
      pred_verlet

CHILDREN

SOURCE

744 subroutine var2hist(acell,hist,natom,rprimd,xred,zDEBUG)
745 
746 
747 !This section has been created automatically by the script Abilint (TD).
748 !Do not modify the following lines by hand.
749 #undef ABI_FUNC
750 #define ABI_FUNC 'var2hist'
751 !End of the abilint section
752 
753  implicit none
754 
755 !Arguments ------------------------------------
756 !scalars
757  integer,intent(in) :: natom
758  type(abihist),intent(inout) :: hist
759  logical,intent(in) :: zDEBUG
760 !arrays
761  real(dp),intent(in) :: acell(3)
762  real(dp),intent(in) :: rprimd(3,3)
763  real(dp),intent(in) :: xred(3,natom)
764 
765 !Local variables-------------------------------
766 !scalars
767  integer :: kk
768 
769 ! *************************************************************
770 
771  hist%xred(:,:,hist%ihist)=xred(:,:)
772  hist%rprimd(:,:,hist%ihist)=rprimd(:,:)
773  hist%acell(:,hist%ihist)=acell(:)
774 
775  if(zDEBUG)then
776    write (std_out,*) 'Atom positions and cell parameters '
777    write (std_out,*) 'ihist: ',hist%ihist
778    write (std_out,*) 'xred:'
779    do kk=1,natom
780      write (std_out,*) xred(:,kk)
781    end do
782    write(std_out,*) 'rprimd:'
783    do kk=1,3
784      write(std_out,*) rprimd(:,kk)
785    end do
786    write(std_out,*) 'acell:'
787    write(std_out,*) acell(:)
788  end if
789 
790 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,

PARENTS

      gstateimg,mover

CHILDREN

SOURCE

 974 subroutine vel2hist(amass,hist,vel,vel_cell)
 975 
 976 
 977 !This section has been created automatically by the script Abilint (TD).
 978 !Do not modify the following lines by hand.
 979 #undef ABI_FUNC
 980 #define ABI_FUNC 'vel2hist'
 981 !End of the abilint section
 982 
 983  implicit none
 984 
 985 !Arguments ------------------------------------
 986 !scalars
 987 type(abihist),intent(inout) :: hist
 988 !arrays
 989 real(dp),intent(in) :: amass(:)
 990 real(dp),intent(in) :: vel(:,:)
 991 real(dp),intent(in) :: vel_cell(:,:)
 992 
 993 !Local variables-------------------------------
 994 !scalars
 995 integer :: ii,jj,natom
 996 real(dp) :: ekin
 997 
 998 ! *************************************************************
 999 
1000  natom=size(vel,2)
1001 
1002  if (hist%isVused) then
1003 
1004 !  Store the velocities
1005    hist%vel(:,:,hist%ihist)=vel(:,:)
1006    hist%vel_cell(:,:,hist%ihist)=vel_cell(:,:)
1007 
1008 !  Compute the Ionic Kinetic energy
1009    ekin=zero
1010    do ii=1,natom
1011      do jj=1,3
1012        ekin=ekin+half*amass(ii)*vel(jj,ii)**2
1013      end do
1014    end do
1015 
1016  else
1017    hist%vel(:,:,hist%ihist)=zero
1018    hist%vel_cell(:,:,hist%ihist)=zero
1019    ekin=zero
1020  end if
1021 
1022 !Store the Ionic Kinetic Energy
1023  hist%ekin(hist%ihist)=ekin
1024 
1025 end subroutine vel2hist

m_abihist/write_csts_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 write_csts_hist

FUNCTION

INPUTS

OUTPUT

PARENTS

      m_abihist

CHILDREN

SOURCE

2265 subroutine write_csts_hist(ncid,dtion,imgmov,typat,znucl,amu,mdtemp)
2266 
2267 
2268 !This section has been created automatically by the script Abilint (TD).
2269 !Do not modify the following lines by hand.
2270 #undef ABI_FUNC
2271 #define ABI_FUNC 'write_csts_hist'
2272 !End of the abilint section
2273 
2274 implicit none
2275 
2276 !Arguments ------------------------------------
2277 !scalars
2278  integer,intent(in) :: ncid
2279  real(dp),intent(in) :: dtion
2280  integer,intent(in) :: imgmov
2281 !arrays
2282  integer,intent(in) :: typat(:)
2283  real(dp),intent(in) :: amu(:),znucl(:), mdtemp(2)
2284 
2285 !Local variables-------------------------------
2286 #if defined HAVE_NETCDF
2287 !scalars
2288  integer :: ncerr
2289  integer :: typat_id,znucl_id,amu_id,dtion_id, imgmov_id, mdtemp_id
2290 #endif
2291 
2292 ! *************************************************************************
2293 
2294 #if defined HAVE_NETCDF
2295 
2296 !1.Get the IDs
2297 
2298  ncerr = nf90_inq_varid(ncid, "typat", typat_id)
2299  NCF_CHECK_MSG(ncerr," get the id for typat")
2300 
2301  ncerr = nf90_inq_varid(ncid, "znucl", znucl_id)
2302  NCF_CHECK_MSG(ncerr," get the id for znucl")
2303 
2304  ncerr = nf90_inq_varid(ncid, "amu", amu_id)
2305  NCF_CHECK_MSG(ncerr," get the id for amu")
2306 
2307  ncerr = nf90_inq_varid(ncid, "dtion", dtion_id)
2308  NCF_CHECK_MSG(ncerr," get the id for dtion")
2309 
2310  if ( nf90_noerr == nf90_inq_varid(ncid, "imgmov", imgmov_id) ) then
2311    ncerr = nf90_put_var(ncid, imgmov_id, imgmov)
2312    NCF_CHECK_MSG(ncerr," write variable imgmov")
2313  end if
2314 
2315  if ( nf90_noerr == nf90_inq_varid(ncid, "mdtemp", mdtemp_id) ) then
2316    ncerr = nf90_put_var(ncid, mdtemp_id, mdtemp)
2317    NCF_CHECK_MSG(ncerr," write variable mdtemp")
2318  end if
2319 
2320 !2.Write the constants
2321 
2322  ncerr = nf90_put_var(ncid, typat_id, typat)
2323  NCF_CHECK_MSG(ncerr," write variable typat")
2324 
2325  ncerr = nf90_put_var(ncid, znucl_id, znucl)
2326  NCF_CHECK_MSG(ncerr," write variable znucl")
2327 
2328  ncerr = nf90_put_var(ncid, amu_id, amu)
2329  NCF_CHECK_MSG(ncerr," write variable amu")
2330 
2331  ncerr = nf90_put_var(ncid, dtion_id, dtion)
2332  NCF_CHECK_MSG(ncerr," write variable dtion")
2333 
2334 #endif
2335 
2336 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 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,
  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)

PARENTS

      mover

CHILDREN

SOURCE

1267 subroutine write_md_hist(hist,filename,ifirst,itime,natom,nctime,ntypat,&
1268 &                        typat,amu,znucl,dtion,mdtemp)
1269 
1270 
1271 !This section has been created automatically by the script Abilint (TD).
1272 !Do not modify the following lines by hand.
1273 #undef ABI_FUNC
1274 #define ABI_FUNC 'write_md_hist'
1275 !End of the abilint section
1276 
1277  implicit none
1278 
1279 !Arguments ------------------------------------
1280 !scalars
1281  integer,intent(in) :: ifirst,itime,natom,nctime,ntypat
1282  real(dp),intent(in) :: dtion
1283  character(len=*),intent(in) :: filename
1284 !arrays
1285  integer,intent(in) :: typat(natom)
1286  real(dp),intent(in) :: amu(ntypat),znucl(:),mdtemp(2)
1287  type(abihist),intent(inout),target :: hist
1288 
1289 !Local variables-------------------------------
1290 #if defined HAVE_NETCDF
1291 !scalars
1292  integer :: itime_file,ncerr,ncid,npsp
1293  integer :: xcart_id,xred_id,fcart_id,fred_id
1294  integer :: vel_id,vel_cell_id,etotal_id,acell_id,rprimd_id,strten_id
1295  integer :: ekin_id,entropy_id,mdtime_id
1296  logical :: has_nimage=.false.,need_to_write
1297  integer, parameter :: imgmov=0
1298 !arrays
1299 #endif
1300 
1301 ! *************************************************************************
1302 
1303 #if defined HAVE_NETCDF
1304 
1305  need_to_write = .FALSE.
1306  if(nctime==0 .or. ifirst==1 .or. (itime > nctime .and.mod(itime,nctime) == 0)) need_to_write = .TRUE.
1307 !Return if we don't need to write the HIST file at this step
1308  if (.not. need_to_write) return
1309 
1310  if (ifirst==1) then
1311 !##### First access: Create NetCDF file and write defs
1312 
1313    write(std_out,*) 'Write iteration in HIST netCDF file'
1314    npsp=size(znucl)
1315 
1316 !  Create netCDF file
1317    ncerr = nf90_create(path=trim(filename),cmode=NF90_CLOBBER,ncid=ncid)
1318    NCF_CHECK_MSG(ncerr," create netcdf history file")
1319 
1320 !  Define all dims and vars
1321    call def_file_hist(ncid,natom,1,ntypat,npsp,has_nimage)
1322 
1323 !  Write variables that do not change
1324 !  (they are not read in a hist structure).
1325    call write_csts_hist(ncid,dtion,imgmov,typat,znucl,amu,mdtemp)
1326 
1327 !  Compute the itime for the hist file
1328    itime_file = 1
1329  else
1330 !##### itime>2 access: just open NetCDF file
1331 
1332    if(need_to_write) then
1333 
1334      write(std_out,*) 'Write iteration in HIST netCDF file'
1335 
1336 !    Open netCDF file
1337      ncerr = nf90_open(path=trim(filename),mode=NF90_WRITE, ncid=ncid)
1338      NCF_CHECK_MSG(ncerr," open netcdf history file")
1339 
1340 !    Compute the itime for the hist file
1341      itime_file = itime
1342      if(nctime > 0) itime_file = int(anint(real(itime / nctime,sp)))
1343 
1344    end if
1345  endif
1346 
1347  if(need_to_write) then
1348   !##### Write variables into the dataset
1349   !Get the IDs
1350    call get_varid_hist(ncid,xcart_id,xred_id,fcart_id,fred_id,vel_id,vel_cell_id,&
1351 &       rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id,has_nimage)
1352 !Write
1353    call write_vars_hist(ncid,hist,natom,has_nimage,1,itime_file,&
1354 &       xcart_id,xred_id,fcart_id,fred_id,vel_id,vel_cell_id,&
1355 &       rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id)
1356 
1357 !##### Close the file
1358    ncerr = nf90_close(ncid)
1359    NCF_CHECK_MSG(ncerr," close netcdf history file")
1360  end if
1361 #endif
1362 
1363 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)

PARENTS

      gstateimg

CHILDREN

SOURCE

1410 subroutine write_md_hist_img(hist,filename,ifirst,itime,natom,ntypat,&
1411 &                            typat,amu,znucl,dtion,&
1412 &                            nimage,imgmov,mdtemp,comm_img,imgtab) ! optional arguments
1413 
1414 
1415 !This section has been created automatically by the script Abilint (TD).
1416 !Do not modify the following lines by hand.
1417 #undef ABI_FUNC
1418 #define ABI_FUNC 'write_md_hist_img'
1419 !End of the abilint section
1420 
1421  implicit none
1422 
1423 !Arguments ------------------------------------
1424 !scalars
1425  integer,intent(in) :: ifirst,itime,natom,ntypat
1426  integer,intent(in),optional :: nimage,imgmov,comm_img
1427  real(dp),intent(in) :: dtion
1428  character(len=*),intent(in) :: filename
1429 !arrays
1430  integer,intent(in) :: typat(natom)
1431  integer,intent(in),optional :: imgtab(:)
1432  real(dp),intent(in) :: amu(ntypat),znucl(:),mdtemp(2)
1433  type(abihist),intent(inout),target :: hist(:)
1434 
1435 !Local variables-------------------------------
1436 #if defined HAVE_NETCDF
1437 !scalars
1438  integer :: ii,iimage,iimg,me_img,my_comm_img,my_nimage,ncerr
1439  integer :: ncid,nimage_,nproc_img,npsp,imgmov_
1440  integer :: xcart_id,xred_id,fcart_id,fred_id
1441  integer :: vel_id,vel_cell_id,etotal_id
1442  integer :: acell_id,rprimd_id,strten_id
1443  integer :: ekin_id,entropy_id,mdtime_id
1444  logical :: has_nimage, has_imgmov
1445  character(len=500) :: msg
1446  type(abihist),pointer :: hist_
1447 !arrays
1448  integer,allocatable :: my_imgtab(:)
1449 #endif
1450 
1451 ! *************************************************************************
1452 
1453 #if defined HAVE_NETCDF
1454 
1455 !Manage multiple images of the cell
1456  has_nimage=present(nimage)
1457  has_imgmov=present(imgmov)
1458  nimage_=merge(nimage,1,has_nimage)
1459  imgmov_=merge(imgmov,0,has_imgmov)
1460  my_nimage=size(hist) ; if (my_nimage==0) return
1461  my_comm_img=xmpi_comm_self;if(present(comm_img)) my_comm_img=comm_img
1462  nproc_img=xmpi_comm_size(my_comm_img)
1463  me_img=xmpi_comm_rank(my_comm_img)
1464  ABI_ALLOCATE(my_imgtab,(my_nimage))
1465  if (present(imgtab)) then
1466   if (size(my_imgtab)/=my_nimage) then
1467     msg='Inconsistency between hist and imgtab!'
1468     MSG_BUG(msg)
1469   end if
1470   my_imgtab(:)=imgtab(:)
1471  else
1472    my_imgtab(:)=(/(iimage,iimage=1,my_nimage)/)
1473  end if
1474 
1475 !Has to access the HIST file sequentially, proc by proc
1476  do ii=0,nproc_img-1
1477    call xmpi_barrier(my_comm_img)
1478    if (me_img==ii) then
1479 
1480 !    ##### First access: Create NetCDF file and write defs
1481      if (ifirst==1.and.me_img==0) then
1482        npsp=size(znucl)
1483 !      Create netCDF file
1484        ncerr = nf90_create(path=trim(filename),cmode=NF90_CLOBBER,ncid=ncid)
1485        NCF_CHECK_MSG(ncerr," create netcdf history file")
1486 !      Define all dims and vars
1487        call def_file_hist(ncid,natom,nimage_,ntypat,npsp,has_nimage)
1488 !      Write variables that do not change
1489 !      (they are not read in a hist structure).
1490        call write_csts_hist(ncid,dtion,imgmov_,typat,znucl,amu,mdtemp)
1491      end if
1492 
1493 !    ##### itime>2 access: just open NetCDF file
1494      if (ifirst/=1.or.me_img/=0) then
1495 !      Open netCDF file
1496        ncerr = nf90_open(path=trim(filename),mode=NF90_WRITE, ncid=ncid)
1497        NCF_CHECK_MSG(ncerr," open netcdf history file")
1498      end if
1499 
1500      write(std_out,*) 'Write iteration in HIST netCDF file'
1501 
1502 !    ##### Write variables into the dataset (loop over images)
1503 !    Get the IDs
1504      call get_varid_hist(ncid,xcart_id,xred_id,fcart_id,fred_id,vel_id,vel_cell_id,&
1505 &         rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id,has_nimage)
1506 
1507 !    Write
1508      do iimage=1,my_nimage
1509        iimg=my_imgtab(iimage)
1510        hist_ => hist(iimage)
1511        call write_vars_hist(ncid,hist_,natom,has_nimage,iimg,itime,&
1512 &           xcart_id,xred_id,fcart_id,fred_id,vel_id,vel_cell_id,&
1513 &           rprimd_id,acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id)
1514      end do
1515 
1516 !    ##### Close the file
1517      ncerr = nf90_close(ncid)
1518      NCF_CHECK_MSG(ncerr," close netcdf history file")
1519      ABI_DEALLOCATE(my_imgtab)
1520 
1521 !  End loop on MPI processes
1522    end if
1523  end do
1524 
1525 #endif
1526 
1527 end subroutine write_md_hist_img

m_abihist/write_vars_hist [ Functions ]

[ Top ] [ m_abihist ] [ Functions ]

NAME

 write_vars_hist

FUNCTION

INPUTS

OUTPUT

PARENTS

      m_abihist

CHILDREN

SOURCE

2359 subroutine write_vars_hist(ncid,hist,natom,has_nimage,iimg,itime,&
2360 &          xcart_id,xred_id,fcart_id,fred_id,vel_id,vel_cell_id,rprimd_id,&
2361 &          acell_id,strten_id,etotal_id,ekin_id,entropy_id,mdtime_id)
2362 
2363 
2364 !This section has been created automatically by the script Abilint (TD).
2365 !Do not modify the following lines by hand.
2366 #undef ABI_FUNC
2367 #define ABI_FUNC 'write_vars_hist'
2368 !End of the abilint section
2369 
2370 implicit none
2371 
2372 !Arguments ------------------------------------
2373 !scalars
2374  integer,intent(in) :: ncid,natom,iimg,itime
2375  integer,intent(in) :: xcart_id,xred_id,fcart_id,fred_id,vel_id
2376  integer,intent(in) :: vel_cell_id,rprimd_id,acell_id,strten_id
2377  integer,intent(in) :: etotal_id,ekin_id,entropy_id,mdtime_id
2378  logical,intent(in) :: has_nimage
2379  type(abihist),intent(inout),target :: hist
2380 
2381 !Local variables-------------------------------
2382 #if defined HAVE_NETCDF
2383 !scalars
2384  integer :: ncerr
2385 !arrays
2386  integer :: count2(2),count3(3),count4(4)
2387  integer :: start1(1),start2(2),start3(3),start4(4)
2388  real(dp),allocatable :: conv(:,:)
2389  real(dp),pointer :: xred(:,:),fcart(:,:),rprimd(:,:),vel(:,:),vel_cell(:,:)
2390 #endif
2391 
2392 ! *************************************************************************
2393 
2394 #if defined HAVE_NETCDF
2395 
2396  xred     => hist%xred(:,:,hist%ihist)
2397  fcart    => hist%fcart(:,:,hist%ihist)
2398  vel      => hist%vel(:,:,hist%ihist)
2399  vel_cell => hist%vel_cell(:,:,hist%ihist)
2400  rprimd   => hist%rprimd(:,:,hist%ihist)
2401 
2402 !Variables not depending on images
2403 
2404 !mdtime
2405  start1=(/itime/)
2406  ncerr = nf90_put_var(ncid,mdtime_id,hist%time(hist%ihist),start=start1)
2407  NCF_CHECK_MSG(ncerr," write variable mdtime")
2408 
2409 !Variables depending on images
2410 
2411  ABI_ALLOCATE(conv,(3,natom))
2412 
2413 !xcart,xred,fcart,fred,vel
2414  if (has_nimage) then
2415    start4=(/1,1,iimg,itime/);count4=(/3,natom,1,1/)
2416    call xred2xcart(natom,rprimd,conv,xred)
2417    ncerr = nf90_put_var(ncid,xcart_id,conv, start = start4,count = count4)
2418    NCF_CHECK_MSG(ncerr," write variable xcart")
2419    ncerr = nf90_put_var(ncid,xred_id,xred, start = start4,count = count4)
2420    NCF_CHECK_MSG(ncerr," write variable xred")
2421    ncerr = nf90_put_var(ncid,fcart_id,fcart,start = start4,count = count4)
2422    NCF_CHECK_MSG(ncerr," write variable fcart")
2423    call fcart2fred(fcart,conv,rprimd,natom)
2424    ncerr = nf90_put_var(ncid,fred_id,conv, start = start4,count = count4)
2425    NCF_CHECK_MSG(ncerr," write variable fred")
2426    ncerr = nf90_put_var(ncid,vel_id,vel, start = start4,count = count4)
2427    NCF_CHECK_MSG(ncerr," write variable vel")
2428  else
2429    start3=(/1,1,itime/);count3=(/3,natom,1/)
2430    call xred2xcart(natom,rprimd,conv,xred)
2431    ncerr = nf90_put_var(ncid,xcart_id,conv, start = start3,count = count3)
2432    NCF_CHECK_MSG(ncerr," write variable xcart")
2433    ncerr = nf90_put_var(ncid,xred_id,xred, start = start3,count = count3)
2434    NCF_CHECK_MSG(ncerr," write variable xred")
2435    ncerr = nf90_put_var(ncid,fcart_id,fcart,start = start3,count = count3)
2436    NCF_CHECK_MSG(ncerr," write variable fcart")
2437    call fcart2fred(fcart,conv,rprimd,natom)
2438    ncerr = nf90_put_var(ncid,fred_id,conv, start = start3,count = count3)
2439    NCF_CHECK_MSG(ncerr," write variable fred")
2440    ncerr = nf90_put_var(ncid,vel_id,vel, start = start3,count = count3)
2441    NCF_CHECK_MSG(ncerr," write variable vel")
2442  end if
2443 
2444  ABI_DEALLOCATE(conv)
2445 
2446 !rprimd,vel_cell
2447  if (has_nimage) then
2448    start4=(/1,1,iimg,itime/);count4=(/3,3,1,1/)
2449    ncerr = nf90_put_var(ncid,rprimd_id,hist%rprimd(:,:,hist%ihist),&
2450 &                       start = start4,count = count4)
2451    NCF_CHECK_MSG(ncerr," write variable rprimd")
2452    ncerr = nf90_put_var(ncid,vel_cell_id,hist%vel_cell(:,:,hist%ihist),&
2453 &                       start = start4,count = count4)
2454    NCF_CHECK_MSG(ncerr," write variable vel_cell")
2455  else
2456    start3=(/1,1,itime/);count3=(/3,3,1/)
2457    ncerr = nf90_put_var(ncid,rprimd_id,hist%rprimd(:,:,hist%ihist),&
2458 &                       start = start3,count = count3)
2459    NCF_CHECK_MSG(ncerr," write variable rprimd")
2460  end if
2461 
2462 !acell
2463  if (has_nimage) then
2464    start3=(/1,iimg,itime/);count3=(/3,1,1/)
2465    ncerr = nf90_put_var(ncid,acell_id,hist%acell(:,hist%ihist),&
2466 &                       start = start3,count = count3)
2467    NCF_CHECK_MSG(ncerr," write variable acell")
2468  else
2469    start2=(/1,itime/);count2=(/3,1/)
2470    ncerr = nf90_put_var(ncid,acell_id,hist%acell(:,hist%ihist),&
2471 &                       start = start2,count = count2)
2472    NCF_CHECK_MSG(ncerr," write variable acell")
2473  end if
2474 
2475 !strten
2476  if (has_nimage) then
2477    start3=(/1,iimg,itime/);count3=(/6,1,1/)
2478    ncerr = nf90_put_var(ncid,strten_id,hist%strten(:,hist%ihist),&
2479 &                       start = start3,count = count3)
2480    NCF_CHECK_MSG(ncerr," write variable strten")
2481  else
2482    start2=(/1,itime/);count2=(/6,1/)
2483    ncerr = nf90_put_var(ncid,strten_id,hist%strten(:,hist%ihist),&
2484 &                       start = start2,count = count2)
2485    NCF_CHECK_MSG(ncerr," write variable strten")
2486  end if
2487 
2488 !etotal,ekin,entropy
2489  if (has_nimage) then
2490    start2=(/iimg,itime/)
2491    ncerr = nf90_put_var(ncid,etotal_id,hist%etot(hist%ihist),start=start2)
2492    NCF_CHECK_MSG(ncerr," write variable etotal")
2493    ncerr = nf90_put_var(ncid,ekin_id,hist%ekin(hist%ihist),start=start2)
2494    NCF_CHECK_MSG(ncerr," write variable ekin")
2495    ncerr = nf90_put_var(ncid,entropy_id,hist%entropy(hist%ihist),start=start2)
2496    NCF_CHECK_MSG(ncerr," write variable entropy")
2497  else
2498    start1=(/itime/)
2499    ncerr = nf90_put_var(ncid,etotal_id,hist%etot(hist%ihist),start=start1)
2500    NCF_CHECK_MSG(ncerr," write variable etotal")
2501    ncerr = nf90_put_var(ncid,ekin_id,hist%ekin(hist%ihist),start=start1)
2502    NCF_CHECK_MSG(ncerr," write variable ekin")
2503    ncerr = nf90_put_var(ncid,entropy_id,hist%entropy(hist%ihist),start=start1)
2504    NCF_CHECK_MSG(ncerr," write variable entropy")
2505  end if
2506 
2507 #endif
2508 
2509 end subroutine write_vars_hist