一般病态方程组初始解的迭代改善
Basic程序
40 INPUT "N,M,E=";N,M,E
50 INPUT "K1=";K1
60 DIM A(N,M),R(N),Z(N),B1(N),X(N),C(N,N)
70 PRINT TAB(3);"FOREXAMPLE":PRINT
80 FOR I=1 TO N
90 FOR J=1 TO M
100 READ A(I,J):PRINT USING"####.######";A(I,J);:PRINT;
110 NEXT J
120 PRINT
130 NEXT I
140 FOR I=1 TO N:B1(I)=A(I,M):NEXT I
150 FOR I=1 TO N
160 FOR J=1 TO N
170 C(I,J)=A(I,J)
180 NEXT J
190 NEXT I
200 PRINT:GOSUB 400
210 PRINT TAB(5);"W=";W:PRINT
220 IF W=1 THEN 340
230 PRINT TAB(3);"MODIFY"
240 FOR G=1 TO N
250 PRINT TAB(5);"X(";G;")=";A(G,M)
260 NEXT G:PRINT
270 GOSUB 830
280 PRINT TAB(3);"PRESENT";TAB(15);"K=";K
290 FOR I=1 TO N
300 PRINT TAB(5);"X(";I;")=";Z(I)
310 NEXT I
320 DATA 5,7,6,5,23,7,10,8,7,32,6,8,10,9,33,5,7,9,10,31
330 DATA 11,10,4,1,12,11,-13,1,14,13,-66,1
340 DATA 1,0.5,0.333333,0.25,1,0.5,0.333333,0.25,0.2,0
350 DATA 0.333333,0.25,0.2,0.2,0.166667,0,0.25,0.2,0.166667
360 DATA 0.142857,0
370 END
400 'SUB(1)
410 IF ABS(A(1,1))<=E THEN 820
420 FOR J=2 TO N
430 A(I,J)=A(I,J)/A(1,1)
440 NEXT J
450 P=0
460 FOR S=2 TO N
470 FOR I=S TO N
480 FOR R=1 TO S-1
490 P=P+A(I,R)*A(R,S)
500 NEXT R
510 A(I,S)=A(I,S)-P
520 P=0
530 NEXT I
540 P=0
550 FOR J=S+1 TO N
560 FOR R=1 TO S-1
570 P=P+A(S,R)*A(R,J)
580 NEXT R
590 IF ABS(A(S,S))<=E THEN 820
600 A(S,J)=(A(S,J)-P)/A(S,S)
610 P=0
620 NEXT J
630 NEXT S
640 W=0
650 '***
660 A(1,M)=(A(1,M))/A(1,1)
670 FOR I=2 TO N
680 Z=0
690 FOR J=1 TO I-1
700 Z=Z+A(I,J)*A(J,M)
710 NEXT J
720 A(I,M)=(A(I,M)-Z)/A(I,I)
730 NEXT I
740 FOR I=N-1 TO 1 STEP -1
750 Z=0
760 FOR J=N TO I+1 STEP -1
770 Z=Z+A(I,J)*A(J,M)
780 NEXT J
790 A(I,M)=A(I,M)-Z
800 NEXT I
810 RETURN
820 W=1:RETURN
830 'SUB(2)
840 K=0
850 FOR I=1 TO N
860 T#=0
870 IF K=0 THEN X(I)=A(I,M) ELSE X(I)=Z(I)
880 FOR J=1 TO N
890 T#=T#+CDBL(C(I,J))*CDBL(A(J,N))
900 NEXT J
910 PRINT TAB(5);T#
920 R(I)=CDBL(B1(I))-T#
930 PRINT R(I)
940 NEXT I
950 FOR I=1 TO N:A(I,M)=R(I):NEXT I
960 K=K+1
970 GOSUB 660
980 FOR I=1 TO N:Z(I)=X(I)+A(I,M):NEXT I
990 T=0
1000 FOR I=1 TO N
1010 T=T+(A(I,M))^2
1020 NEXT I
1030 T=SQR(T)
1040 IF T>E THEN 1060
1050 GOTO 1090
1060 IF K>K1 THEN 1090
1070 FOR I=1 TO N:A(I,M)=Z(I):NEXT I
1080 GOTO 850
1090 PRINT TAB(5);T
1100 RETURN
