mirror of
https://github.com/LCPQ/quantum_package
synced 2024-12-22 20:35:19 +01:00
Thomas Applencourt
commit
3d37e2895d
@ -389,16 +389,17 @@ subroutine H_u_0_mrcc(v_0,u_0,H_jj,n,keys_tmp,Nint,istate)
|
||||
Vt = 0.d0
|
||||
!$OMP DO SCHEDULE(guided)
|
||||
do i=1,n
|
||||
idx(0) = i
|
||||
call filter_connected_davidson(keys_tmp,keys_tmp(1,1,i),Nint,i-1,idx)
|
||||
do jj=1,idx(0)
|
||||
j = idx(jj)
|
||||
if ( (dabs(u_0(j)) > 1.d-7).or.((dabs(u_0(i)) > 1.d-7)) ) then
|
||||
! idx(0) = i
|
||||
! call filter_connected_davidson(keys_tmp,keys_tmp(1,1,i),Nint,i-1,idx)
|
||||
! do jj=1,idx(0)
|
||||
! j = idx(jj)
|
||||
! if ( (dabs(u_0(j)) > 1.d-7).or.((dabs(u_0(i)) > 1.d-7)) ) then
|
||||
do j = 1, i-1
|
||||
call i_H_j(keys_tmp(1,1,j),keys_tmp(1,1,i),Nint,hij)
|
||||
hij = hij + delta_ij(j,i,istate)
|
||||
vt (i) = vt (i) + hij*u_0(j)
|
||||
vt (j) = vt (j) + hij*u_0(i)
|
||||
endif
|
||||
! endif
|
||||
enddo
|
||||
enddo
|
||||
!$OMP END DO
|
||||
|
||||
@ -48,7 +48,7 @@ subroutine run_mrcc
|
||||
E_new = 0.d0
|
||||
delta_E = 1.d0
|
||||
iteration = 0
|
||||
do while (delta_E > 1.d-8)
|
||||
do while (delta_E > 1.d-10)
|
||||
iteration += 1
|
||||
print *, '==========================='
|
||||
print *, 'MRCC Iteration', iteration
|
||||
|
||||
@ -138,7 +138,11 @@ subroutine mrcc_dress(delta_ij_,Ndet,i_generator,n_selected,det_buffer,Nint,ipro
|
||||
do i_state=1,N_states
|
||||
delta_ij_(idx_non_cas(k_sd),idx_cas(i_I),i_state) += dIa_hla(i_state,k_sd)
|
||||
delta_ij_(idx_cas(i_I),idx_non_cas(k_sd),i_state) += dIa_hla(i_state,k_sd)
|
||||
delta_ij_(idx_cas(i_I),idx_cas(i_I),i_state) -= dIa_hla(i_state,k_sd) * ci_inv(i_state) * psi_non_cas_coef(k_sd,i_state)
|
||||
if(dabs(psi_cas_coef(i_I,i_state)).ge.5.d-5)then
|
||||
delta_ij_(idx_cas(i_I),idx_cas(i_I),i_state) -= dIa_hla(i_state,k_sd) * ci_inv(i_state) * psi_non_cas_coef(k_sd,i_state)
|
||||
else
|
||||
delta_ij_(idx_cas(i_I),idx_cas(i_I),i_state) = 0.d0
|
||||
endif
|
||||
enddo
|
||||
enddo
|
||||
call omp_unset_lock( psi_cas_lock(i_I) )
|
||||
|
||||
@ -1,35 +1,28 @@
|
||||
|
||||
BEGIN_PROVIDER [ double precision, lambda_mrcc, (N_states,psi_det_size) ]
|
||||
BEGIN_PROVIDER [ double precision, lambda_mrcc, (N_states,psi_det_size) ]
|
||||
&BEGIN_PROVIDER [ double precision, lambda_pert, (N_states,psi_det_size) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! cm/<Psi_0|H|D_m>
|
||||
! cm/<Psi_0|H|D_m> or perturbative 1/Delta_E(m)
|
||||
END_DOC
|
||||
integer :: i,k
|
||||
double precision :: ihpsi(N_states), hij(N_states)
|
||||
|
||||
double precision :: ihpsi(N_states), hii
|
||||
|
||||
do i=1,N_det_non_cas
|
||||
call i_h_psi(psi_non_cas(1,1,i), psi_cas, psi_cas_coef, N_int, N_det_cas, &
|
||||
size(psi_cas_coef,1), n_states, ihpsi)
|
||||
call i_h_j(psi_non_cas(1,1,i),psi_non_cas(1,1,i),N_int,hij)
|
||||
call i_h_j(psi_non_cas(1,1,i),psi_non_cas(1,1,i),N_int,hii)
|
||||
do k=1,N_states
|
||||
lambda_pert(k,i) = 1d0 / (CI_electronic_energy(k)-hij(k))
|
||||
|
||||
lambda_pert(k,i) = 1.d0 / (psi_cas_energy_diagonalized(k)-hii)
|
||||
lambda_mrcc(k,i) = psi_non_cas_coef(i,k)/ihpsi(k)
|
||||
if ((lambda_mrcc(k,i)/lambda_pert(k,i))<0.d0 .or. (lambda_mrcc(k,i)/lambda_pert(k,i))>4.d0) then
|
||||
lambda_mrcc(k,i) = lambda_pert(k,i)
|
||||
else
|
||||
if ((lambda_mrcc(k,i)/lambda_pert(k,i))<0.1d0 .or. (lambda_mrcc(k,i)/lambda_pert(k,i))>=0d0) then
|
||||
lambda_mrcc(k,i) = lambda_mrcc(k,i)*((cos((lambda_mrcc(k,i)/lambda_pert(k,i))*3.141592653589793d0/0.1d0+3.141592653589793d0)+1d0)/2.d0) &
|
||||
+ lambda_pert(k,i)*(1.d0-((cos((lambda_mrcc(k,i)/lambda_pert(k,i))*3.141592653589793d0/0.1d0+3.141592653589793d0)+1.d0)/2.d0))
|
||||
elseif ((lambda_mrcc(k,i)/lambda_pert(k,i))<=4.0d0 .or. (lambda_mrcc(k,i)/lambda_pert(k,i))>2.0d0) then
|
||||
lambda_mrcc(k,i) = lambda_mrcc(k,i)*(1.d0-(cos(abs(2.d0-(lambda_mrcc(k,i)/lambda_pert(k,i)))*3.141592653589793d0/2.0d0+3.141592653589793d0)+1.d0)/2d0) &
|
||||
+ lambda_pert(k,i)*((cos(abs(2.d0-(lambda_mrcc(k,i)/lambda_pert(k,i)))*3.141592653589793d0/2.0d0+3.141592653589793d0)+1.d0)/2.d0)
|
||||
else
|
||||
lambda_mrcc(k,i) = lambda_mrcc(k,i)
|
||||
endif
|
||||
endif
|
||||
if (dabs(ihpsi(k)).le.1.d-3) then
|
||||
lambda_mrcc(k,i) = 1.d0 / (psi_cas_energy_diagonalized(k)-hii)
|
||||
icount_manu = icount_manu+1
|
||||
cycle
|
||||
endif
|
||||
enddo
|
||||
enddo
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
@ -71,6 +64,16 @@ BEGIN_PROVIDER [ double precision, delta_ij, (N_det,N_det,N_states) ]
|
||||
enddo
|
||||
enddo
|
||||
endif
|
||||
do i = 1, N_det
|
||||
do j = 1, N_det
|
||||
do m = 1, N_states
|
||||
if(isnan(delta_ij(j,i,m)))then
|
||||
delta_ij(j,i,m) = 0.d0
|
||||
endif
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ double precision, h_matrix_dressed, (N_det,N_det) ]
|
||||
|
||||
1
plugins/loc_cele/NEEDED_CHILDREN_MODULES
Normal file
1
plugins/loc_cele/NEEDED_CHILDREN_MODULES
Normal file
@ -0,0 +1 @@
|
||||
AO_Basis Electrons Ezfio_files MO_Basis Nuclei Utils
|
||||
28
plugins/loc_cele/README.rst
Normal file
28
plugins/loc_cele/README.rst
Normal file
@ -0,0 +1,28 @@
|
||||
===============
|
||||
loc_cele Module
|
||||
===============
|
||||
|
||||
Documentation
|
||||
=============
|
||||
|
||||
.. Do not edit this section. It was auto-generated from the
|
||||
.. by the `update_README.py` script.
|
||||
|
||||
`loc_rasorb <http://github.com/LCPQ/quantum_package/tree/master/src/loc_cele/loc_cele.irp.f#L1>`_
|
||||
Undocumented
|
||||
|
||||
Needed Modules
|
||||
==============
|
||||
|
||||
.. Do not edit this section. It was auto-generated from the
|
||||
.. by the `update_README.py` script.
|
||||
|
||||
.. image:: tree_dependency.png
|
||||
|
||||
* `AO_Basis <http://github.com/LCPQ/quantum_package/tree/master/src/AO_Basis>`_
|
||||
* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
|
||||
* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
|
||||
* `MO_Basis <http://github.com/LCPQ/quantum_package/tree/master/src/MO_Basis>`_
|
||||
* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
|
||||
* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
|
||||
|
||||
163
plugins/loc_cele/loc.f
Normal file
163
plugins/loc_cele/loc.f
Normal file
@ -0,0 +1,163 @@
|
||||
c************************************************************************
|
||||
subroutine maxovl(n,m,s,t,w)
|
||||
C
|
||||
C This subprogram contains an iterative procedure to find the
|
||||
C unitary transformation of a set of n vectors which maximizes
|
||||
C the sum of their square overlaps with a set of m reference
|
||||
C vectors (m.le.n)
|
||||
C
|
||||
C S: overlap matrix <ref|vec>
|
||||
C T: rotation matrix
|
||||
C W: new overlap matrix
|
||||
C
|
||||
C
|
||||
implicit real*8(a-h,o-y),logical*1(z)
|
||||
parameter (id1=300)
|
||||
dimension s(id1,id1),t(id1,id1),w(id1,id1)
|
||||
data small/1.d-6/
|
||||
|
||||
zprt=.true.
|
||||
niter=100
|
||||
conv=1.d-8
|
||||
|
||||
write (6,5) n,m,conv
|
||||
5 format (//5x,'Unitary transformation of',i3,' vectors'/
|
||||
* 5x,'following the principle of maximum overlap with a set of',
|
||||
* i3,' reference vectors'/5x,'required convergence on rotation ',
|
||||
* 'angle =',f13.10///5x,'Starting overlap matrix'/)
|
||||
do 6 i=1,m
|
||||
write (6,145) i
|
||||
6 write (6,150) (s(i,j),j=1,n)
|
||||
8 mm=m-1
|
||||
if (m.lt.n) mm=m
|
||||
iter=0
|
||||
do 20 j=1,n
|
||||
do 16 i=1,n
|
||||
t(i,j)=0.d0
|
||||
16 continue
|
||||
do 18 i=1,m
|
||||
18 w(i,j)=s(i,j)
|
||||
20 t(j,j)=1.d0
|
||||
sum=0.d0
|
||||
do 10 i=1,m
|
||||
sum=sum+s(i,i)*s(i,i)
|
||||
10 continue
|
||||
sum=sum/m
|
||||
if (zprt) write (6,12) sum
|
||||
12 format (//5x,'Average square overlap =',f10.6)
|
||||
if (n.eq.1) goto 100
|
||||
last=n
|
||||
j=1
|
||||
21 if (j.ge.last) goto 30
|
||||
sum=0.d0
|
||||
|
||||
do 22 i=1,n
|
||||
22 sum=sum+s(i,j)*s(i,j)
|
||||
if (sum.gt.small) goto 28
|
||||
do 24 i=1,n
|
||||
sij=s(i,j)
|
||||
s(i,j)=-s(i,last)
|
||||
s(i,last)=sij
|
||||
tij=t(i,j)
|
||||
t(i,j)=-t(i,last)
|
||||
t(i,last)=tij
|
||||
24 continue
|
||||
last=last-1
|
||||
goto 21
|
||||
28 j=j+1
|
||||
goto 21
|
||||
30 iter=iter+1
|
||||
imax=0
|
||||
jmax=0
|
||||
dmax=0.d0
|
||||
amax=0.d0
|
||||
do 60 i=1,mm
|
||||
ip=i+1
|
||||
do 50 j=ip,n
|
||||
a=s(i,j)*s(i,j)-s(i,i)*s(i,i)
|
||||
b=-s(i,i)*s(i,j)
|
||||
if (j.gt.m) goto 31
|
||||
a=a+s(j,i)*s(j,i)-s(j,j)*s(j,j)
|
||||
b=b+s(j,i)*s(j,j)
|
||||
31 b=b+b
|
||||
if (a.eq.0.d0) goto 32
|
||||
ba=b/a
|
||||
if (dabs(ba).gt.small) goto 32
|
||||
if (a.gt.0.d0) goto 33
|
||||
tang=-0.5d0*ba
|
||||
cosine=1.d0/dsqrt(1.d0+tang*tang)
|
||||
sine=tang*cosine
|
||||
goto 34
|
||||
32 tang=0.d0
|
||||
if (b.ne.0.d0) tang=(a+dsqrt(a*a+b*b))/b
|
||||
cosine=1.d0/dsqrt(1.d0+tang*tang)
|
||||
sine=tang*cosine
|
||||
goto 34
|
||||
33 cosine=0.d0
|
||||
sine=1.d0
|
||||
34 delta=sine*(a*sine+b*cosine)
|
||||
if (zprt.and.delta.lt.0.d0) write (6,71) i,j,a,b,sine,cosine,delta
|
||||
do 35 k=1,m
|
||||
p=s(k,i)*cosine-s(k,j)*sine
|
||||
q=s(k,i)*sine+s(k,j)*cosine
|
||||
s(k,i)=p
|
||||
35 s(k,j)=q
|
||||
do 40 k=1,n
|
||||
p=t(k,i)*cosine-t(k,j)*sine
|
||||
q=t(k,i)*sine+t(k,j)*cosine
|
||||
t(k,i)=p
|
||||
t(k,j)=q
|
||||
40 continue
|
||||
45 d=dabs(sine)
|
||||
if (d.le.amax) goto 50
|
||||
imax=i
|
||||
jmax=j
|
||||
amax=d
|
||||
dmax=delta
|
||||
50 continue
|
||||
60 continue
|
||||
if (zprt) write (6,70) iter,amax,imax,jmax,dmax
|
||||
70 format (' iter=',i4,' largest rotation=',f12.8,
|
||||
* ', vectors',i3,' and',i3,', incr. of diag. squares=',g12.5)
|
||||
71 format (' i,j,a,b,sin,cos,delta =',2i3,5f10.5)
|
||||
if (amax.lt.conv) goto 100
|
||||
if (iter.lt.niter) goto 30
|
||||
write (6,80)
|
||||
write (6,*) 'niter=',niter
|
||||
80 format (//5x,'*** maximum number of cycles exceeded ',
|
||||
* 'in subroutine maxovl ***'//)
|
||||
stop
|
||||
100 continue
|
||||
do 120 j=1,n
|
||||
if (s(j,j).gt.0.d0) goto 120
|
||||
do 105 i=1,m
|
||||
105 s(i,j)=-s(i,j)
|
||||
do 110 i=1,n
|
||||
110 t(i,j)=-t(i,j)
|
||||
120 continue
|
||||
sum=0.d0
|
||||
do 125 i=1,m
|
||||
125 sum=sum+s(i,i)*s(i,i)
|
||||
sum=sum/m
|
||||
do 122 i=1,m
|
||||
do 122 j=1,n
|
||||
sw=s(i,j)
|
||||
s(i,j)=w(i,j)
|
||||
122 w(i,j)=sw
|
||||
if (.not.zprt) return
|
||||
write (6,12) sum
|
||||
write (6,130)
|
||||
130 format (//5x,'transformation matrix')
|
||||
do 140 i=1,n
|
||||
write (6,145) i
|
||||
140 write (6,150) (t(i,j),j=1,n)
|
||||
145 format (i8)
|
||||
150 format (2x,10f12.8)
|
||||
write (6,160)
|
||||
160 format (//5x,'new overlap matrix'/)
|
||||
do 170 i=1,m
|
||||
write (6,145) i
|
||||
170 write (6,150) (w(i,j),j=1,n)
|
||||
return
|
||||
end
|
||||
|
||||
333
plugins/loc_cele/loc_cele.irp.f
Normal file
333
plugins/loc_cele/loc_cele.irp.f
Normal file
@ -0,0 +1,333 @@
|
||||
program loc_rasorb
|
||||
|
||||
implicit none
|
||||
|
||||
!
|
||||
|
||||
! This program performs a localization of the active orbitals
|
||||
|
||||
! of a CASSCF wavefunction, reading the orbitals from a RASORB
|
||||
|
||||
! file of molcas.
|
||||
|
||||
!
|
||||
|
||||
! id1=max number of MO in a given symmetry.
|
||||
|
||||
!
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
integer id1
|
||||
|
||||
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(1,1) = 20 ! the first mo to be rotated is the 19 th MO
|
||||
irot(2,1) = 21 ! the first mo to be rotated is the 20 th MO
|
||||
irot(3,1) = 22 ! etc....
|
||||
irot(4,1) = 23 !
|
||||
irot(5,1) = 24 !
|
||||
irot(6,1) = 25 !
|
||||
|
||||
! 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
|
||||
cmoref(3,1,1) = 1.d0 !
|
||||
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 isym=1,nsym
|
||||
|
||||
if (nrot(isym).eq.0) cycle
|
||||
|
||||
do i=1,ao_num
|
||||
|
||||
s(i,i,isym)=1.d0
|
||||
|
||||
do j=1,ao_num
|
||||
|
||||
if (i.ne.j) s(i,j,isym)=0.d0
|
||||
|
||||
ddum(i,j)=0.d0
|
||||
|
||||
do k=1,nmo(isym)
|
||||
|
||||
ddum(i,j)=ddum(i,j)+cmo(i,k,isym)*cmo(j,k,isym)
|
||||
|
||||
enddo
|
||||
|
||||
enddo
|
||||
|
||||
enddo
|
||||
|
||||
call dgesv(ao_num,ao_num,ddum,id1,ipiv,s(1,1,isym),id1,info)
|
||||
|
||||
if (info.ne.0) then
|
||||
|
||||
write (6,*) 'Something wrong in dgsev',isym
|
||||
|
||||
stop
|
||||
|
||||
endif
|
||||
|
||||
|
||||
|
||||
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 (4E18.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
|
||||
BIN
plugins/loc_cele/tree_dependency.png
Normal file
BIN
plugins/loc_cele/tree_dependency.png
Normal file
Binary file not shown.
|
After Width: | Height: | Size: 48 KiB |
@ -376,7 +376,7 @@ end
|
||||
! Can be : [ energy | residual | both | wall_time | cpu_time | iterations ]
|
||||
END_DOC
|
||||
davidson_criterion = 'residual'
|
||||
davidson_threshold = 1.d-6
|
||||
davidson_threshold = 1.d-9
|
||||
END_PROVIDER
|
||||
|
||||
subroutine davidson_converged(energy,residual,wall,iterations,cpu,N_st,converged)
|
||||
|
||||
@ -256,6 +256,7 @@ subroutine make_s2_eigenfunction
|
||||
integer :: N_det_new
|
||||
integer, parameter :: bufsze = 1000
|
||||
logical, external :: is_in_wavefunction
|
||||
return
|
||||
|
||||
! !TODO DEBUG
|
||||
! do i=1,N_det
|
||||
|
||||
@ -13,6 +13,9 @@ use bitmasks
|
||||
logical :: good
|
||||
N_det_cas = 0
|
||||
do i=1,N_det
|
||||
do l = 1, N_states
|
||||
psi_cas_coef(i,l) = 0.d0
|
||||
enddo
|
||||
do l=1,n_cas_bitmask
|
||||
good = .True.
|
||||
do k=1,N_int
|
||||
@ -109,6 +112,57 @@ END_PROVIDER
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
BEGIN_PROVIDER [double precision, H_matrix_cas, (N_det_cas,N_det_cas)]
|
||||
implicit none
|
||||
integer :: i,j
|
||||
double precision :: hij
|
||||
do i = 1, N_det_cas
|
||||
do j = 1, N_det_cas
|
||||
call i_H_j(psi_cas(1,1,i),psi_cas(1,1,j),N_int,hij)
|
||||
H_matrix_cas(i,j) = hij
|
||||
enddo
|
||||
enddo
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [double precision, psi_coef_cas_diagonalized, (N_det_cas,N_states)]
|
||||
&BEGIN_PROVIDER [double precision, psi_cas_energy_diagonalized, (N_states)]
|
||||
implicit none
|
||||
integer :: i,j
|
||||
double precision, allocatable :: eigenvectors(:,:), eigenvalues(:)
|
||||
allocate (eigenvectors(size(H_matrix_cas,1),N_det_cas))
|
||||
allocate (eigenvalues(N_det_cas))
|
||||
call lapack_diag(eigenvalues,eigenvectors, &
|
||||
H_matrix_cas,size(H_matrix_cas,1),N_det_cas)
|
||||
do i = 1, N_states
|
||||
psi_cas_energy_diagonalized(i) = eigenvalues(i)
|
||||
do j = 1, N_det_cas
|
||||
psi_coef_cas_diagonalized(j,i) = eigenvectors(j,i)
|
||||
enddo
|
||||
enddo
|
||||
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [double precision, psi_cas_energy, (N_states)]
|
||||
implicit none
|
||||
integer :: i,j,k
|
||||
double precision :: hij,norm,u_dot_v
|
||||
psi_cas_energy = 0.d0
|
||||
|
||||
|
||||
do k = 1, N_states
|
||||
norm = 0.d0
|
||||
do i = 1, N_det_cas
|
||||
norm += psi_cas_coef(i,k) * psi_cas_coef(i,k)
|
||||
do j = 1, N_det_cas
|
||||
psi_cas_energy(k) += psi_cas_coef(i,k) * psi_cas_coef(j,k) * H_matrix_cas(i,j)
|
||||
enddo
|
||||
enddo
|
||||
psi_cas_energy(k) = psi_cas_energy(k) /norm
|
||||
enddo
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
Loading…
Reference in New Issue
Block a user