1! Lawrence University Department of Chemistry Appleton, Wisconsin 54911 (414) 739-3681, Ext. 456 2! Program name, version, date, authorship: HYBRID, Version 1A, 14-Dec-74 (Slave program for SPPLT) Written by Stephen L. Holmgren, class of 1973, and maintained by Dr. James S. Evans, Department of Chemistry, Lawrence University. 3! Purpose of program: SPPLT prepares contour plots of hybrid atomic orbitals with specified contour levels, effective nuclear charge, amounts of 2s and 2p character, and geometric operations (scaling, translation, rotation). 4! Included are subroutines for several plotting devices: Tektronix 4010 graphics terminal; TSP-212 plotter system; Hewlett-Packard 7200A/7202A graphic plotter. 5! Language and operating system: SPPLT, written in BASIC-PLUS, uses virtual core files and several non-privileged SYS functions, and runs in 8K under RSTS Versions 05-21 and 05B-24 (PDP-11/45). 6! Availability: SPPLT is a non-proprietary program product of Lawrence University. It is distributed at cost for educational use, on the condition that it not be sold, rented, or leased for profit. 7! Limitations: Variables in NMRSIM are dimensioned to handle up to 6 spin-1/2 nuclei, or fewer with larger spins. Larger systems might be feasible on an 11/45 configured for 16K of user core. 8! Further documentation: User's guide: 3 pages of run instructions and operational features for students who will use SPPLT for class assignments. 9! Programmer's guide: discussion of algorithms and instructions for system implementation. Literature reference: S. L. Holmgren and J. S. Evans, Journal of Chemical Education, vol. 51, pp. 189-191 (1974). 10! Disclaimer: Neither the author of this program nor Lawrence University assumes any liability, expressed or implied, with respect to the correctness or performance of this program. 40 PRINT "CAN'T RUN HYBRID" : GOTO 32000 50 GOSUB 5600 : M$=CHR$(K0%) 60 D=0.05 : A0=0.52915 : S=1.E-20 : A1(K1%)=F0 : C,S2=F1 65 INPUT "EFFECTIVE NUCLEAR CHARGE"; Z0 : M6%=K0% : GOTO 6040 UNLESS Z0>=F1 75 INPUT "AMOUNT OF P CHARACTER";A ! ALPHA IS ENTERED 77 IF A>=20 THEN S2=F0 ! DEFAULT TO PURE "P" IF A > 20 79 A=SQR(A) 86 INPUT "# OF CONTOURS";I5% : IF I5%>K3% THEN I5%=K3% 90 DIM P(3%),A1(38%),R1(37%,3%),R2(37%,3%),R3(37%,3%),R4(37%,3%) 100 DIM#1,Z1(37%,3%),Z2(37%,3%),Z3(37%,3%),Z4(37%,3%) : Z4(37%,K3%)=K0% ! OPEN AND EXTEND VIRT. CORE FILE 110 DEF FNF(S)=C1*E*(S2*(F2-S)+S1*A*C*S) ! PSI 120 DEF FNF1(S)=C2*E*(S2*(-F2+S)+S1*A*C*(F1-S)) ! D/DR(PSI) 140 Z1=Z0/A0 : Z9=F1/SQR(Z0) : C0=0.25/SQR(F2*PI) : C1=(C0*Z1^1.5)/SQR(S2+A*A) : C2=C1*Z1 ! NORMALIZATION CONST. 150 ! **********BEGINNING OF FIRST LOOP****** 160 FOR J%=K1% TO N% : J1%,J2%=K1% ! J1%,J2%=2 WHEN RESPECTIVE ROOTS FOUND 180 GOTO 400 UNLESS J%=K1% OR J%>11% ! MAX & MIN CHANGE LITTLE AT FIRST 205 IF A<=F1 THEN J1%=K2% : M1=F0 ! MAX PSI AT R=0 IF 0<=A<=1 209 R0=F0 ! ***********FIND MAX & MIN******* 210 FOR I%=K1% UNTIL J1%+J2%=K4% OR I%=100% 220 E=EXP(-Z1*R0/F2) 230 ON J1% GOTO 240,280 ! HAS MAX BEEN FOUND 240 S1=F1 250 IF FNF1(Z1*R0/F2)>=F0 GOTO 280 260 M1=(R0-0.025) : J1%=K2% ! ROOT SET HALFWAY BETWEEN PRESENT&PAST R0 280 ON J2% GOTO 290,320 ! HAS OTHER MAX BEEN FOUND 290 S1=-F1 ! NEGATIVE LOBE 300 IF FNF1(Z1*R0/F2)<=F0 GOTO 320 310 M2=(R0-.025) : J2%=K2% 320 R0=R0+.05 ! RADIUS INCREMENTED BY .05 ANGSTROMS 330 NEXT I% ! ******END OF MAX-MIN PART*** 334 GOTO 400 UNLESS J%=K1% 335 S1=F1 : E=EXP(-Z1*M1/F2) : PRINT "APPROX PSI MAX ="; FNF(Z1*M1) : PRINT " AT R="; M1 : INPUT "CHOOSE CONTOURS"; C$ 340 IF ASCII(C$)=89% THEN 350 ! TEST FOR Y[ES] 345 P(K1%)=0.03 : P(K2%)=0.06 : P(K3%)=0.12 ! DEFAULT VALUES OF CONTOURS 346 GOTO 354 350 PRINT "ENTER SMALLEST FIRST" : INPUT P(I%) FOR I%=K1% TO I5% 354 PRINT "CONTOURS USED:" 355 G$="##.### ##.### ##.###" 360 PRINT USING G$,P(I%); FOR I%=K1% TO I5% : PRINT 400 C=COS(A1(J%)) ! CALCULATE COSINE AT THIS ANGLE 410! PSIMAX,MIN,AND NODE POSITION ARE CALCULATED FOR EACH ANGLE 420 S1=+F1 : E=EXP(-Z1*M1/F2) : P2=FNF(Z1*M1) ! MAX OF POS. LOBE 430 S1=-F1 : E=EXP(-Z1*M2/F2) : P1=FNF(Z1*M2) ! MIN OF NEG. LOBE 440 M5,P4=F0 : M6=M2 : C$=" NODE " : R8=0.2 445 IF S2=F0 THEN R9=F0 ELSE GOSUB 5001 ! SOLVE FOR NODE IF NOT PURE 'P' 447 N0=R9 ! NODE IS ALWAYS 0 FOR PURE "P" ORBITAL 455 !.....INITIAL GUESSES FOR ROOTS ARE MADE 460 R1(K0%,K1%)=(M1+6.0*Z9)/F2 : R1(K0%,K2%),R1(K0%,K3%)=M1+.75 461 FOR I2%=K1% TO K3% 462 R2(K0%,I2%)=M1/F2 + S ! ADDED TO PREVENT DIV BY ZERO 463 R3(K0%,I2%)=(M2+N0)/F2 : R4(K0%,I2%)=(M2+6.0*Z9)/F2 464 NEXT I2% 470 FOR I%=K1% TO I5% ! *************BEGINNING OF SECOND LOOP**** 480 P4=-ABS(P(I%)) ! SELECT CONTOUR LEVEL 490 IF P419% GOTO 750 ! ROOTS 1&2 ARE REDUNDANT FOR THETA > 90 605 P4=ABS(P(I%)) 610 IF P4>=P2 THEN 750 ! CHECK TO SEE IF ROOTS 1&2 ARE POSSIBLE 620 ! WHAT SIDE OF ORIGIN IS ROOT 2 ON?? 625 E=F1 630 Q=P4-FNF(F0) ! COMPARE CONTOUR LEVEL TO FUNCTION AT ZERO 640 IF Q>F0 THEN 649 ELSE IF Q=30% 5025 E=EXP(-R9*Z1/F2) 5030 R3=(FNF(Z1*R9)-P4)/(FNF1(Z1*R9/F2)+S) 5035 IF ABS(R3)>10 THEN 5055 ! CHECK FOR DIVERGENCE 5040 R8=R9-R3 5045 D1=ABS(R3/(R9+S)) ! CHECK CONVERGENCE CRITERION (0.05) 5050 R9=R8 : GO TO 5070 5055 R9=-F1 5060 K%=30% 5070 NEXT K% 5080 IF R9M5 THEN RETURN ! IF ROOT IN PROPER RANGE,RETURN 5085 T=T+F2 ! ...T IS INCREMENTED AFTER EACH TRY FOR ROOT 5087 IF T>F3 THEN 5089 ELSE 5100 5089 ! ...TRY OTHER METHOD IF T>3 5090 IF Q=F0 THEN 5225 ELSE 5091 5091 R9=-F2 : RETURN 5100 R9=(M6-M5)*T/4+M5 ! TRY ANOTHER STARTING VALUE WITHIN THE INTERVAL 5150 GO TO 5015 5200! ADDITIONAL SUBROUTINE FOR STUBBORN ROOTS; USES THE SIMPLE INTERVAL-HALVING TECHNIQUE 5225 Q=F1 5230 R6=M6 : R5=M5 : R8=(R6+R5)/F2 5235 E=EXP(-R5*Z1/F2) : R3=FNF(Z1*R5)-P4 5240 E=EXP(-R8*Z1/F2) : IF R3*(FNF(Z1*R8)-P4)>F0 THEN 5270 5250 R9=(R8+R5)/F2 : R6=R8 5260 GOTO 5280 5270 R9=(R8+R6)/F2 : R5=R8 5280 D1=ABS((R9-R8)/R8) 5290 IF D1<.001 THEN 5300 5295 R8=R9 : GOTO 5235 5300 E=EXP(-R6*Z1/F2) : IF R3*(FNF(Z1*R6)-P4)