Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- C THE ATOM HARTREE FOCK
- C 1-3-2016 Eric Lovejoy removes impicit declaration
- C organized code, and fixed typos
- PROGRAM HFOCK
- implicit none
- INTEGER I
- REAL A, B, C, E, F, T
- REAL C1, C2, ET
- REAL F11,F12,F22,H11,H12,H22,R11,R12,R22
- REAL S12,T11,T12,T22,V11,V12,V22
- REAL DIFF,EOLD,PLUS,PROD,C1C2
- REAL ZETA1,ZETA2
- REAL ER1111,ER1112,ER1122,ER1212,ER1222,ER2222
- ZETA1 = 2.173171
- ZETA2 = 1.188530
- WRITE (*, '( "Zeta1:", F7.4 / "Zeta2:", F7.4 )' ) ZETA1, ZETA2
- C COMPUTE INTEGRALS
- PLUS=ZETA1+ZETA2
- PROD=ZETA1*ZETA2
- DIFF=ZETA1-ZETA2
- S12=8.*PROD**1.5/PLUS**3
- T11=0.5*ZETA1**2
- T22=0.5*ZETA2**2
- T12=4.*PROD**2.5/PLUS**3
- V11=-2.*ZETA1
- V22=-2.*ZETA2
- V12=-8.*PROD**1.5/PLUS**2
- H11=T11+V11
- H12=T12+V12
- H22=T22+V22
- ER1111=5.*ZETA1/8.
- ER2222=5.*ZETA2/8.
- T=(DIFF/PLUS)**2
- ER1122=PLUS*(1.-T)*(5.-T)/16.
- ER1212=4.*PROD**3/PLUS**5
- F=0.25*PROD**1.5/PLUS**3
- T=(DIFF/(3.*ZETA1+ZETA2)**2)
- ER1112=F*(3.*ZETA1+ZETA2)*(1.-T)*(5.-T)
- T=(DIFF/(3.*ZETA2+ZETA1)**2)
- ER1222=F*(3.*ZETA2+ZETA1)*(1.-T)*(5.-T)
- R11=1.
- R12=0.
- R22=0.
- EOLD=0
- c variables
- c S12 = 0.874115, T11 = 2.361336, T12 = 1.128867, T22 = 0.706302
- c V11 = -4.346342, V12 = -2.938515, V22 = -2.377060
- c ER1111, = 1.358232, ER1222 = 0.756250
- c ER1112 = 1.032210 ER2222 = 0.742831
- c ER1122, = 0.943928 ER1212 = 0.802688
- WRITE(*,*) '***********************************************'
- WRITE(*,*) '*****System Variables: Actual vs Expected*****'
- WRITE(*,*) '***********************************************'
- WRITE(*,*) 'S12',S12, ' Expected:0.874115'
- WRITE(*,*) 'T11',T11, ' Expected:2.361336'
- WRITE(*,*) 'T12',T12, ' Expected:1.128867'
- WRITE(*,*) 'T22',T22, ' Expected:0.706302'
- WRITE(*,*) 'V11',V11, ' Expected:-4.346342'
- WRITE(*,*) 'V12',V12, ' Expected:-2.938515'
- WRITE(*,*) 'V22',V22, ' Expected:-2.377060'
- WRITE(*,*) 'ER1111',ER1111, ' Expected:1.358232'
- WRITE(*,*) 'ER1222',ER1222, ' Expected:0.756250'
- WRITE(*,*) 'ER1112',ER1112, ' Expected:1.032210'
- WRITE(*,*) 'ER2222',ER2222, ' Expected:0.742831'
- WRITE(*,*) 'ER1122',ER1122, ' Expected:0.943928'
- WRITE(*,*) 'ER1212',ER1212, ' Expected:0.802688'
- WRITE(*,*) ' '
- WRITE(*,*) '***********************************************'
- WRITE(*,*) '* ITERATION ORBITAL ENERGY TOTAL E EIGENVECTOR*'
- WRITE(*,*) '***********************************************'
- C START HARTREE-FOCK LOOP
- DO 3 I=1,15
- F11=H11+R11*ER1111+2.*R12*ER1112+R22*(2.*ER1122-ER1212)
- F12=H12+R11*ER1112+R12*(3.*ER1212-ER1122)+R22*ER1222
- F22=H22+R11*(2.*ER1122-ER1212)+2.*R12*ER1222+R22*ER2222
- A=1.-S12**2
- B=2.*F12*S12-F11-F22
- C=F11*F22-F12**2
- E=(-B-(B**2-4.*A*C)**0.5)/(2.*A)
- C1C2=-(F12-E*S12)/(F11-E)
- C2=1./(C1C2**2+2.*C1C2*S12+1.)**0.5
- C2=1./(C1C2**2+2.*C1C2*S12+1.)**0.5
- C1=C1C2*C2
- R11=C1**2
- R12=C1*C2
- R22=C2**2
- ET=E+R11*H11+2.*R12+R22*H22
- WRITE(*,*) I,E,ET,C1,C2
- IF((ABS(ET-EOLD)).LT.1.E-6) STOP
- EOLD=ET
- 3 CONTINUE
- STOP
- END PROGRAM HFOCK
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement