mirror of
https://github.com/LCPQ/quantum_package
synced 2024-11-19 12:32:34 +01:00
500 lines
12 KiB
Fortran
500 lines
12 KiB
Fortran
program loc_rasorb
|
|
|
|
implicit none
|
|
|
|
BEGIN_DOC
|
|
! This program performs a localization of the active orbitals
|
|
! of a CASSCF wavefunction, reading the orbitals from a RASORB
|
|
! file of molcas.
|
|
! id1=max is the number of MO in a given symmetry.
|
|
END_DOC
|
|
|
|
integer id1,i_atom,shift,shift_h
|
|
|
|
parameter (id1=300)
|
|
|
|
|
|
|
|
character*1 jobz,uplo
|
|
|
|
character*64 file1,file2
|
|
|
|
character*72 string(id1,8),cdum
|
|
|
|
double precision :: cmo(id1,id1,1),cmoref(id1,id1,1),newcmo(id1,id1,1)
|
|
|
|
double precision ::s(id1,id1,1),dum,ddum(id1,id1),ovl(id1,id1)
|
|
|
|
double precision :: w(id1),work(3*id1),t(id1,id1),wi(id1,id1)
|
|
|
|
integer n,i,j,k,l,nmo(8),isym,nsym,idum,nrot(8),irot(id1,8)
|
|
|
|
integer ipiv(id1),info,lwork
|
|
|
|
logical *1 z54
|
|
print*,'passed the first copy'
|
|
|
|
z54=.false.
|
|
|
|
|
|
|
|
!Read the name of the RasOrb file
|
|
|
|
|
|
print*,'Entering in the loc program'
|
|
|
|
! read(5,*) z54
|
|
print*,'before = '
|
|
accu_norm = 0.d0
|
|
do i =1,mo_tot_num
|
|
accu_norm += dabs(mo_overlap(i,i))
|
|
enddo
|
|
print*,'accu_norm = ',accu_norm
|
|
|
|
nsym = 1
|
|
|
|
nmo(1) = mo_tot_num
|
|
|
|
print*,'nmo(1) = ',nmo(1)
|
|
|
|
cmo = 0.d0
|
|
do isym=1,nsym
|
|
|
|
do i=1,nmo(isym)
|
|
|
|
do j = 1, ao_num
|
|
|
|
cmo(j,i,isym) = mo_coef(j,i)
|
|
|
|
enddo
|
|
|
|
enddo
|
|
|
|
enddo
|
|
print*,'passed the first copy'
|
|
|
|
|
|
|
|
do isym=1,nsym
|
|
|
|
do j=1,mo_tot_num
|
|
|
|
do i=1,ao_num
|
|
|
|
newcmo(i,j,isym)=cmo(i,j,isym)
|
|
|
|
enddo
|
|
|
|
enddo
|
|
|
|
enddo
|
|
print*,'passed the copy'
|
|
|
|
|
|
|
|
nrot(1) = 6 ! number of orbitals to be localized
|
|
|
|
|
|
integer :: index_rot(1000,1)
|
|
|
|
|
|
cmoref = 0.d0
|
|
|
|
! Definition of the index of the MO to be rotated
|
|
! irot(2,1) = 21 ! the first mo to be rotated is the 21 th MO
|
|
! irot(3,1) = 22 ! etc....
|
|
! irot(4,1) = 23 !
|
|
! irot(5,1) = 24 !
|
|
! irot(6,1) = 25 !
|
|
! do i = 1,12
|
|
! irot(i,1) = i+6
|
|
! enddo
|
|
irot(1,1) = 5
|
|
irot(2,1) = 6
|
|
irot(3,1) = 7
|
|
irot(4,1) = 8
|
|
irot(5,1) = 9
|
|
irot(6,1) = 10
|
|
do i = 1, nrot(1)
|
|
print*,'irot(i,1) = ',irot(i,1)
|
|
enddo
|
|
pause
|
|
cmoref(4,1,1) = 1.d0 ! 2S function
|
|
cmoref(5,2,1) = 1.d0 ! 2S function
|
|
cmoref(6,3,1) = 1.d0 ! 2S function
|
|
cmoref(19,4,1) = 1.d0 ! 2S function
|
|
cmoref(20,5,1) = 1.d0 ! 2S function
|
|
cmoref(21,6,1) = 1.d0 ! 2S function
|
|
|
|
! you define the guess vectors that you want
|
|
! the new MO to be close to
|
|
! cmore(i,j,1) = < AO_i | guess_vector_MO(j) >
|
|
! i goes from 1 to ao_num
|
|
! j goes from 1 to nrot(1)
|
|
|
|
! Here you must go to the GAMESS output file
|
|
! where the AOs are listed and explicited
|
|
! From the basis of this knowledge you can build your
|
|
! own guess vectors for the MOs
|
|
! The new MOs are provided in output
|
|
! in the same order than the guess MOs
|
|
|
|
! C-C bonds
|
|
! 1-2
|
|
! i_atom = 1
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,1,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,1,1) = 0.18d0 !
|
|
! cmoref(3+shift,1,1) = 0.1d0 !
|
|
|
|
! cmoref(5+shift,1,1) = -0.1d0 ! 2pX function
|
|
! cmoref(6+shift,1,1) = -0.1d0 ! 2pZ function
|
|
|
|
! i_atom = 2
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,1,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,1,1) = 0.18d0 !
|
|
! cmoref(3+shift,1,1) = 0.1d0 !
|
|
|
|
! cmoref(5+shift,1,1) = 0.1d0 ! 2pX function
|
|
! cmoref(6+shift,1,1) = 0.1d0 ! 2pZ function
|
|
|
|
|
|
! ! 1-3
|
|
! i_atom = 1
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,2,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,2,1) = 0.18d0 !
|
|
! cmoref(3+shift,2,1) = 0.1d0 !
|
|
|
|
! cmoref(5+shift,2,1) = 0.1d0 ! 2pX function
|
|
! cmoref(6+shift,2,1) = -0.1d0 ! 2pZ function
|
|
|
|
! i_atom = 3
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,2,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,2,1) = 0.18d0 !
|
|
! cmoref(3+shift,2,1) = 0.1d0 !
|
|
|
|
! cmoref(5+shift,2,1) = -0.1d0 ! 2pX function
|
|
! cmoref(6+shift,2,1) = 0.1d0 ! 2pZ function
|
|
|
|
! ! 4-6
|
|
! i_atom = 4
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,3,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,3,1) = 0.18d0 !
|
|
! cmoref(3+shift,3,1) = 0.1d0 !
|
|
|
|
! cmoref(5+shift,3,1) = 0.1d0 ! 2pX function
|
|
! cmoref(6+shift,3,1) = -0.1d0 ! 2pZ function
|
|
|
|
! i_atom = 6
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,3,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,3,1) = 0.18d0 !
|
|
! cmoref(3+shift,3,1) = 0.1d0 !
|
|
|
|
! cmoref(5+shift,3,1) = -0.1d0 ! 2pX function
|
|
! cmoref(6+shift,3,1) = 0.1d0 ! 2pZ function
|
|
|
|
|
|
! ! 6-5
|
|
! i_atom = 6
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,4,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,4,1) = 0.18d0 !
|
|
! cmoref(3+shift,4,1) = 0.1d0 !
|
|
|
|
! cmoref(5+shift,4,1) = 0.1d0 ! 2pX function
|
|
! cmoref(6+shift,4,1) = 0.1d0 ! 2pZ function
|
|
|
|
! i_atom = 5
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,4,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,4,1) = 0.18d0 !
|
|
! cmoref(3+shift,4,1) = 0.1d0 !
|
|
|
|
! cmoref(5+shift,4,1) = -0.1d0 ! 2pX function
|
|
! cmoref(6+shift,4,1) = -0.1d0 ! 2pZ function
|
|
|
|
|
|
! ! 2-4
|
|
! i_atom = 2
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,5,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,5,1) = 0.18d0 !
|
|
! cmoref(3+shift,5,1) = 0.1d0 !
|
|
|
|
! cmoref(6+shift,5,1) = 0.1d0 ! 2pZ function
|
|
|
|
! i_atom = 4
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,5,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,5,1) = 0.18d0 !
|
|
! cmoref(3+shift,5,1) = 0.1d0 !
|
|
|
|
! cmoref(6+shift,5,1) = -0.1d0 ! 2pZ function
|
|
|
|
|
|
! ! 3-5
|
|
! i_atom = 3
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,6,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,6,1) = 0.18d0 !
|
|
! cmoref(3+shift,6,1) = 0.1d0 !
|
|
|
|
! cmoref(6+shift,6,1) = 0.1d0 ! 2pZ function
|
|
|
|
! i_atom = 5
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,6,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,6,1) = 0.18d0 !
|
|
! cmoref(3+shift,6,1) = 0.1d0 !
|
|
|
|
! cmoref(6+shift,6,1) = -0.1d0 ! 2pZ function
|
|
|
|
! ! C-H bonds
|
|
! ! 2-7
|
|
! i_atom = 2
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,7,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,7,1) = 0.18d0 !
|
|
! cmoref(3+shift,7,1) = 0.1d0 !
|
|
|
|
! cmoref(5+shift,7,1) = -0.1d0 ! 2pX function
|
|
! cmoref(6+shift,7,1) = 0.1d0 ! 2pZ function
|
|
!
|
|
! i_atom = 7
|
|
! shift_h = (6-1) * 15 + (i_atom - 6)*5
|
|
! cmoref(1+shift_h,7,1) = 0.12d0 ! 1S function
|
|
|
|
! ! 4-10
|
|
! i_atom = 4
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,8,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,8,1) = 0.18d0 !
|
|
! cmoref(3+shift,8,1) = 0.1d0 !
|
|
|
|
! cmoref(5+shift,8,1) = -0.1d0 ! 2pX function
|
|
! cmoref(6+shift,8,1) = -0.1d0 ! 2pZ function
|
|
!
|
|
! i_atom = 10
|
|
! shift_h = (6-1) * 15 + (i_atom - 6)*5
|
|
! cmoref(1+shift_h,8,1) = 0.12d0 ! 1S function
|
|
|
|
! ! 5-11
|
|
! i_atom = 5
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,9,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,9,1) = 0.18d0 !
|
|
! cmoref(3+shift,9,1) = 0.1d0 !
|
|
|
|
! cmoref(5+shift,9,1) = 0.1d0 ! 2pX function
|
|
! cmoref(6+shift,9,1) = -0.1d0 ! 2pZ function
|
|
!
|
|
! i_atom = 11
|
|
! shift_h = (6-1) * 15 + (i_atom - 6)*5
|
|
! cmoref(1+shift_h,9,1) = 0.12d0 ! 1S function
|
|
|
|
! ! 3-8
|
|
! i_atom = 3
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,10,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,10,1) = 0.18d0 !
|
|
! cmoref(3+shift,10,1) = 0.1d0 !
|
|
!
|
|
! cmoref(5+shift,10,1) = 0.1d0 ! 2pX function
|
|
! cmoref(6+shift,10,1) = 0.1d0 ! 2pZ function
|
|
!
|
|
! i_atom = 8
|
|
! shift_h = (6-1) * 15 + (i_atom - 6)*5
|
|
! cmoref(1+shift_h,10,1) = 0.12d0 ! 1S function
|
|
|
|
! ! 1-9
|
|
! i_atom = 1
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,11,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,11,1) = 0.18d0 !
|
|
! cmoref(3+shift,11,1) = 0.1d0 !
|
|
!
|
|
! cmoref(6+shift,11,1) = 0.1d0 ! 2pZ function
|
|
|
|
! i_atom = 9
|
|
! shift_h = (6-1) * 15 + (i_atom - 6)*5
|
|
! cmoref(1+shift_h,11,1) = 0.12d0 ! 1S function
|
|
|
|
!
|
|
! ! 6-12
|
|
! i_atom = 6
|
|
! shift = (i_atom -1) * 15
|
|
! cmoref(1+shift,12,1) = -0.012d0 ! 2S function
|
|
! cmoref(2+shift,12,1) = 0.18d0 !
|
|
! cmoref(3+shift,12,1) = 0.1d0 !
|
|
!
|
|
! cmoref(6+shift,12,1) = -0.1d0 ! 2pZ function
|
|
|
|
! i_atom = 12
|
|
! shift_h = (6-1) * 15 + (i_atom - 6)*5
|
|
! cmoref(1+shift_h,12,1) = 0.12d0 ! 1S function
|
|
! cmoref(12,1,1) = 1.d0 !
|
|
|
|
! cmoref(21,2,1) = 1.d0 !
|
|
! cmoref(30,2,1) = 1.d0 !
|
|
|
|
! cmoref(39,3,1) = 1.d0 !
|
|
! cmoref(48,3,1) = 1.d0 !
|
|
|
|
! cmoref(3,4,1) = 1.d0 !
|
|
! cmoref(12,4,1) =-1.d0 !
|
|
|
|
! cmoref(21,5,1) = 1.d0 !
|
|
! cmoref(30,5,1) =-1.d0 !
|
|
|
|
! cmoref(39,6,1) = 1.d0 !
|
|
! cmoref(48,6,1) =-1.d0 !
|
|
|
|
|
|
|
|
|
|
print*,'passed the definition of the referent vectors '
|
|
!Building the S (overlap) matrix in the AO basis.
|
|
|
|
|
|
|
|
do i = 1, ao_num
|
|
do j = 1, ao_num
|
|
s(i,j,1) = ao_overlap(i,j)
|
|
enddo
|
|
enddo
|
|
!Now big loop over symmetry
|
|
|
|
|
|
|
|
do isym=1,nsym
|
|
|
|
if (nrot(isym).eq.0) cycle
|
|
|
|
|
|
|
|
write (6,*)
|
|
|
|
write (6,*)
|
|
|
|
write (6,*)
|
|
|
|
write (6,*) 'WORKING ON SYMMETRY',isym
|
|
|
|
write (6,*)
|
|
|
|
|
|
|
|
|
|
|
|
!Compute the overlap matrix <ref|vec>
|
|
|
|
|
|
|
|
|
|
|
|
! do i=1,nmo(isym)
|
|
do i=1,ao_num
|
|
|
|
do j=1,nrot(isym)
|
|
|
|
ddum(i,j)=0.d0
|
|
|
|
do k=1,ao_num
|
|
|
|
ddum(i,j)=ddum(i,j)+s(i,k,isym)*cmo(k,irot(j,isym),isym)
|
|
|
|
enddo
|
|
|
|
enddo
|
|
|
|
enddo
|
|
|
|
|
|
|
|
do i=1,nrot(isym)
|
|
|
|
do j=1,nrot(isym)
|
|
|
|
ovl(i,j)=0.d0
|
|
|
|
do k=1,ao_num
|
|
! do k=1,mo_tot_num
|
|
|
|
ovl(i,j)=ovl(i,j)+cmoref(k,i,isym)*ddum(k,j)
|
|
|
|
enddo
|
|
|
|
enddo
|
|
|
|
enddo
|
|
|
|
|
|
|
|
call maxovl(nrot(isym),nrot(isym),ovl,t,wi)
|
|
|
|
|
|
|
|
do i=1,nrot(isym)
|
|
do j=1,ao_num
|
|
write (6,*) 'isym,',isym,nrot(isym),nmo(isym)
|
|
newcmo(j,irot(i,isym),isym)=0.d0
|
|
do k=1,nrot(isym)
|
|
newcmo(j,irot(i,isym),isym)=newcmo(j,irot(i,isym),isym) + cmo(j,irot(k,isym),isym)*t(k,i)
|
|
enddo
|
|
enddo
|
|
enddo
|
|
! if(dabs(newcmo(3,19,1) - mo_coef(3,19)) .gt.1.d-10 )then
|
|
! print*,'Something wrong bitch !!'
|
|
! print*,'newcmo(3,19,1) = ',newcmo(3,19,1)
|
|
! print*,'mo_coef(3,19) = ',mo_coef(3,19)
|
|
! stop
|
|
! endif
|
|
|
|
|
|
|
|
enddo !big loop over symmetry
|
|
|
|
10 format (4E20.12)
|
|
|
|
|
|
! Now we copyt the newcmo into the mo_coef
|
|
|
|
mo_coef = 0.d0
|
|
do isym=1,nsym
|
|
do i=1,nmo(isym)
|
|
do j = 1, ao_num
|
|
mo_coef(j,i) = newcmo(j,i,isym)
|
|
enddo
|
|
enddo
|
|
enddo
|
|
! if(dabs(newcmo(3,19,1) - mo_coef(3,19)) .gt.1.d-10 )then
|
|
print*,'mo_coef(3,19)',mo_coef(3,19)
|
|
pause
|
|
|
|
|
|
! we say that it hase been touched, and valid and that everything that
|
|
! depends on mo_coef must not be reprovided
|
|
double precision :: accu_norm
|
|
touch mo_coef
|
|
print*,'after = '
|
|
accu_norm = 0.d0
|
|
do i =1,mo_tot_num
|
|
accu_norm += dabs(mo_overlap(i,i))
|
|
enddo
|
|
print*,'accu_norm = ',accu_norm
|
|
! We call the routine that saves mo_coef in the ezfio format
|
|
call save_mos
|
|
|
|
|
|
|
|
|
|
|
|
stop
|
|
|
|
end
|