100'  NAME--DEBYE
110'
120'  DESCRIPTION--COMPUTES THE DEBYE OR THE EINSTEIN FUNCTION.
130'  GIVEN TWO OF THE THREE VALUES (TEMPERATURE,SPECIFIC HEAT,
140'  AND THETA) THE PROGRAM WILL CALCULATE THE THIRD AND IN 
150'  ADDITION THE NORMALIZED ENERGY FUNCTION AT THE GIVEN TEMPERATURE.
160'
170'  SOURCE--FORTRAN CODE BY MRS. MARTIN REILLY OF THE NATIONAL 
180'  BUREAU OF STANDARDS, WASHINGTON D.C.  
190'  BASIC CODE BY MRS. CARLA MESSINA 221.04, NBS, 3/31/65.
200'
210'  INSTRUCTIONS--DATA IS ENTERED IN LINE 1350 AND FOLLOWING.
220'  DATA IS ENTERED IN SETS OF 3 NUMBERS, FIRST IS THE 
230'  TEMPERATURE, SECOND THE SPECIFIC HEAT, AND THIRD IS THETA.
240'  ENTER 0 FOR THE UNKNOWN VALUE.
250'  SAMPLE DATA IS IN LINE 1350.
260'
270'
280'  *  *  *  *  *  *   MAIN PROGRAM   *  *  *  *  *  *  *  *  *
290'
300 PRINT "DO YOU WANT THE DEBYE OR THE EINSTEIN FUNCTION";
310 INPUT T$
320 IF T$="DEBYE" THEN 410
330 IF T$="EINSTEIN" THEN 360
340 PRINT "INCORRECT ANSWER.  PLEASE TYPE 'DEBYE' OR 'EINSTEIN'";
350 GOTO 310
360 PRINT "EINSTEIN FUNCTION"
370 LET T4=1
380 PRINT
390 PRINT " TEMPERATURE","   CV","   OMEGA","   Q","  CV/3R"
400 GO TO 440
410 PRINT "   DEBYE FUNCTION"
420 PRINT
430 PRINT " TEMPERATURE","   CV","   THETA","   Q","   CV/3R"
440 PRINT
450 READ T,C,O
460 LET M = 0
470 LET F = C/5.96151
480 IF T4 <>0 THEN 990
490 IF F <= 1 THEN 520
500 PRINT "NOT DEFINED FOR  T = "T;" C/3R = "F;"THETA = "O
510 GO TO 450
520 IF F< 1 THEN 540
530 GO TO 500
540 IF F < 0 THEN 500
550 IF F <> 0 THEN 720
560 LET M = 1
570 IF T < 0 THEN 500
580 IF T > 0 THEN 620
590 LET F4 = 0
600 LET F1 = 0
610 GO TO 950
620 LET X = O/T
630 IF X <= 0 THEN 500
640 IF X<=176 THEN 680
650 LET P = 0
660 LET S = 0
670 GO TO 790
680 LET R = EXP(X)
690 LET S = 3*X/(R-1)
700 IF X >=8 THEN 780
710 GO TO 740
720 LET X = 4.0
730 GO TO 680
740 LET D = (((X+15.121491)*X+143.155337)*X+682.0012)*X+3953.632
750 LET U =((((0.0946173*X)-4.432582)*X+85.07724)*X-800.6087)*X+3953.632
760 LET F4 = U/D
770 GO TO 800
780 LET P=((((6.0/X+6.0)/X+3.0)/X+1.0)/R)*3.0
790 LET F4 = (19.48182/X^3)-P
800 LET F1 = 4*F4 - S
810 IF M<> 0 THEN 890
820 LET F2 = F - F1
830 LET T1 = F2/F
840 IF ABS(T1)<=.000001 THEN 890
850 LET F3 = (((R/3.0)*S)/3.0)*S - F1
860 LET X = (((F2/F3)*X)/3.0)+X
870 IF X <= 0 THEN 500
880 GO TO 640
890 IF T >0 THEN 910
900 LET  T = O/X
910 IF O > 0 THEN 940
920 LET O = T*X
930 GO TO 950
940 LET C = F1*5.96151
950 LET Q = F4
960 LET V = F1
970 PRINT T,C,O,Q,V
980 GO TO 450
990 LET O=O
1000 IF F-1 <=0 THEN 1020
1010 GO TO 500
1020 IF F >=0 THEN 1040
1030 GO TO 500
1040 IF F =0 THEN 1080
1050 LET X = 2.35
1060 LET M = 1
1070 GO TO 1160
1080 IF T <0 THEN 500
1090 IF T = 0 THEN 1130
1100 LET X = O/T
1110 IF X <=0 THEN500
1120 IF X <=176 THEN 1160
1130 LET S = 0
1140 LET P = 0
1150 GO TO 1250
1160 LET R = EXP(X)
1170 LET B = R-1.0
1180 LET S = X/B
1190 LET P = S*S*R
1200 IF M =0 THEN 1250
1210 LET Z =(1-P/F)
1220 IF ABS(Z)<=.000001 THEN 1250
1230 LET X = X + Z*X*B/(2.0*B-X*(B+2.0))
1240 GO TO 1110
1250 LET Q = S
1260 LET V = P
1270 IF T >0 THEN 1300
1280 LET T = O/X
1290 GO TO 970
1300 IF O > 0 THEN 1330
1310 LET O = T*X
1320 GO TO 970
1330 LET C = V*5.96151
1340 GO TO 970
1350 DATA 1000,0,1,100,0,1,10,0,1
1360 END