TABLE OF CONTENTS


ABINIT/m_per_cond [ Modules ]

[ Top ] [ Modules ]

NAME

  m_per_cond

FUNCTION

  This module contains basic tools for periodic conditions traitement.

COPYRIGHT

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

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

25 #if defined HAVE_CONFIG_H
26 #include "config.h"
27 #endif
28 
29 #include "abi_common.h"
30 
31 MODULE m_per_cond
32 
33  use defs_basis
34  use m_abicore
35 
36  implicit none
37 
38  private
39 
40  public :: per_cond   ! Return an arithmetic progression
41  public :: per_dist   ! Return the distance on a periodic grid
42 
43 
44  interface per_cond
45   module procedure per_cond_re
46   module procedure per_cond_int
47   module procedure per_cond_int3
48  end interface per_cond
49 
50  interface per_dist
51   module procedure per_dist_int
52   module procedure per_dist_int1
53   module procedure per_dist_int3
54  end interface per_dist
55 
56 CONTAINS  !===========================================================

m_per_cond/per_cond_int [ Functions ]

[ Top ] [ m_per_cond ] [ Functions ]

NAME

  per_cond_int

FUNCTION

 Given a 2d-array of integer initial(3,nb), it calulates the values
 of any of the array in the periodic orthogonal discretized
 grid begining in 0 and of lengths dim_grid(3)

INPUTS

  initial(1:3,0:nb-1)=initial point
  dim_grid(1:3)= box lengths
  nb= dimension
  metric=is the factor scale of the box axes

OUTPUT

  per_cond_int= initial in the box

SOURCE

132 function per_cond_int(nb,initial,dim_grid)
133 
134 !Arguments ------------------------------------
135 !scalars
136 
137 !This section has been created automatically by the script Abilint (TD).
138 !Do not modify the following lines by hand.
139 #undef ABI_FUNC
140 #define ABI_FUNC 'per_cond_int'
141 !End of the abilint section
142 
143  integer,intent(in) :: nb
144  integer,intent(in) :: dim_grid(3)
145  integer :: per_cond_int(3,0:nb-1)
146  integer, intent(in) :: initial(3,0:nb-1)
147 
148 !Local variables-------------------------------
149  integer :: ii,jj
150 
151 ! *********************************************************************
152  do jj=0,nb-1
153   do ii=1,3
154    per_cond_int(ii,jj)=modulo(initial(ii,jj),dim_grid(ii))
155   end do
156  end do
157 
158 
159 end function per_cond_int

m_per_cond/per_cond_int3 [ Functions ]

[ Top ] [ m_per_cond ] [ Functions ]

NAME

  per_cond_int3

FUNCTION

 Given a 2d-array of integer initial(3,nb), it calulates the values
 of any of the array in the periodic orthogonal discretized
 grid begining in 0 and of lengths dim_grid(3)

INPUTS

  initial(1:3,0:nb-1)=initial point
  dim_grid(1:3)= box lengths
  nb= dimension
  metric=is the factor scale of the box axes

OUTPUT

  per_cond_int3= initial in the box

SOURCE

183 function per_cond_int3(initial,dim_grid)
184 
185 !Arguments ------------------------------------
186 !scalars
187 
188 !This section has been created automatically by the script Abilint (TD).
189 !Do not modify the following lines by hand.
190 #undef ABI_FUNC
191 #define ABI_FUNC 'per_cond_int3'
192 !End of the abilint section
193 
194  integer,intent(in) :: dim_grid(3)
195  integer :: per_cond_int3(3)
196  integer,intent(in) :: initial(3)
197 
198 
199 !Local variables-------------------------------
200  integer :: ii
201 
202 ! *********************************************************************
203 
204  do ii=1,3
205   per_cond_int3(ii) = modulo(initial(ii),dim_grid(ii))
206  end do
207 
208 end function per_cond_int3

m_per_cond/per_cond_re [ Functions ]

[ Top ] [ m_per_cond ] [ Functions ]

NAME

  per_cond_re

FUNCTION

 Given a 2d-array of integer initial(3,nb), it calulates the values
 of any of the array in the periodic orthogonal discretized
 grid begining in 0 and of lengths dim_grid(3)

INPUTS

  initial(1:3,0:nb-1)=initial point
  dim_grid(1:3)= box lengths
  nb= dimension
  metric=is the factor scale of the box axes

OUTPUT

  per_cond_re= initial in the box

SOURCE

 80 function per_cond_re(nb,initial,dim_grid,metric)
 81 
 82 !Arguments ------------------------------------
 83 !scalars
 84 
 85 !This section has been created automatically by the script Abilint (TD).
 86 !Do not modify the following lines by hand.
 87 #undef ABI_FUNC
 88 #define ABI_FUNC 'per_cond_re'
 89 !End of the abilint section
 90 
 91  integer,intent(in) :: nb
 92  integer,intent(in) :: dim_grid(3)
 93  integer :: per_cond_re(3,0:nb-1)
 94  real(dp),intent(in) :: initial(3,0:nb-1)
 95  real(dp),intent(in) :: metric(1:3)
 96 
 97 !Local variables-------------------------------
 98  integer :: ii,jj
 99 
100 ! *********************************************************************
101  do jj=0,nb-1
102   do ii=1,3
103    per_cond_re(ii,jj)=modulo(nint(initial(ii,jj)/metric(ii)),dim_grid(ii))
104   end do
105  end do
106 
107 
108 end function per_cond_re

m_per_cond/per_dist_int [ Functions ]

[ Top ] [ m_per_cond ] [ Functions ]

NAME

  per_dist_int

FUNCTION

 Given a two 2d-array of integer initA(3,nb),initB(3,nb) in the
 periodic grid, it calulates the values
 of any of the distances in the periodic orthogonal discretized
 grid begining in 0 and of lengths dim_grid(3)

INPUTS

  initA(1:3,0:nb-1)=initial point
  initB(1:3,0:nb-1)=initial point
  dim_grid(1:3)= box lengths
  nb= dimension

OUTPUT

  per_dist_int= initial in the box

SOURCE

233 function per_dist_int(nb,initA,initB,dim_grid)
234 
235 !Arguments ------------------------------------
236 !scalars
237 
238 !This section has been created automatically by the script Abilint (TD).
239 !Do not modify the following lines by hand.
240 #undef ABI_FUNC
241 #define ABI_FUNC 'per_dist_int'
242 !End of the abilint section
243 
244  integer,intent(in) :: nb
245  integer,intent(in) :: dim_grid(3)
246  integer :: per_dist_int(3,0:nb-1)
247  integer,intent(in) :: initA(3,0:nb-1),initB(3,0:nb-1)
248 
249 
250 !Local variables-------------------------------
251  integer :: ii,jj, bo
252 
253 ! *********************************************************************
254  do jj=0,nb-1
255   do ii=1,3
256    bo = abs(initA(ii,jj)-initB(ii,jj))
257    per_dist_int(ii,jj) = minval((/ bo,dim_grid(ii)-bo/))
258   end do
259  end do
260 
261 
262 end function per_dist_int

m_per_cond/per_dist_int1 [ Functions ]

[ Top ] [ m_per_cond ] [ Functions ]

NAME

  per_dist_int1

FUNCTION

 Given a two scalars of integer initA,initB in the
 periodic grid, it calulates the values
 of any of the distances in the periodic orthogonal discretized
 grid begining in 0 and of lengths dim_grid

INPUTS

  initA=initial point
  initB=initial point
  dim_grid= box lengths

OUTPUT

  per_dist_int1= initial in the box

SOURCE

287 function per_dist_int1(initA,initB,dim_grid)
288 
289 !Arguments ------------------------------------
290 !scalars
291 
292 !This section has been created automatically by the script Abilint (TD).
293 !Do not modify the following lines by hand.
294 #undef ABI_FUNC
295 #define ABI_FUNC 'per_dist_int1'
296 !End of the abilint section
297 
298  integer,intent(in) :: dim_grid,initA,initB
299  integer :: per_dist_int1
300 
301 
302 !Local variables-------------------------------
303  integer :: bo
304 
305 ! *********************************************************************
306  bo = abs(initA-initB)
307  per_dist_int1 = minval((/ bo,dim_grid-bo/))
308 
309 end function per_dist_int1

m_per_cond/per_dist_int3 [ Functions ]

[ Top ] [ m_per_cond ] [ Functions ]

NAME

  per_dist_int3

FUNCTION

 Given a two 3d-vector of integer initA(3),initB(3) in the
 periodic grid, it calulates the values
 of any of the distances in the periodic orthogonal discretized
 grid begining in 0 and of lengths dim_grid(3)

INPUTS

  initA(1:3)=initial point
  initB(1:3)=initial point
  dim_grid(1:3)= box lengths

OUTPUT

  per_dist_int3= initial in the box

SOURCE

334 function per_dist_int3(initA,initB,dim_grid)
335 
336 !Arguments ------------------------------------
337 !scalars
338 
339 !This section has been created automatically by the script Abilint (TD).
340 !Do not modify the following lines by hand.
341 #undef ABI_FUNC
342 #define ABI_FUNC 'per_dist_int3'
343 !End of the abilint section
344 
345  integer,intent(in) :: dim_grid(3)
346  integer :: per_dist_int3(3)
347  integer,intent(in) :: initA(3),initB(3)
348 
349 
350 !Local variables-------------------------------
351  integer :: ii, bo
352 
353 ! *********************************************************************
354  do ii=1,3
355   bo = abs(initA(ii)-initB(ii))
356   per_dist_int3(ii) = minval((/ bo,dim_grid(ii)-bo/))
357  end do
358 
359 end function per_dist_int3