TABLE OF CONTENTS


ABINIT/m_frskerker1 [ Modules ]

[ Top ] [ Modules ]

NAME

 m_frskerker1

FUNCTION

 provide the ability to compute the
 penalty function and its first derivative associated
 with some residuals and a real space dielectric function

COPYRIGHT

 Copyright (C) 1998-2018 ABINIT group (DCA, XG, MT)
 This file is distributed under the terms of the
 GNU General Public License, see ~ABINIT/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~ABINIT/Infos/contributors .

NOTES

 this is neither a function nor a subroutine. This is a module
 It is made of two functions and one init subroutine

PARENTS

CHILDREN

SOURCE

28 #if defined HAVE_CONFIG_H
29 #include "config.h"
30 #endif
31 
32 #include "abi_common.h"
33 
34 module m_frskerker1
35 
36   use defs_basis
37   use defs_abitypes
38   use m_abicore
39 
40   use m_spacepar,  only : laplacian
41   use m_numeric_tools, only : dotproduct
42 
43   implicit none
44   !! common variables copied from input
45   integer,save,private                  :: nfft,nspden,ngfft(18)
46   real(dp),save,allocatable,private    :: deltaW(:,:),mat(:,:),g2cart(:)
47   real(dp),save,private                :: gprimd(3,3),dielng
48   type(dataset_type),pointer,save,private  :: dtset_ptr
49   type(MPI_type),save,private,pointer :: mpi_enreg_ptr
50   !! common variables computed
51   logical,save,private :: ok=.false.
52 
53 contains

m_frskerker1/frskerker1__dpf [ Functions ]

[ Top ] [ m_frskerker1 ] [ Functions ]

NAME

 frskerker1__dpf

FUNCTION

 derivative of the penalty function
 actually not the derivative but something allowing minimization
 at constant density
 formula from the work of rackowski,canning and wang
 H*phi - int(phi**2H d3r)phi

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

280 function frskerker1__dpf(nv1,nv2,vrespc)
281 
282 
283 !This section has been created automatically by the script Abilint (TD).
284 !Do not modify the following lines by hand.
285 #undef ABI_FUNC
286 #define ABI_FUNC 'frskerker1__dpf'
287 !End of the abilint section
288 
289  implicit none
290 !Arguments ------------------------------------
291  integer,intent(in) :: nv1,nv2
292  real(dp),intent(in):: vrespc(nv1,nv2)
293  real(dp) ::frskerker1__dpf(nv1,nv2)
294 
295 !Local variables-------------------------------
296  real(dp) :: buffer1(nv1,nv2),buffer2(nv1,nv2)
297 ! *************************************************************************
298 
299   if(ok) then
300    buffer1=vrespc
301    call laplacian(gprimd,mpi_enreg_ptr,nfft,nspden,ngfft,dtset_ptr%paral_kgb,&
302 &    rdfuncr=buffer1,laplacerdfuncr=buffer2,g2cart_in=g2cart)
303    frskerker1__dpf(:,:)= vrespc(:,:)-deltaW-((dielng)**2)*buffer2(:,:)
304   else
305    frskerker1__dpf = zero
306   end if
307 
308 end function frskerker1__dpf

m_frskerker1/frskerker1__end [ Functions ]

[ Top ] [ m_frskerker1 ] [ Functions ]

NAME

 frskerker1__end

FUNCTION

 ending subroutine
 deallocate memory areas

INPUTS

OUTPUT

PARENTS

      prcrskerker1

CHILDREN

      laplacian

SOURCE

140 subroutine frskerker1__end()
141 
142 
143 !This section has been created automatically by the script Abilint (TD).
144 !Do not modify the following lines by hand.
145 #undef ABI_FUNC
146 #define ABI_FUNC 'frskerker1__end'
147 !End of the abilint section
148 
149   implicit none
150 
151 ! *************************************************************************
152   if(ok) then
153 !  ! set ok to false which prevent using the pf and dpf
154    ok = .false.
155 !  ! free memory
156    ABI_DEALLOCATE(deltaW)
157    ABI_DEALLOCATE(mat)
158    ABI_DEALLOCATE(g2cart)
159   end if
160 
161  end subroutine frskerker1__end

m_frskerker1/frskerker1__init [ Functions ]

[ Top ] [ m_frskerker1 ] [ Functions ]

NAME

 frskerker1__init

FUNCTION

 initialisation subroutine
 Copy every variables required for the energy calculation
 Allocate the required memory

INPUTS

OUTPUT

PARENTS

      prcrskerker1

CHILDREN

      laplacian

SOURCE

 77 subroutine frskerker1__init(dtset_in,mpi_enreg_in,nfft_in,ngfft_in,nspden_in,dielng_in,deltaW_in,gprimd_in,mat_in,g2cart_in )
 78 
 79 
 80 !This section has been created automatically by the script Abilint (TD).
 81 !Do not modify the following lines by hand.
 82 #undef ABI_FUNC
 83 #define ABI_FUNC 'frskerker1__init'
 84 !End of the abilint section
 85 
 86  implicit none
 87 
 88 !Arguments ------------------------------------
 89  integer,intent(in)                  :: nfft_in,ngfft_in(18),nspden_in
 90  real(dp),intent(in)                 :: deltaW_in(nfft_in,nspden_in),mat_in(nfft_in,nspden_in),g2cart_in(nfft_in)
 91  real(dp),dimension(3,3),intent(in)  :: gprimd_in
 92  real(dp),intent(in)                 :: dielng_in
 93  type(dataset_type),target,intent(in) :: dtset_in
 94  type(MPI_type),target,intent(in)  :: mpi_enreg_in
 95 
 96 ! *************************************************************************
 97 
 98 ! !allocation and data transfer
 99 ! !Thought it would have been more logical to use the privates intrinsic of the module as
100 ! !input variables it seems that it is not possible...
101   if(.not.ok) then
102    dtset_ptr  => dtset_in
103    mpi_enreg_ptr => mpi_enreg_in
104    nspden=nspden_in
105    ngfft=ngfft_in
106    nfft=nfft_in
107    ABI_ALLOCATE(deltaW,(size(deltaW_in,1),size(deltaW_in,2)))
108    ABI_ALLOCATE(mat,(size(mat_in,1),size(mat_in,2)))
109    ABI_ALLOCATE(g2cart,(size(g2cart_in,1)))
110    deltaW=deltaW_in
111    dielng=dielng_in
112    mat=mat_in
113    gprimd=gprimd_in
114    g2cart=g2cart_in
115    ok = .true.
116   end if
117  end subroutine frskerker1__init

m_frskerker1/frskerker1__newvres [ Functions ]

[ Top ] [ m_frskerker1 ] [ Functions ]

NAME

 frskerker1__newvres

FUNCTION

 affectation subroutine
 do the required renormalisation when providing a new value for
 the density after application of the gradient

INPUTS

OUTPUT

PARENTS

CHILDREN

      laplacian

SOURCE

184 subroutine frskerker1__newvres(nv1,nv2,x, grad, vrespc)
185 
186 
187 !This section has been created automatically by the script Abilint (TD).
188 !Do not modify the following lines by hand.
189 #undef ABI_FUNC
190 #define ABI_FUNC 'frskerker1__newvres'
191 !End of the abilint section
192 
193  implicit none
194 
195 !Arguments ------------------------------------
196  integer,intent(in) :: nv1,nv2
197  real(dp),intent(in):: x
198  real(dp),intent(inout)::grad(nv1,nv2)
199  real(dp),intent(inout)::vrespc(nv1,nv2)
200 
201 ! *************************************************************************
202 
203   grad(:,:)=x*grad(:,:)
204   vrespc(:,:)=vrespc(:,:)+grad(:,:)
205 
206 end subroutine frskerker1__newvres

m_frskerker1/frskerker1__pf [ Functions ]

[ Top ] [ m_frskerker1 ] [ Functions ]

NAME

 frskerker1__pf

FUNCTION

 penalty function associated with the preconditionned residuals

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

226 function frskerker1__pf(nv1,nv2,vrespc)
227 
228 
229 !This section has been created automatically by the script Abilint (TD).
230 !Do not modify the following lines by hand.
231 #undef ABI_FUNC
232 #define ABI_FUNC 'frskerker1__pf'
233 !End of the abilint section
234 
235  implicit none
236 
237 !Arguments ------------------------------------
238  integer,intent(in) :: nv1,nv2
239  real(dp),intent(in)::vrespc(nv1,nv2)
240  real(dp)           ::frskerker1__pf
241 
242 !Local variables-------------------------------
243  real(dp):: buffer1(nv1,nv2),buffer2(nv1,nv2)
244 ! *************************************************************************
245 
246   if(ok) then
247    buffer1=vrespc
248    call laplacian(gprimd,mpi_enreg_ptr,nfft,nspden,ngfft,dtset_ptr%paral_kgb,&
249 &    rdfuncr=buffer1,laplacerdfuncr=buffer2,g2cart_in=g2cart)
250    buffer2(:,:)=(vrespc(:,:)-((dielng)**2)*buffer2(:,:)) * half -deltaW
251    frskerker1__pf=dotproduct(nv1,nv2,vrespc,buffer2) !*half-dotproduct(nv1,nv2,vrespc,deltaW)
252   else
253    frskerker1__pf=zero
254   end if
255 
256  end function frskerker1__pf