TABLE OF CONTENTS
ABINIT/addout [ Functions ]
NAME
addout
FUNCTION
Output density and laplacian (see input variables denout and lapout)
COPYRIGHT
Copyright (C) 2002-2017 ABINIT group (PCasek,FF,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 . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
aim_dtset=the structured entity containing all input variables also, uses the variables saved in the module "defs_aimprom"
OUTPUT
(print) WARNING This file does not follow the ABINIT coding rules (yet) : the use of a module to transfer data should be avoided
PARENTS
drvaim
CHILDREN
vgh_rho
SOURCE
35 #if defined HAVE_CONFIG_H 36 #include "config.h" 37 #endif 38 39 #include "abi_common.h" 40 41 subroutine addout(aim_dtset) 42 43 use m_profiling_abi 44 45 use defs_basis 46 use defs_aimprom 47 use defs_abitypes 48 49 !This section has been created automatically by the script Abilint (TD). 50 !Do not modify the following lines by hand. 51 #undef ABI_FUNC 52 #define ABI_FUNC 'addout' 53 use interfaces_63_bader, except_this_one => addout 54 !End of the abilint section 55 56 implicit none 57 58 !Arguments ------------------------------------ 59 !scalars 60 type(aim_dataset_type),intent(in) :: aim_dtset 61 62 !Local variables ------------------------------ 63 !scalars 64 integer :: cod,dims,iat,ii,ipos,jj,nn,tgrd 65 real(dp) :: alfa,rho,rr,xx,yy 66 !arrays 67 real(dp) :: grho(3),hrho(3,3),orig(3),vv(3) 68 real(dp),allocatable :: dfld(:),lfld(:),nr(:),stp(:),uu(:,:) 69 70 !************************************************************************ 71 orig(:)=aim_dtset%vpts(:,1) 72 if (aim_dtset%denout > 0) then 73 dims=aim_dtset%denout 74 elseif (aim_dtset%lapout > 0) then 75 dims=aim_dtset%lapout 76 end if 77 78 select case (aim_dtset%dltyp) 79 case (1) 80 cod=1 81 case (2) 82 cod=2 83 case default 84 cod=0 85 end select 86 87 ABI_ALLOCATE(uu,(3,dims)) 88 ABI_ALLOCATE(nr,(dims)) 89 ABI_ALLOCATE(stp,(dims)) 90 91 write(std_out,*) 'grid:', aim_dtset%ngrid(1:dims) 92 write(std_out,*) 'kod :', cod 93 tgrd=1 94 do ii=1,dims 95 tgrd=tgrd*aim_dtset%ngrid(ii) 96 uu(:,ii)=aim_dtset%vpts(:,ii+1)-aim_dtset%vpts(:,1) 97 nr(ii)=vnorm(uu(:,ii),0) 98 stp(ii)=nr(ii)/(aim_dtset%ngrid(ii)-1) 99 uu(:,ii)=uu(:,ii)/nr(ii) 100 end do 101 write(std_out,*) 'tgrd :', tgrd 102 do ii=1,dims 103 write(std_out,*) 'uu :', uu(1:3,ii) 104 end do 105 106 if (aim_dtset%denout > 0) then 107 ABI_ALLOCATE(dfld,(tgrd+1)) 108 dfld(:)=0._dp 109 end if 110 if (aim_dtset%lapout > 0) then 111 ABI_ALLOCATE(lfld,(tgrd+1)) 112 end if 113 114 select case (dims) 115 case (1) 116 nn=0 117 do ii=0,aim_dtset%ngrid(1)-1 118 nn=nn+1 119 vv(:)=orig(:)+ii*stp(1)*uu(:,1) 120 call vgh_rho(vv,rho,grho,hrho,rr,iat,ipos,cod) 121 if (aim_dtset%denout > 0) dfld(nn)=rho 122 if (aim_dtset%lapout > 0) lfld(nn)=hrho(1,1)+hrho(2,2)+hrho(3,3) 123 end do 124 if (aim_dtset%denout==1) then 125 do ii=0,aim_dtset%ngrid(1)-1 126 xx=ii*stp(1) 127 write(untd,'(2E16.8)') xx, dfld(ii+1) 128 end do 129 end if 130 if (aim_dtset%lapout==1) then 131 do ii=0,aim_dtset%ngrid(1)-1 132 xx=ii*stp(1) 133 write(untl,'(2E16.8)') xx, lfld(ii+1) 134 end do 135 end if 136 case (2) 137 nn=0 138 alfa=dot_product(uu(:,1),uu(:,2)) 139 alfa=acos(alfa) 140 do ii=0,aim_dtset%ngrid(2)-1 141 do jj=0,aim_dtset%ngrid(1)-1 142 nn=nn+1 143 vv(:)=orig(:)+jj*uu(:,2)*stp(2)+ii*stp(1)*uu(:,1) 144 call vgh_rho(vv,rho,grho,hrho,rr,iat,ipos,cod) 145 if (aim_dtset%denout > 0) dfld(nn)=rho 146 if (aim_dtset%lapout > 0) lfld(nn)=hrho(1,1)+hrho(2,2)+hrho(3,3) 147 end do 148 end do 149 write(std_out,*) 'generace hotova', nn 150 nn=0 151 if (aim_dtset%denout==2) then 152 do ii=0,aim_dtset%ngrid(2)-1 153 do jj=0,aim_dtset%ngrid(1)-1 154 nn=nn+1 155 xx=jj*stp(1)+cos(alfa)*ii*stp(2) 156 yy=sin(alfa)*ii*stp(2) 157 write(untd,'(3E16.8)') xx, yy, dfld(nn) 158 end do 159 write(untd,*) ' ' 160 end do 161 end if 162 nn=0 163 if (aim_dtset%lapout==2) then 164 write(std_out,*) 'lezes sem?' 165 do ii=0,aim_dtset%ngrid(2)-1 166 do jj=0,aim_dtset%ngrid(1)-1 167 nn=nn+1 168 xx=jj*stp(1)+cos(alfa)*ii*stp(2) 169 yy=sin(alfa)*ii*stp(2) 170 write(untl,'(3E16.8)') xx, yy, lfld(nn) 171 end do 172 write(untl,*) ' ' 173 end do 174 end if 175 end select 176 ABI_DEALLOCATE(uu) 177 ABI_DEALLOCATE(stp) 178 ABI_DEALLOCATE(nr) 179 if(aim_dtset%denout>0) then 180 ABI_DEALLOCATE(dfld) 181 end if 182 if(aim_dtset%lapout>0) then 183 ABI_DEALLOCATE(lfld) 184 end if 185 186 end subroutine addout