PTEROSOR/QuAcK
10
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2024-12-22 12:23:42 +01:00
Code Releases Activity

Fixing up xcDFT

Browse Source
This commit is contained in:
Pierre-Francois Loos 2019-03-13 12:06:10 +01:00
parent 247bc25db5
commit 1d5ea165b1
3 changed files with 0 additions and 239 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

119
src/xcDFT/read_basis.f90
View File

@ -1,119 +0,0 @@
subroutine read_basis(nAt,rAt,nBas,nO,nV,nShell,TotAngMomShell,CenterShell,KShell,DShell,ExpShell)
! Read basis set information
implicit none
include 'parameters.h'
! Input variables
integer,intent(in) :: nAt,nO(nspin)
double precision,intent(in) :: rAt(nAt,ncart)
! Local variables
integer :: nShAt,iAt,iShell
integer :: i,j,k
character :: shelltype
! Output variables
integer,intent(out) :: nShell,nBas
double precision,intent(out) :: CenterShell(maxShell,ncart)
integer,intent(out) :: TotAngMomShell(maxShell),KShell(maxShell)
double precision,intent(out) :: DShell(maxShell,maxK),ExpShell(maxShell,maxK)
integer,intent(out) :: nV(nspin)
!------------------------------------------------------------------------
! Primary basis set information
!------------------------------------------------------------------------
! Open file with basis set specification
open(unit=2,file='input/basis')
! Read basis information
write(*,'(A28)') 'Gaussian basis set'
write(*,'(A28)') '------------------'
nShell = 0
do i=1,nAt
read(2,*) iAt,nShAt
write(*,'(A28,1X,I16)') 'Atom n. ',iAt
write(*,'(A28,1X,I16)') 'number of shells ',nShAt
write(*,'(A28)') '------------------'
! Basis function centers
do j=1,nShAt
nShell = nShell + 1
do k=1,ncart
CenterShell(nShell,k) = rAt(iAt,k)
enddo
! Shell type and contraction degree
read(2,*) shelltype,KShell(nShell)
if(shelltype == "S") then
TotAngMomShell(nShell) = 0
write(*,'(A28,1X,I16)') 's-type shell with K = ',KShell(nShell)
elseif(shelltype == "P") then
TotAngMomShell(nShell) = 1
write(*,'(A28,1X,I16)') 'p-type shell with K = ',KShell(nShell)
elseif(shelltype == "D") then
TotAngMomShell(nShell) = 2
write(*,'(A28,1X,I16)') 'd-type shell with K = ',KShell(nShell)
elseif(shelltype == "F") then
TotAngMomShell(nShell) = 3
write(*,'(A28,1X,I16)') 'f-type shell with K = ',KShell(nShell)
elseif(shelltype == "G") then
TotAngMomShell(nShell) = 4
write(*,'(A28,1X,I16)') 'g-type shell with K = ',KShell(nShell)
elseif(shelltype == "H") then
TotAngMomShell(nShell) = 5
write(*,'(A28,1X,I16)') 'h-type shell with K = ',KShell(nShell)
elseif(shelltype == "I") then
TotAngMomShell(nShell) = 6
write(*,'(A28,1X,I16)') 'i-type shell with K = ',KShell(nShell)
endif
! Read exponents and contraction coefficients
write(*,'(A28,1X,A16,A16)') '','Exponents','Contraction'
do k=1,Kshell(nShell)
read(2,*) ExpShell(nShell,k),DShell(nShell,k)
write(*,'(A28,1X,F16.10,F16.10)') '',ExpShell(nShell,k),DShell(nShell,k)
enddo
enddo
write(*,'(A28)') '------------------'
enddo
! Total number of shells
write(*,'(A28,1X,I16)') 'Number of shells',nShell
write(*,'(A28)') '------------------'
write(*,*)
! Close file with basis set specification
close(unit=2)
! Calculate number of basis functions
nBas = 0
do iShell=1,nShell
nBas = nBas + (TotAngMomShell(iShell)*TotAngMomShell(iShell) + 3*TotAngMomShell(iShell) + 2)/2
enddo
write(*,'(A28)') '------------------'
write(*,'(A28,1X,I16)') 'Number of basis functions',NBas
write(*,'(A28)') '------------------'
write(*,*)
! Number of virtual orbitals
nV(1) = nBas - nO(1)
nV(2) = nBas - nO(2)
end subroutine read_basis

68
src/xcDFT/read_geometry.f90
View File

@ -1,68 +0,0 @@
subroutine read_geometry(nAt,ZNuc,rA,ENuc)
! Read molecular geometry
implicit none
include 'parameters.h'
! Ouput variables
integer,intent(in) :: nAt
! Local variables
integer :: i,j
double precision :: RAB
character(len=2) :: El
integer,external :: element_number
! Ouput variables
double precision,intent(out) :: ZNuc(NAt),rA(nAt,ncart),ENuc
! Open file with geometry specification
open(unit=1,file='input/molecule')
! Read geometry
read(1,*)
read(1,*)
read(1,*)
do i=1,nAt
read(1,*) El,rA(i,1),rA(i,2),rA(i,3)
ZNuc(i) = element_number(El)
enddo
! Compute nuclear repulsion energy
ENuc = 0
do i=1,nAt-1
do j=i+1,nAt
RAB = (rA(i,1)-rA(j,1))**2 + (rA(i,2)-rA(j,2))**2 + (rA(i,3)-rA(j,3))**2
ENuc = ENuc + ZNuc(i)*ZNuc(j)/sqrt(RAB)
enddo
enddo
! Close file with geometry specification
close(unit=1)
! Print geometry
write(*,'(A28)') '------------------'
write(*,'(A28)') 'Molecular geometry'
write(*,'(A28)') '------------------'
do i=1,NAt
write(*,'(A28,1X,I16)') 'Atom n. ',i
write(*,'(A28,1X,F16.10)') 'Z = ',ZNuc(i)
write(*,'(A28,1X,F16.10,F16.10,F16.10)') 'Atom coordinates:',(rA(i,j),j=1,ncart)
enddo
write(*,*)
write(*,'(A28)') '------------------'
write(*,'(A28,1X,F16.10)') 'Nuclear repulsion energy = ',ENuc
write(*,'(A28)') '------------------'
write(*,*)
end subroutine read_geometry

52
src/xcDFT/read_molecule.f90
View File

@ -1,52 +0,0 @@
subroutine read_molecule(nAt,nEl,nO)
! Read number of atoms nAt and number of electrons nEl
implicit none
include 'parameters.h'
! Input variables
integer,intent(out) :: nAt,nEl(nspin),nO(nspin)
! Local variables
integer :: n
! Open file with geometry specification
open(unit=1,file='input/molecule')
! Read number of atoms and number of electrons
read(1,*)
read(1,*) nAt,n,nEl(1),nEl(2)
! Check imput consistency
if(n /= sum(nEl(:))) then
write(*,*) 'number of electrons inconsistent'
stop
endif
nO(:) = nEl(:)
! Print results
write(*,'(A28)') '----------------------'
write(*,'(A28,1X,I16)') 'Number of atoms',nAt
write(*,'(A28)') '----------------------'
write(*,*)
write(*,'(A28)') '----------------------'
write(*,'(A28,1X,I16)') 'Number of spin-up electrons',nEl(1)
write(*,'(A28,1X,I16)') 'Number of spin-down electrons',nEl(2)
write(*,'(A28,1X,I16)') 'Total number of electrons',sum(nEl(:))
write(*,'(A28)') '----------------------'
write(*,*)
! Close file with geometry specification
close(unit=1)
end subroutine read_molecule