C C $Id: wrtvec.F,v 1.3 1998/07/16 16:40:57 jjv5 Exp $ C C------------------------------------------------------------------------ SUBROUTINE WRTVECSUBF C C-RDC C WRITES EIGENVECTORS (MOs) TO UNIT IOUT. C for frozen calculation: eigenvectors can be printed at C end of SCF iteration C C By Arjan van der Vaart IMPLICIT DOUBLE PRECISION (A-H,O-Z) #include "divcon.dim" #include "divcon.h" C C C LOCAL: C CHARACTER LABEL(4)*2 DATA LABEL /' S','PX','PY','PZ'/ SAVE LABEL character*1 occ dimension occ(6) logical found C C-RDC WRITE(IOUT,'(///" FINAL EIGENVECTORS:")') c c c MO -> IIII IIII IIII IIII IIII IIII c ENERGY -> -FFFF.FFFF -FFFF.FFFF -FFFF.FFFF -FFFF.FFFF -FFFF.FFFF -FFFF.FFFF C OCC. -> Y Y N N N N c cIIII BR PX -FF.FFFFFF -FF.FFFFFF -FF.FFFFFF -FF.FFFFFF -FF.FFFFFF -FF.FFFFFF ifrpn = 0 ipr = 1 do 1000 ksub=1,nsub C . find out if we need to print the eigenvectors of this subsystem 10 if (ipr.gt.nprsub) return if (iprsub(ipr).gt.ksub) then norbs = iorbpt(ksub+1)-iorbpt(ksub) nk2 = norbs*norbs ifrpn = ifrpn + nk2 goto 1000 elseif (iprsub(ipr).lt.ksub) then ipr = ipr + 1 goto 10 endif C . find the Fermi energy: base occupation on the Fermi energy, C . since we don't know the charge of this subsystem do icr=1,ncores i1=isubend(icr-1)+1 i2=isubend(icr) if ((ksub.ge.i1).and.(ksub.le.i2)) then ef = efermi(icr) goto 11 endif enddo 11 continue C-RDC 11 write(iout,'(//" SUBSYSTEM ",I4," (EFERMI=",F12.5," eV)", C-RDC & /" -----------------------------------")') ksub,ef k1 = iorbpt(ksub)-1 k2 = iorbpt(ksub+1)-1 ipr2 = 2*ipr ipr1 = ipr2-1 j1 = 0 found = .false. do k=k1+1,k2 j1 = j1 + 1 if (eval(k).ge.eprsub(ipr1)) then found = .true. kk=k goto 20 endif enddo 20 if (found) then maxj2 = j1 do k=kk,k2 if (eval(k).gt.eprsub(ipr2)) goto 30 maxj2 = maxj2 + 1 enddo else C . did not find any eigenvector within the requested energies C-RDC write(iout,'(/" NONE FOUND WITHIN THE ", C-RDC $ "REQUESTED ENERGY INTERVAL")') norbs = iorbpt(ksub+1)-iorbpt(ksub) nk2 = norbs*norbs ifrpn = ifrpn + nk2 goto 1000 endif 30 maxj2=maxj2-1 norbs = maxj2-j1+1 if (norbs.eq.0) then C-RDC write(iout,'(/" NONE FOUND WITHIN THE ", C-RDC $ "REQUESTED ENERGY INTERVAL")') norbs = iorbpt(ksub+1)-iorbpt(ksub) nk2 = norbs*norbs ifrpn = ifrpn + nk2 goto 1000 endif NBLOCK = NORBS/6 + MIN(1,MOD(NORBS,6)) NORBS = IORBPT(KSUB+1)-IORBPT(KSUB) nk2 = norbs*norbs k1 = iorbpt(ksub)-1 DO 500 IBLOCK=1,NBLOCK IAO = 0 J2 = MIN(J1+5,maxj2) kk = 0 do j=j1,j2 kk = kk + 1 if (eval(j+k1).lt.ef) then occ(kk) = 'Y' else occ(kk) = 'N' endif enddo C-RDC WRITE(IOUT,'(/" MO -> ",I4,5(7X,I4))') C-RDC & (j,j=j1,j2) C-RDC write(iout,'(" ENERGY ->",6(1X,F10.4))') C-RDC . (EVAL(J+k1),J=J1,J2) C-RDC write(iout,'(" OCC. ->",6(6x,a1,4x))') C-RDC $ (occ(j-j1+1),j=j1,j2) C-RDC WRITE(IOUT,'("")') NATMS = IATOM1(KSUB+1) - IATOM1(KSUB) iatom0 = iatom1(ksub)-1 DO 200 IIATM=1,NATMS IATM = IATOMS(IATOM0+iiatm) IAI = IATNUM(IATM) NORBSI = NATORB(IAI) DO 100 IORB=1,NORBSI IAO = IAO + 1 C-RDC WRITE(IOUT,'(i4,1x,I4,1X,A2,1X,A2,6(1X,F10.7))') C-RDC . IAO,iatm,SYMBOL(IAI), C-RDC . LABEL(IORB), C-RDC . (EVEC(ifrpn+IAO+(J-1)*norbs),J=J1,J2) 100 CONTINUE 200 CONTINUE C J1 = J1 + 6 500 CONTINUE ifrpn = ifrpn + nk2 1000 continue RETURN END CC---------------------------------------------------------------CC SUBROUTINE WRTVECSUB(ksub) C C-RDC C WRITES EIGENVECTORS (MOs) TO UNIT IOUT. C for a DC calculation: one extra SCF iteration will be performed C to print the eigenvectors during the last cycle C C By Arjan van der Vaart C IMPLICIT DOUBLE PRECISION (A-H,O-Z) #include "divcon.dim" #include "divcon.h" C C C LOCAL: C data ipr /0/ CHARACTER LABEL(4)*2 DATA LABEL /' S','PX','PY','PZ'/ SAVE LABEL save ipr character*1 occ dimension occ(6) logical found C if (ksub.eq.1) then C-RDC WRITE(IOUT,'(///" FINAL EIGENVECTORS:")') ipr = 1 endif C find out if we need to print the eigenvectors of this subsystem 10 if (ipr.gt.nprsub) return if (iprsub(ipr).gt.ksub) then return elseif (iprsub(ipr).lt.ksub) then ipr = ipr + 1 goto 10 endif C find the Fermi energy: base occupation on the Fermi energy, C since we don't know the charge of this subsystem do icr=1,ncores i1=isubend(icr-1)+1 i2=isubend(icr) if ((ksub.ge.i1).and.(ksub.le.i2)) then ef = efermi(icr) goto 11 endif enddo c c c MO -> IIII IIII IIII IIII IIII IIII c ENERGY -> -FFFF.FFFF -FFFF.FFFF -FFFF.FFFF -FFFF.FFFF -FFFF.FFFF -FFFF.FFFF C OCC. -> Y Y N N N N c cIIII BR PX -FF.FFFFFF -FF.FFFFFF -FF.FFFFFF -FF.FFFFFF -FF.FFFFFF -FF.FFFFFF 11 continue C-RDC 11 write(iout,'(//" SUBSYSTEM ",I4," (EFERMI=",F12.5,")", C-RDC & /" -----------------------------------")') ksub,ef k1 = iorbpt(ksub)-1 k2 = iorbpt(ksub+1)-1 ipr2 = 2*ipr ipr1 = ipr2-1 j1 = 0 found = .false. do k=k1+1,k2 j1 = j1 + 1 if (eval(k).ge.eprsub(ipr1)) then found = .true. kk=k goto 20 endif enddo 20 if (found) then maxj2 = j1 do k=kk,k2 if (eval(k).gt.eprsub(ipr2)) goto 30 maxj2 = maxj2 + 1 enddo else C . did not find any eigenvector within the requested energies C-RDC write(iout, C-RDC $ '(/" NONE FOUND WITHIN THE REQUESTED ENERGY INTERVAL")') return endif 30 maxj2=maxj2-1 norbs = maxj2-j1+1 if (norbs.eq.0) then C-RDC write(iout, C-RDC $ '(/" NONE FOUND WITHIN THE REQUESTED ENERGY INTERVAL")') return endif NBLOCK = NORBS/6 + MIN(1,MOD(NORBS,6)) C-RDC C now we have nblock, overwrite norbs with the number of orbitals C in this subsystem, needed as offset when we print the eigenvectors NORBS = IORBPT(KSUB+1)-IORBPT(KSUB) DO 500 IBLOCK=1,NBLOCK IAO = 0 J2 = MIN(J1+5,maxj2) kk = 0 do j=j1,j2 kk = kk + 1 if (eval(j+k1).lt.ef) then occ(kk) = 'Y' else occ(kk) = 'N' endif enddo C-RDC WRITE(IOUT,'(/" MO -> ",I4,5(7X,I4))') C-RDC & (j,j=j1,j2) C-RDC write(iout,'(" ENERGY ->",6(1X,F10.4))') C-RDC . (EVAL(J+k1),J=J1,J2) C-RDC write(iout,'(" OCC. ->",6(6x,a1,4x))') C-RDC $ (occ(j-j1+1),j=j1,j2) C-RDC WRITE(IOUT,'("")') NATMS = IATOM1(KSUB+1) - IATOM1(KSUB) iatom0 = iatom1(ksub)-1 DO 200 IIATM=1,NATMS IATM = IATOMS(IATOM0+iiatm) IAI = IATNUM(IATM) NORBSI = NATORB(IAI) DO 100 IORB=1,NORBSI IAO = IAO + 1 C-RDC WRITE(IOUT,'(i4,1x,I4,1X,A2,1X,A2,6(1X,F10.7))') C-RDC . IAO,iatm,SYMBOL(IAI), C-RDC . LABEL(IORB), C-RDC . (EVEC(IAO+(J-1)*norbs),J=J1,J2) 100 CONTINUE 200 CONTINUE C J1 = J1 + 6 500 CONTINUE RETURN END