TABLE OF CONTENTS
- ABINIT/m_frskerker2
- m_frskerker2/frskerker2__dpf
- m_frskerker2/frskerker2__end
- m_frskerker2/frskerker2__init
- m_frskerker2/frskerker2__newvres2
- m_frskerker2/frskerker2__pf
ABINIT/m_frskerker2 [ Modules ]
NAME
m_frskerker2
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-2024 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
SOURCE
23 #if defined HAVE_CONFIG_H 24 #include "config.h" 25 #endif 26 27 #include "abi_common.h" 28 29 module m_frskerker2 30 31 use defs_basis 32 use m_abicore 33 use m_dtset 34 35 use defs_abitypes, only : MPI_type 36 use m_spacepar, only : laplacian 37 use m_numeric_tools, only : dotproduct 38 39 implicit none 40 41 !! common variables copied from input 42 integer,save,private :: nfft,nspden,ngfft(18) 43 real(dp),save,allocatable,private :: deltaW(:,:),mat(:,:),rdielng(:) 44 real(dp),save,private :: gprimd(3,3) 45 type(dataset_type),pointer,save,private :: dtset_ptr 46 type(MPI_type),pointer,save,private :: mpi_enreg_ptr 47 !! common variables computed 48 logical,save,private :: ok=.false. 49 ! ************************************************************************* 50 51 contains
m_frskerker2/frskerker2__dpf [ Functions ]
[ Top ] [ m_frskerker2 ] [ Functions ]
NAME
frskerker2__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 that is the simple projection of the change on a direction normal to the density changes
INPUTS
OUTPUT
SOURCE
220 function frskerker2__dpf(nv1,nv2,vrespc) 221 222 !Arguments ------------------------------------ 223 integer,intent(in) :: nv1,nv2 224 real(dp),intent(in)::vrespc(nv1,nv2) 225 real(dp) :: frskerker2__dpf(nv1,nv2) 226 227 !Local variables------------------------------- 228 real(dp):: buffer1(nv1,nv2),buffer2(nv1,nv2) 229 integer :: ispden 230 231 ! ************************************************************************* 232 233 if(ok) then 234 buffer1=vrespc 235 call laplacian(gprimd,mpi_enreg_ptr,nfft,nspden,ngfft,rdfuncr=buffer1,laplacerdfuncr=buffer2) 236 do ispden=1,nspden 237 frskerker2__dpf(:,ispden)= vrespc(:,ispden)-deltaW(:,ispden)-((rdielng(:))**2)*buffer2(:,ispden) 238 end do 239 else 240 frskerker2__dpf = zero 241 end if 242 243 end function frskerker2__dpf
m_frskerker2/frskerker2__end [ Functions ]
[ Top ] [ m_frskerker2 ] [ Functions ]
NAME
frskerker2__end
FUNCTION
ending subroutine deallocate memory areas
INPUTS
OUTPUT
SOURCE
116 subroutine frskerker2__end() 117 118 ! ************************************************************************* 119 if(ok) then 120 ! ! set ok to false which prevent using the pf and dpf 121 ok = .false. 122 ! ! free memory 123 ABI_FREE(deltaW) 124 ABI_FREE(mat) 125 ABI_FREE(rdielng) 126 end if 127 128 end subroutine frskerker2__end
m_frskerker2/frskerker2__init [ Functions ]
[ Top ] [ m_frskerker2 ] [ Functions ]
NAME
frskerker2__init
FUNCTION
initialisation subroutine Copy every variables required for the energy calculation Allocate the required memory
INPUTS
OUTPUT
SOURCE
69 subroutine frskerker2__init(dtset_in,mpi_enreg_in,nfft_in,ngfft_in,nspden_in,rdielng_in,deltaW_in,gprimd_in,mat_in ) 70 71 !Arguments ------------------------------------ 72 type(dataset_type),target,intent(in) :: dtset_in 73 integer,intent(in) :: nfft_in,ngfft_in(18),nspden_in 74 real(dp),intent(in) :: deltaW_in(nfft_in,nspden_in),mat_in(nfft_in,nspden_in) 75 real(dp),intent(in) :: rdielng_in(nfft_in) 76 real(dp),intent(in) :: gprimd_in(3,3) 77 type(MPI_type),target,intent(in) :: mpi_enreg_in 78 79 ! ************************************************************************* 80 ! !allocation and data transfer 81 ! !Thought it would have been more logical to use the privates intrinsic of the module as 82 ! !input variables it seems that it is not possible... 83 if(.not.ok) then 84 dtset_ptr => dtset_in 85 mpi_enreg_ptr => mpi_enreg_in 86 nspden=nspden_in 87 ngfft=ngfft_in 88 nfft=nfft_in 89 ABI_MALLOC(deltaW,(size(deltaW_in,1),size(deltaW_in,2))) 90 ABI_MALLOC(mat,(size(mat_in,1),size(mat_in,2))) 91 ABI_MALLOC(rdielng,(size(rdielng_in))) 92 deltaW=deltaW_in 93 rdielng=rdielng_in 94 mat=mat_in 95 gprimd=gprimd_in 96 ok = .true. 97 end if 98 99 end subroutine frskerker2__init
m_frskerker2/frskerker2__newvres2 [ Functions ]
[ Top ] [ m_frskerker2 ] [ Functions ]
NAME
frskerker2__newvres2
FUNCTION
affectation subroutine do the required renormalisation when providing a new value for the density after application of the gradient
INPUTS
OUTPUT
SOURCE
146 subroutine frskerker2__newvres2(nv1,nv2,x, grad, vrespc) 147 148 !Arguments ------------------------------------ 149 integer,intent(in) :: nv1,nv2 150 real(dp),intent(in) :: x 151 real(dp),intent(inout)::grad(nv1,nv2) 152 real(dp),intent(inout)::vrespc(nv1,nv2) 153 154 ! ************************************************************************* 155 grad(:,:)=x*grad(:,:) 156 vrespc(:,:)=vrespc(:,:)+grad(:,:) 157 158 end subroutine frskerker2__newvres2
m_frskerker2/frskerker2__pf [ Functions ]
[ Top ] [ m_frskerker2 ] [ Functions ]
NAME
frskerker2__pf
FUNCTION
penalty function associated with the preconditionned residuals
INPUTS
OUTPUT
SOURCE
174 function frskerker2__pf(nv1,nv2,vrespc) 175 176 !Arguments ------------------------------------ 177 integer,intent(in) :: nv1,nv2 178 real(dp),intent(in) ::vrespc(nv1,nv2) 179 real(dp) ::frskerker2__pf 180 181 !Local variables------------------------------- 182 real(dp) :: buffer1(nv1,nv2),buffer2(nv1,nv2) 183 integer :: ispden 184 ! ************************************************************************* 185 186 if(ok) then 187 buffer1=vrespc 188 call laplacian(gprimd,mpi_enreg_ptr,nfft,nspden,ngfft,rdfuncr=buffer1,laplacerdfuncr=buffer2) 189 do ispden=1,nspden 190 buffer2(:,ispden)=(vrespc(:,ispden)-((rdielng(:))**2)*buffer2(:,ispden)) & 191 & *half - deltaW(:,ispden) 192 end do 193 ! pf_rscgres=dotproduct(vrespc,buffer2)*half-dotproduct(vrespc,deltaW) 194 frskerker2__pf=dotproduct(nv1,nv2,vrespc,buffer2) 195 else 196 frskerker2__pf=zero 197 end if 198 end function frskerker2__pf