TABLE OF CONTENTS


ABINIT/tools_lotf [ Modules ]

[ Top ] [ Modules ]

NAME

 tools_lotf

FUNCTION

  Contains simple functions for LOTF

COPYRIGHT

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

16 #if defined HAVE_CONFIG_H
17 #include "config.h"
18 #endif
19 
20 #include "abi_common.h"
21 
22 module tools_lotf
23 
24  use defs_basis
25  use m_errors
26 
27  implicit none
28 
29  private
30 
31  public ::             &
32    pinterp,            &
33    pinterp_nolinear,   &
34    dlvsum,             &
35    icf
36 
37 
38 contains

tools_lotf/dlvsum [ Functions ]

[ Top ] [ tools_lotf ] [ Functions ]

NAME

 dlvsum

FUNCTION

INPUTS

PARENTS

      m_lotf

CHILDREN

SOURCE

144  subroutine dlvsum(local_pe,npes,sum_local_pe,ndimve)
145   !------------------------------------------------------------
146   ! (ADV) GLOBAL VECTOR SUM, REAL 28/5/94
147 
148 !This section has been created automatically by the script Abilint (TD).
149 !Do not modify the following lines by hand.
150 #undef ABI_FUNC
151 #define ABI_FUNC 'dlvsum'
152 !End of the abilint section
153 
154   implicit none
155 
156   !Arguments ------------------------
157   integer :: local_pe,npes,ndimve
158   real(dp) :: sum_local_pe(ndimve)
159 
160 ! *************************************************************************
161 
162   ! do nothing here !!! 
163   ! better remove dlvsum  everywhere in this version.
164   return
165  end subroutine dlvsum

tools_lotf/icf [ Functions ]

[ Top ] [ tools_lotf ] [ Functions ]

NAME

 icf

FUNCTION

INPUTS

PARENTS

CHILDREN

SOURCE

181  FUNCTION icf(ix,iy,iz,ic1,ic2,ic3) 
182 
183 
184 !This section has been created automatically by the script Abilint (TD).
185 !Do not modify the following lines by hand.
186 #undef ABI_FUNC
187 #define ABI_FUNC 'icf'
188 !End of the abilint section
189 
190   implicit none
191   !Arguments ------------------------
192   integer :: ix,iy,iz,ic1,ic2,ic3,icf  
193   icf = 1 &
194     + mod( ix - 1 + ic1, ic1 ) &
195     + mod( iy - 1 + ic2, ic2 ) * ic1 &
196     + mod( iz - 1 + ic3, ic3 ) * ic1 * ic2
197  end FUNCTION icf
198 
199 end module tools_lotf

tools_lotf/pinterp [ Functions ]

[ Top ] [ tools_lotf ] [ Functions ]

NAME

 pinterp

FUNCTION

INPUTS

PARENTS

      m_lotf

CHILDREN

SOURCE

54  subroutine  pinterp(a0,a1,ainterpoled,ndim,nitex,n)
55   ! interpolator : a0 for n=0
56   !                a1 for n=nx
57   !                a1 for nx = 0
58 
59 !This section has been created automatically by the script Abilint (TD).
60 !Do not modify the following lines by hand.
61 #undef ABI_FUNC
62 #define ABI_FUNC 'pinterp'
63 !End of the abilint section
64 
65   implicit none
66 
67   !Arguments ------------------------
68   integer,intent(in) :: ndim,nitex,n 
69   real(dp),intent(in) :: a0(ndim),a1(ndim)
70   real(dp),intent(out) :: ainterpoled(ndim)
71   !Local --------------------------- 
72   character(len=500) :: message
73 
74 ! *************************************************************************
75 
76   !--Control
77   if(nitex < 0) then
78     write(message,'(a,i4,a)') '  LOTF: pinterp nitex =',nitex,'smaller than 0'
79     MSG_ERROR(message)
80 
81   elseif(nitex >= 1) then
82     ainterpoled = a0 + (n/real(nitex,dp))*(a1-a0) 
83   elseif(nitex == 0) then
84     ainterpoled = a1 
85   endif
86  end subroutine pinterp

tools_lotf/pinterp_nolinear [ Functions ]

[ Top ] [ tools_lotf ] [ Functions ]

NAME

 pinterp_nolinear

FUNCTION

INPUTS

PARENTS

CHILDREN

SOURCE

102  subroutine  pinterp_nolinear(a0,a1,ainterpoled,ndim,nitex,n)
103   ! interpolator : a0 for n=0
104   !                a1 for n=nx
105   !                a1 for nx = 0
106 
107 !This section has been created automatically by the script Abilint (TD).
108 !Do not modify the following lines by hand.
109 #undef ABI_FUNC
110 #define ABI_FUNC 'pinterp_nolinear'
111 !End of the abilint section
112 
113   implicit none
114 
115   !Arguments ------------------------
116   integer,intent(in) :: ndim,nitex,n 
117   real(dp),intent(in) :: a0(ndim),a1(ndim)
118   real(dp),intent(out) :: ainterpoled(ndim)
119   !Local --------------------------- 
120   real(dp) :: lambda
121 
122 ! *************************************************************************
123 
124   !--Control
125   lambda = n/real(nitex,dp)
126   ainterpoled = a0 - six*(a1-a0)*(lambda**3/three-lambda**2/two) 
127  end subroutine pinterp_nolinear