Trailing-Edge
-
PDP-10 Archives
-
decuslib10-01
-
43,50144/rooter.bas
There are 2 other files named rooter.bas in the archive. Click here to see a list.
10 DATA 1,3,6
11 DATA 2,1,-6,5
12 DATA 3,3,-7,0,4
305 DIM A(25),B(25),X(25)
310 LET G8=0
315 PRINT
320 LET G8=G8+1
325 PRINT
330 READ N
335 IF N=0 THEN 999
340 PRINT
345 PRINT "POLYNOMIAL NUMBER";G8;"IS OF ORDER";N
350 PRINT
355 PRINT " COEFFICIENTS (IN DESCENDING ORDER) ARE:"
360 PRINT
365 PRINT " ";
370 FOR I=0 TO N
375 READ A(I)
380 PRINT A(I);
385 LET B(I)=A(I)
390 NEXT I
395 PRINT
400 PRINT
405 PRINT " THE ROOTS ARE:"
410 PRINT
415 IF N<=2 THEN 665
420 IF A(N) = 0 THEN 705
425 IF (N/2-INT(N/2))=0 THEN 440
430 GOSUB 745
435 GOTO 415
440 IF ABS(A(N-2))< 1.E-25 THEN 460
445 LET P=A(N-1)/A(N-2)
450 LET Q=A(N)/A(N-2)
455 GOTO 470
460 LET P=A(N-1)
465 LET Q=A(N)
470 FOR I = 0 TO N
475 LET X(I)=A(I)
480 NEXT I
485 GOSUB 720
490 FOR I = 0 TO N-2
495 LET B(I)=X(I)
500 NEXT I
505 LET R=X(N-1)
510 LET S=A(N)-P*X(N-1)-Q*X(N-2)
515 GOSUB 720
520 LET X(N-1)=-P*X(N-2)-Q*X(N-3)
525 LET D=X(N-2)^2-X(N-1)*X(N-3)
530 IF ABS(D ) > 1.E-25 THEN 545
535 PRINT " UNOBTAINABLE WITH THIS PROGRAM. SORRY ABOUT THAT."
540 GOTO 315
545 LET P1=P+(R*X(N-2)-S*X(N-3))/D
550 LET Q1=Q+(S*X(N-2)-R*X(N-1))/D
555 IF ABS(P ) > 1.E-25 THEN 575
560 IF ABS(P1) > 1.E-25 THEN 575
565 IF ABS( Q) > 1.E-25 THEN 580
570 GOTO 585
575 IF ABS(P1/P-1)>.000001 THEN 585
580 IF ABS(Q1/Q-1)<.000001 THEN 600
585 LET P=P1
590 LET Q=Q1
595 GOTO 470
600 FOR I=1 TO N-2
605 LET A(I)=B(I)
610 NEXT I
615 LET N=N-2
620 LET D=P*P-4*Q
625 IF D<0 THEN 650
630 LET D=SQR(D)
635 PRINT " ";(-P+D)/2;" AND ";(-P-D)/2
640 IF N-2>0 THEN 420
645 GOTO 665
650 LET D=SQR(-D)
655 PRINT " ";-P/2;"+ J *";D/2;" AND ";-P/2;"- J *";D/2
660 IF N-2>0 THEN 420
665 IF N=1 THEN 695
670 IF N=0 THEN 315
675 LET P=B(1)/B(0)
680 LET Q=B(2)/B(0)
685 LET N=0
690 GOTO 620
695 PRINT " ";-B(1)/B(0)
700 GOTO 315
705 PRINT " 0.00000"
710 LET N=N-1
715 GOTO 415
720 LET X(1)=X(1)-P*X(0)
725 FOR I=2 TO N-1
730 LET X(I)=X(I)-P*X(I-1)-Q*X(I-2)
735 NEXT I
740 RETURN
745 IF B(1) = 0 THEN 760
750 LET X=-B(1)/B(0)
755 GOTO 765
760 LET X=-B(N)/B(0)
765 LET F=0
770 LET F1=0
775 FOR I=0 TO N
780 LET J=N-I
785 IF B(J)=0 THEN 805
790 LET F=B(J)*X^I+F
795 IF I = 0 THEN 805
800 LET F1=I*B(J)*X^(I-1)+F1
805 NEXT I
810 LET X1=X-F/F1
815 IF ABS(X/X1-1) < .000001 THEN 830
820 LET X=X1
825 GOTO 765
830 PRINT " ";X1
835 LET N=N-1
840 FOR I=1 TO N
845 LET A(I)=B(I)+X1*A(I-1)
850 LET B(I)=A(I)
855 NEXT I
860 RETURN
865 STOP
870 PRINT
875 PRINT "THIS PROGRAM FINDS THE ROOTS OF POLYNOMIALS"
880 PRINT "USING BAIRSTOW'S METHOD. TO USE, ENTER DATA"
885 PRINT "IN DATA STATEMENTS AS FOLLOWS:"
890 PRINT
895 PRINT " 10 DATA N, A(N), A(N-1),......, A(1), A(0)"
900 PRINT
905 PRINT "WHERE N IS THE ORDER OF THE POLYNOMIAL AND"
910 PRINT "A(I) IS THE COEFFICIENT OF THE I-TH DEGREE"
915 PRINT "TERM. MORE THAN ONE DATA LINE MAY BE USED"
920 PRINT "TO SUPPLY COEFFICIENTS FOR ONE POLYNOMIAL,"
925 PRINT "AND ADDITIONAL POLYNOMIALS MAY BE SOLVED ON"
930 PRINT "A SINGLE RUN BY SUPPLYING DATA FOR THEM ON"
935 PRINT "SUBSEQUENT DATA LINES (NOT BEYOND LINE 299)."
940 PRINT
945 PRINT "CAUTION: THERE ARE A FEW TYPES OF POLYNOMIALS"
950 PRINT "WHICH THIS PROGRAM IS UNABLE TO SOLVE. PROGRAM"
955 PRINT "WILL SO INDICATE AND GO ON TO NEXT CASE. ALSO,"
960 PRINT "FOR HIGH ORDER POLYNOMIALS, PROGRAM MAY REQUIRE"
965 PRINT "MANY ITERATIONS, SO RUNNING TIME CAN BE HIGH."
975 PRINT"ENTER DATA BEGINNING ON LINE 10 THEN TYPE 'RUN'"
998 DATA 0
999 END