TABLE OF CONTENTS


ABINIT/m_m1geo [ Modules ]

[ Top ] [ Modules ]

NAME

 m_m1geo

FUNCTION

 This module contains definition the m1geo type (="move unique geometry")
 and its related init and destroy routines

COPYRIGHT

 Copyright (C) 2018 ABINIT group (XG)
 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

18 #if defined HAVE_CONFIG_H
19 #include "config.h"
20 #endif
21 
22 #include "abi_common.h"
23 
24 module m_m1geo
25 
26  use defs_basis
27  use defs_abitypes
28  use m_abimover
29  use m_abihist
30 
31  implicit none
32 
33  private
34 
35  public :: m1geo_init
36  public :: m1geo_destroy

m_gstateimg/m1geo_destroy [ Functions ]

[ Top ] [ m_gstateimg ] [ Functions ]

NAME

 m1geo_destroy

FUNCTION

 Destroy the m1geo structure

INPUTS

OUTPUT

SOURCE

225  subroutine m1geo_destroy(m1geo_param)
226 
227 
228 !This section has been created automatically by the script Abilint (TD).
229 !Do not modify the following lines by hand.
230 #undef ABI_FUNC
231 #define ABI_FUNC 'm1geo_destroy'
232 !End of the abilint section
233 
234  implicit none
235 
236 !Arguments ------------------------------------
237 !scalars
238  type(m1geo_type),intent(inout) :: m1geo_param
239 
240 !Local variables
241 !integer
242  integer :: iatom,ionmov,natom,nimage,ntimimage,ntypat
243 
244 !************************************************************************
245 
246  if(allocated(m1geo_param%mixesimgf))then
247    ABI_FREE(m1geo_param%mixesimgf)
248  endif
249 
250  if (allocated(m1geo_param%ab_xfh_1geo%xfhist))then
251    ABI_FREE(m1geo_param%ab_xfh_1geo%xfhist)
252  endif
253 
254  call abimover_destroy(m1geo_param%ab_mover)
255  call delocint_fin(m1geo_param%deloc)
256  call abihist_free(m1geo_param%hist_1geo)
257  call mttk_fin(m1geo_param%mttk_vars)
258 
259  end subroutine m1geo_destroy
260 
261 !----------------------------------------------------------------------
262 
263 end module m_m1geo

m_gstateimg/m1geo_init [ Functions ]

[ Top ] [ m_gstateimg ] [ Functions ]

NAME

 m1geo_init

FUNCTION

 Initializes the m1geo structure

INPUTS

OUTPUT

SOURCE

101  subroutine m1geo_init(dtfil,dtset,m1geo_param)
102 
103 
104 !This section has been created automatically by the script Abilint (TD).
105 !Do not modify the following lines by hand.
106 #undef ABI_FUNC
107 #define ABI_FUNC 'm1geo_init'
108 !End of the abilint section
109 
110  implicit none
111 
112 !Arguments ------------------------------------
113 !scalars
114  type(datafiles_type),target,intent(in) :: dtfil
115  type(dataset_type),target,intent(in) :: dtset
116  type(m1geo_type),intent(inout) :: m1geo_param
117 
118 !Local variables
119 !integer
120  integer :: iatom,ionmov,natom,ncycle,nhisttot,nimage,ntimimage,ntypat
121  real(dp), allocatable :: amu_curr(:)
122 
123 !************************************************************************
124 
125 
126  ionmov=dtset%ionmov
127  natom=dtset%natom
128  nimage=dtset%nimage
129  ntimimage=dtset%ntimimage
130  ntypat=dtset%ntypat
131 
132 !Straight initialization from dtset
133  m1geo_param%dt_chkdilatmx  =dtset%chkdilatmx
134  m1geo_param%dilatmx        =dtset%dilatmx
135  m1geo_param%hmctt          =dtset%hmctt
136  m1geo_param%nctime         =dtset%nctime
137  m1geo_param%nimage         =dtset%nimage
138  m1geo_param%npsp           =dtset%npsp
139  m1geo_param%usewvl         =dtset%usewvl
140 
141  ABI_MALLOC(m1geo_param%mixesimgf,(nimage))
142  m1geo_param%mixesimgf(:)=dtset%mixesimgf(1:nimage)
143 
144 
145 !From dtset, with change or adjustment of meaning
146  m1geo_param%ntime          =dtset%ntimimage  ! Beware ntime vs ntimimage
147  m1geo_param%rprimd_orig    =dtset%rprimd_orig(:,:,1)
148 
149  ABI_MALLOC(amu_curr,(ntypat))
150  amu_curr(:)=dtset%amu_orig(:,1)
151 
152 !From dtfil
153  m1geo_param%filnam_ds4     =dtfil%filnam_ds(4)
154 
155 !###########################################################
156 !Init sub-datastructures
157 
158 !1) ab_mover
159  call abimover_ini(m1geo_param%ab_mover,amu_curr,dtfil,dtset,m1geo_param%specs)
160  ABI_DEALLOCATE(amu_curr)
161 
162 !2) deloc
163  if(ionmov==10 .or. ionmov==11)then
164    call delocint_ini(m1geo_param%deloc)
165  end if
166 
167 !3) hist
168 !In a simple approach, hist_1geo is initialized here.
169 !Likely one should also make an interface with the hist(nimage) from the calling routine, and average it over the images...
170 
171  ncycle=m1geo_param%specs%ncycle
172  if(ionmov==25.and.dtset%hmctt>=0) then
173    ncycle=dtset%hmctt
174    if(dtset%hmcsst>0.and.dtset%optcell/=0) then
175      ncycle=ncycle+dtset%hmcsst
176    endif
177  endif
178  nhisttot=ncycle*ntimimage  ! WARNING, here ntime is used instead of ntimimage
179  if (dtset%nctime>0) nhisttot=nhisttot+1
180 
181 !We just store the needed history step not all of them.
182  if(m1geo_param%specs%nhist/=-1) nhisttot = m1geo_param%specs%nhist ! We don't need to store all the history
183 
184  call abihist_init(m1geo_param%hist_1geo,natom,ntimimage,m1geo_param%specs%isVused,m1geo_param%specs%isARused)
185 
186 !4) ab_xfh
187  m1geo_param%ab_xfh_1geo%mxfh=ntimimage
188  m1geo_param%ab_xfh_1geo%nxfh=0
189  m1geo_param%ab_xfh_1geo%nxfhr=0
190  ABI_MALLOC(m1geo_param%ab_xfh_1geo%xfhist,(3,natom+4,2,ntimimage))
191 
192 !5) mttk
193  if (ionmov==13)then
194    call mttk_ini(m1geo_param%mttk_vars,m1geo_param%ab_mover%nnos)
195  end if
196 
197 !###########################################################
198 !Init other values
199 
200  m1geo_param%ncycle=m1geo_param%specs%ncycle
201  m1geo_param%icycle=0
202  m1geo_param%iexit=0
203  m1geo_param%itime=0
204  m1geo_param%nerr_dilatmx=0
205  m1geo_param%skipcycle=.FALSE.
206 
207  end subroutine m1geo_init

m_m1geo/m1geo [ Types ]

[ Top ] [ m_m1geo ] [ Types ]

NAME

 m1geo type

FUNCTION

 This datatype has the purpose to store all the data taken
 usually from dtset (but not only) needed to move a single geometry
 opbtained by the weighting of different images,
 using the algorithms defined by ionmov.

SOURCE

54 type, public :: m1geo_type
55 ! Scalars
56   integer  :: dt_chkdilatmx  ! chkdilatmx from dtset
57   integer  :: hmctt          ! hmctt from dtset
58   integer  :: icycle         ! cycle index
59   integer  :: iexit          ! if nonzero, will have to exit
60   integer  :: itime          ! time index for the move1geo algorithms
61   integer  :: nctime         !
62   integer  :: ncycle         ! foreseen number of cycles
63   integer  :: nerr_dilatmx   ! number of dilatmx trespass errors, if chkdilatmx
64   integer  :: nimage         ! nimage from dtset
65   integer  :: npsp           ! npsp from dtset
66   integer  :: ntime          ! number of time steps for the move1geo algorithms
67   integer  :: usewvl         ! usewvl from dtset
68   real(dp) :: dilatmx        ! dilatmx from dtset
69   logical  :: skipcycle      ! .TRUE. when the remaining of the cycle has to be skipped. .FALSE. otherwise.
70 ! Arrays
71   real(dp) :: rprimd_orig(3,3)   ! rprimd from dtset
72   real(dp),allocatable :: mixesimgf(:)   ! mixesimgf from dtset 
73 ! Character string
74   character(len=fnlen) :: filnam_ds4
75 ! Datatypes
76   type(abimover) :: ab_mover
77   type(abimover_specs) :: specs
78   type(ab_xfh_type) :: ab_xfh_1geo
79   type(delocint) :: deloc
80   type(abihist) :: hist_1geo
81   type(mttk_type) :: mttk_vars
82 
83  end type m1geo_type
84 
85 contains  !=============================================================