Bairstow in matlab
Aiuto... Devo implementare il metodo di Bairstow in Matlab senza utilizzare funzioni proprie di matlab... Come faccio? Sto esame mi fà diventare matta... Aiutatemi...
Risposte
Riporto qui il testo di un post da me scritto sul forum generale un paio di anni fa… spero possa esserti utile…
cari amici
ho notato che alcuni di voi hanno mostrato interesse verso il procedimento di soluzione di un’equazione algebrica noto come ‘algoritmo di Bayrstow’. Dal momento che una dei primi programmi per computer che ho creato [tanti ma tanti anni fa, ahimè!…] implementava appunto questo algoritmo perso di fare cosa gradita presentandovi l’algoritmo insieme al programma da me creato nel lontano passato…
Un’equazione algebrica a coefficienti reali può essere scritta nella forma…
$p(x)=0$ (1)
… dove è…
$p(x)= a_n x^n+a_(n-1) x^(n-1)+…+a_0$ (2)
Il metodo ora descritto consente di isolare un fattore quadratico di p(x), dal quale è possibile ottenere una coppia di radici [reali o complesse e coniugate…] dell’equazione (1) e operare poi con un polinomio di grado n-2 fino al calcolo completo di tutte le radici. Per ottenere questo ci serviamo della seguente identità…
$p(x)= a_n x^n+ a_(n-1) x^(n-1)+…+a_0 = (x^2+u*x+v)*q(x)+r(x)$ (3)
… dove poniamo…
$q(x)= b_nx^(n-2)+b_(n-1)*x^(n-3)+…+b_2$
$r(x)= b_1*(x-u)+b_0$ (4)
Nel caso in cui fosse $b_1=b_0=0$ il polinomio q(x) divide esattamente p(x), il polinomio $x^2+u x+v$ è un fattore quadratico di p(x) e questo ci consente di isolare due radici e continuare la ricerca su q(x). Se così non è, considerato il fatto che sia $b_1$ sia $b_0$ sono funzioni di u e v, per trovare i valori di u e v è necessario risolvere il sistema [non lineare] di due equazioni in due incognite…
$b_1(u,v)=0, b_0(u,v)=0$ (5)
Ciò può essere fatto in maniera iterativa nel modo seguente. Supponiamo di scegliere per u e v arbitrariamente due valori $u_0$ e $v_0$ e di calcolare con essi $b_1 (u_0,v_0)$ e $b_0(u_0,v_0)$. Se è $b_1=b_0=0$ siamo a posto. Se no proveremo stimare i valori Du e Dv che occorre aggiungere a $u_0$ e $v_0$ per ottenere tale condizione. Se sviluppiamo $b_1$ e $b_2$ in serie di Taylor nell’intorno di $u_0$ e $v_0$ e tronchiamo la serie ai termini di grado 1 avremo…
$b_1(u,v)= b_1(u_0,v_0)+ Du (db_1)/(du) (u_0,v_0)+Dv (db_1)/(dv) (u_0,v_0)$
$b_0(u,v)=b_0(u_0,v_0)+Du (db_0)/(du) (u_0,v_0)+ Dv (db_0)/(dv) (u_0,v_0)$ (6)
… in cui compaiono le derivate parziali di $b_1$ e $b_0$ rispetto a u e v calcolate in $u_0$ e $v_0$. I valori Du e Dv che, nell’ipotesi di sufficiente accuratezza dell’approssimazione ‘lineare’ di $b_1$ e $b_0$ consentono di annullare entrambi saranno dunque la soluzione del sistema lineare…
$Du (db_1)/(du) (u_0,v_0)+Dv (db_1)/(dv) (u_0,v_0)=-b_1 (u_0,v_0)$
$Du (db_0)/(du) (u_0,v_0)+ Dv (db_0)/(dv) (u_0,v_0)=-b_0 (u_0,v_0)$ (7)
Se l’approssimazione raggiunta dalla prima iterazione sarà accettabile sin terranno buoni i valori di u e v trovati. Viceversa partendo da essi si opererà una nuova iterazione e cosi via…
Il primo passo consiste, dati u, v e le $a_i$, nel calcolare le $b_i$ per i=0,1,…,n. A tale scopo eseguiamo la moltiplicazione membro a membro illustrate nella (3) ed eguagliamo le potenze di x. Otteniamo…
$b_n=a_n$
$b_(n-1)=a_(n-1)+u b_n$
$b_(n-k)=a_(n-k)+u b_(n-k+1)+v b_(n-k+2)$ k=2,3,…,n (8)
Se poniamo $b_(n+1)=b_(n+2)=0$ l’ultima formula vale per k da 0 ad n. Compiuto il primo passo calcoliamo ora le derivate parziali che compaiono nella (6). Per fare questo procediamo in maniera simile a quella seguita per trovare le bi nel punto precedente. Per il calcolo delle derivate rispetto ad u poniamo $c_k=d b_k/du$. Otteniamo…
$c_n=(d b_n)/(du)=0$
$c_(n-1)=(d b_(n-1))/(du)=b_n+ u c_n
$c_(n-k)=(d b_(n-k))/(du)=b_(n-k+1)+u c_(n-k+1)+v c_(n-k+2)$ k=2,3,…,n (9)
Un semplice raffronto tra le (8) e le (9) ci rivela che le ci si ricavano dalle bi esattamente nello steso modo in cui le $b_i$ sono ricavate dalle $a_i$. Poniamo ora $d_k=(d b_k) /dv$. Otteniamo…
$d_n=(d b_n)/(dv)=0$
$d_(n-1)=(d b_(n-1))/(dv)=0$
$d_(n-k)=(d b_(n-k))/(dv)= b_(n-k+2)+u d_(n-k+1)+v d_(n-k+2)$ k=2,3,...,n (10)
Un semplice raffronto tra le (9) e le (10) ci dice che è $d_(n-k)=c_(n-k+1)$ per cui alla fine abbiamo tutto ciò che serve per impostare il sistema (7). Otteniamo dunque…
$(d b_1)/(du)= c_1, (db_1)/(dv)=c_2, (db_0)/(du)=c_0, (d b_0)/(dv)=c_1$ (11)
… per cui il sistema (7) diviene…
$c_1 Du+c_2*Dv=-b_1$
$c_0 Du+c_1 Dv=-b_0$ (12)
… e può essere risolto facilmente, ad esempio con la regola di Kramer.
Questa l’impostazione generale dell’algoritmo. Alcune precauzioni debbono essere seguite nell’applicarlo se non ci si vuole imbattere in spiacevoli sorprese. I due punti critici sono i seguenti, entrambi legati alla scelta dei valori iniziali [arbitrari] di u e v in ciascuna iterazione…
a) ad un certo punto può accadere che si incappi in un sistema lineare [espresso dalla (12)] mal condizionato. In tal caso è opportuno diagnosticare questo e ripartire con diversi valori di u e v
b) la convergenza dell’algoritmo non è garantita. Se una iterazione dura troppo tempo senza convergere ad una soluzione o addirittura diverge è necessario, come nel caso a), ripartire con diversi valori di u e v
Nel seguito è riportato un programma da me scritto oltre vent’anni fa in HT Basic che implementa l’algoritmo ora descritto. Ai giovani il compito di tradurlo in un linguaggio più moderno [eliminando in particolare tutte le istruzioni GO TO, ‘messe all’indice’ dagli ideatori del linguaggio C…]
cordiali saluti!…
lupo grigio
10 PRINTER IS 1
20 CLEAR SCREEN
30 PRINT “* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *”
40 PRINT “* *”
50 PRINT “* PENSIERO - Ricerca delle radici di una equazione algebrica *”
60 PRINT “* *”
70 PRINT “* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *”
80 PRINT
90 OPTION BASE 0
100 Err$=”NO”
110 INTEGER Gr,Nfq, Nfu
120 REAL A(256),B(256),C(256),P(256),Q(256),T(256),H(256,1)
130 REAL Ao,A1,B1,B2,C1,C2
140 DISP “Grado del polinomio (>1)”;
150 INPUT Gr
160 IF Gr<2 THEN 140
170 PRINT “Grado del polinomio:”;Gr
180 PRINT
190 PRINT “Coefficienti del polinomio”;
200 PRINT
210 PRINT
220 FOR I=Gr TO 0 STEP –1
230 DISP “Coefficiente di grado:”;I;
240 INPUT P(I)
250 PRINT “a(“;I;”):”;TAB(20);P(I)
260 NEXT I
270 PRINT
280 PRINT
290 Ko=1/P(Gr)
300 FOR I=0 TO Gr
310 P(I)=P(I)*Ko
320 NEXT I
330 A(Gr)=1
340 Nfq=INT(Gr/2)!Numero fattori quadratici
350 Nfu=Gr-2*Nfq!Numero fattori unitari
360 MAT A=P
370 GOSUB Bairstow
380 IF Err$=”SI” THEN 810
390 PRINT “Fattori quadratici trovati”;
400 PRINT
410 PRINT
420 FOR I=1 TO Nfq
430 PRINT “p^2”;
440 IF H(I,1)<0 THEN 460
450 PRINT “+”;
460 PRINT H(I,1);”p”;
470 IF H(I,0)<0 THEN 490
480 PRINT “+”;
490 PRINT H(I,0)
500 NEXT I
510 IF Nfu=0 THEN 560
520 PRINT “p”;
530 IF Ao<0 THEN 550
540 PRINT “+”
550 PRINT Ao
560 PRINT
570 PRINT
580 PRINT “Radici del polinomio”
590 PRINT
600 PRINT TAB(25);”Parte reale”;TAB(50);”Parte immaginaria”
610 PRINT
620 FOR I=1 TO Nfq
630 Delta=H(I,1)^2-4*H(I,0)
640 I1=2*I
650 I2=I1+1
660 IF Delta<0 THEN 700
670 Q(I1)=-.5*(H(I,1)+SQR(Delta))
680 Q(I2)=-.5*(H(I,1)-SQR(Delta))
690 GOTO 740
700 Q(I1)=-.5*H(I,1)
710 Q(I2)=Q(I1)
720 T(I1)=-.5*SQR(-Delta)
730 T(I2)=-T(I1)
740 PRINT TAB(20);Q(I1);TAB(45);T(I1)
750 PRINT TAB(20);Q(I2);TAB(45);T(I2)
760 NEXT I
770 IF Nfq=0 THEN 800
780 PRINT TAB(20);-Ao;TAB(45);”.00000000”
790 PRINT
800 GOTO 840
810 PRINT “Algoritmo fermatosi, non è possibile fornire radici”
820 PRINT
830 PRINT
840 STOP
850 Bairstow: Ip=Gr
860 FOR I=1 TO Nfq
870 Nit=0
880 Ao=RND!L’istruzione RND genera un numero casuale
890 A1=RND!Idem come sopra
900 DISP “I:”;I;”Nit:”;Nit
910 Ip1=Ip-1
920 Ip2=Ip-2
930 B(Ip)=0
940 B(Ip1)=0
950 FOR J=Ip2 TO 0 STEP –1
960 J1=J+1
970 J2=J+2
980 B(J)=A(J2)-A1*B(J1)-Ao*B(J2)
990 NEXT J
1000 B1= A(1)-A1*B(0)-Ao*B(1)
1010 B2=A(0)-A1*B1-Ao*B(0)!Abbiamo b1 e b2
1020 IF B1=0 AND B2=0 THEN 1730
1030 C(Ip)=0
1040 C(Ip1)=0
1050 FOR J=Ip2 TO 0 STEP –1
1060 J1=J+1
1070 J2=J+2
1080 C(J)=-B(J1)-A1*C(J1)-Ao*C(J2)
1090 NEXT J
1100 C1=-B(0)-A1*C(0)-Ao*C(1)
1110 C2=-B1-A1*C1-Ao*C(0)
1120 Delta= C1^2-C(0)*C2
1130 Md=SQR(2*C1^2+C(0)^2+C2^2)
1140 IF ABS(Delta/Md)<1.E-8 THEN 1320
1150 Da1=(C(0)*B2-C1*B2)/Delta
1160 Dao=(C2*B1-C1*B2)/Delta
1170 A1=A1+Da1
1180 Ao=Ao+Dao
1190 IF A1=0 OR Ao=0 THEN 1250
1200 Amax=MAX(A1,Ao)
1210 IF ABS(Da1/Amax)<1.E-12 AND ABS(Dao/Amax)<1.E-12 THEN 1250
1220 GOTO 910
1230 IF ABS(Da1)<1.E-12 AND ABS(Dao)<1.E-12 THEN 1250
1240 GOTO 910
1250 H(I,1)=A1
1260 H(I,0)=Ao
1270 Ip=Ip2
1280 FOR J=0 TO Ip2
1290 A(J)=B(J)
1300 NEXT J
1310 GOTO 1350
1320 IF Nit=1000 THEN 1400
1330 Nit=Nit+1
1340 GOTO 880
1350 NEXT I
1360 IF Nfu=0 THEN 1400
1370 Ao=A(0)
1380 GOTO 1400
1390 Err$=”SI”
1400 RETURN
1410 STOP
1420 END
cari amici
ho notato che alcuni di voi hanno mostrato interesse verso il procedimento di soluzione di un’equazione algebrica noto come ‘algoritmo di Bayrstow’. Dal momento che una dei primi programmi per computer che ho creato [tanti ma tanti anni fa, ahimè!…] implementava appunto questo algoritmo perso di fare cosa gradita presentandovi l’algoritmo insieme al programma da me creato nel lontano passato…
Un’equazione algebrica a coefficienti reali può essere scritta nella forma…
$p(x)=0$ (1)
… dove è…
$p(x)= a_n x^n+a_(n-1) x^(n-1)+…+a_0$ (2)
Il metodo ora descritto consente di isolare un fattore quadratico di p(x), dal quale è possibile ottenere una coppia di radici [reali o complesse e coniugate…] dell’equazione (1) e operare poi con un polinomio di grado n-2 fino al calcolo completo di tutte le radici. Per ottenere questo ci serviamo della seguente identità…
$p(x)= a_n x^n+ a_(n-1) x^(n-1)+…+a_0 = (x^2+u*x+v)*q(x)+r(x)$ (3)
… dove poniamo…
$q(x)= b_nx^(n-2)+b_(n-1)*x^(n-3)+…+b_2$
$r(x)= b_1*(x-u)+b_0$ (4)
Nel caso in cui fosse $b_1=b_0=0$ il polinomio q(x) divide esattamente p(x), il polinomio $x^2+u x+v$ è un fattore quadratico di p(x) e questo ci consente di isolare due radici e continuare la ricerca su q(x). Se così non è, considerato il fatto che sia $b_1$ sia $b_0$ sono funzioni di u e v, per trovare i valori di u e v è necessario risolvere il sistema [non lineare] di due equazioni in due incognite…
$b_1(u,v)=0, b_0(u,v)=0$ (5)
Ciò può essere fatto in maniera iterativa nel modo seguente. Supponiamo di scegliere per u e v arbitrariamente due valori $u_0$ e $v_0$ e di calcolare con essi $b_1 (u_0,v_0)$ e $b_0(u_0,v_0)$. Se è $b_1=b_0=0$ siamo a posto. Se no proveremo stimare i valori Du e Dv che occorre aggiungere a $u_0$ e $v_0$ per ottenere tale condizione. Se sviluppiamo $b_1$ e $b_2$ in serie di Taylor nell’intorno di $u_0$ e $v_0$ e tronchiamo la serie ai termini di grado 1 avremo…
$b_1(u,v)= b_1(u_0,v_0)+ Du (db_1)/(du) (u_0,v_0)+Dv (db_1)/(dv) (u_0,v_0)$
$b_0(u,v)=b_0(u_0,v_0)+Du (db_0)/(du) (u_0,v_0)+ Dv (db_0)/(dv) (u_0,v_0)$ (6)
… in cui compaiono le derivate parziali di $b_1$ e $b_0$ rispetto a u e v calcolate in $u_0$ e $v_0$. I valori Du e Dv che, nell’ipotesi di sufficiente accuratezza dell’approssimazione ‘lineare’ di $b_1$ e $b_0$ consentono di annullare entrambi saranno dunque la soluzione del sistema lineare…
$Du (db_1)/(du) (u_0,v_0)+Dv (db_1)/(dv) (u_0,v_0)=-b_1 (u_0,v_0)$
$Du (db_0)/(du) (u_0,v_0)+ Dv (db_0)/(dv) (u_0,v_0)=-b_0 (u_0,v_0)$ (7)
Se l’approssimazione raggiunta dalla prima iterazione sarà accettabile sin terranno buoni i valori di u e v trovati. Viceversa partendo da essi si opererà una nuova iterazione e cosi via…
Il primo passo consiste, dati u, v e le $a_i$, nel calcolare le $b_i$ per i=0,1,…,n. A tale scopo eseguiamo la moltiplicazione membro a membro illustrate nella (3) ed eguagliamo le potenze di x. Otteniamo…
$b_n=a_n$
$b_(n-1)=a_(n-1)+u b_n$
$b_(n-k)=a_(n-k)+u b_(n-k+1)+v b_(n-k+2)$ k=2,3,…,n (8)
Se poniamo $b_(n+1)=b_(n+2)=0$ l’ultima formula vale per k da 0 ad n. Compiuto il primo passo calcoliamo ora le derivate parziali che compaiono nella (6). Per fare questo procediamo in maniera simile a quella seguita per trovare le bi nel punto precedente. Per il calcolo delle derivate rispetto ad u poniamo $c_k=d b_k/du$. Otteniamo…
$c_n=(d b_n)/(du)=0$
$c_(n-1)=(d b_(n-1))/(du)=b_n+ u c_n
$c_(n-k)=(d b_(n-k))/(du)=b_(n-k+1)+u c_(n-k+1)+v c_(n-k+2)$ k=2,3,…,n (9)
Un semplice raffronto tra le (8) e le (9) ci rivela che le ci si ricavano dalle bi esattamente nello steso modo in cui le $b_i$ sono ricavate dalle $a_i$. Poniamo ora $d_k=(d b_k) /dv$. Otteniamo…
$d_n=(d b_n)/(dv)=0$
$d_(n-1)=(d b_(n-1))/(dv)=0$
$d_(n-k)=(d b_(n-k))/(dv)= b_(n-k+2)+u d_(n-k+1)+v d_(n-k+2)$ k=2,3,...,n (10)
Un semplice raffronto tra le (9) e le (10) ci dice che è $d_(n-k)=c_(n-k+1)$ per cui alla fine abbiamo tutto ciò che serve per impostare il sistema (7). Otteniamo dunque…
$(d b_1)/(du)= c_1, (db_1)/(dv)=c_2, (db_0)/(du)=c_0, (d b_0)/(dv)=c_1$ (11)
… per cui il sistema (7) diviene…
$c_1 Du+c_2*Dv=-b_1$
$c_0 Du+c_1 Dv=-b_0$ (12)
… e può essere risolto facilmente, ad esempio con la regola di Kramer.
Questa l’impostazione generale dell’algoritmo. Alcune precauzioni debbono essere seguite nell’applicarlo se non ci si vuole imbattere in spiacevoli sorprese. I due punti critici sono i seguenti, entrambi legati alla scelta dei valori iniziali [arbitrari] di u e v in ciascuna iterazione…
a) ad un certo punto può accadere che si incappi in un sistema lineare [espresso dalla (12)] mal condizionato. In tal caso è opportuno diagnosticare questo e ripartire con diversi valori di u e v
b) la convergenza dell’algoritmo non è garantita. Se una iterazione dura troppo tempo senza convergere ad una soluzione o addirittura diverge è necessario, come nel caso a), ripartire con diversi valori di u e v
Nel seguito è riportato un programma da me scritto oltre vent’anni fa in HT Basic che implementa l’algoritmo ora descritto. Ai giovani il compito di tradurlo in un linguaggio più moderno [eliminando in particolare tutte le istruzioni GO TO, ‘messe all’indice’ dagli ideatori del linguaggio C…]
cordiali saluti!…
lupo grigio
10 PRINTER IS 1
20 CLEAR SCREEN
30 PRINT “* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *”
40 PRINT “* *”
50 PRINT “* PENSIERO - Ricerca delle radici di una equazione algebrica *”
60 PRINT “* *”
70 PRINT “* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *”
80 PRINT
90 OPTION BASE 0
100 Err$=”NO”
110 INTEGER Gr,Nfq, Nfu
120 REAL A(256),B(256),C(256),P(256),Q(256),T(256),H(256,1)
130 REAL Ao,A1,B1,B2,C1,C2
140 DISP “Grado del polinomio (>1)”;
150 INPUT Gr
160 IF Gr<2 THEN 140
170 PRINT “Grado del polinomio:”;Gr
180 PRINT
190 PRINT “Coefficienti del polinomio”;
200 PRINT
210 PRINT
220 FOR I=Gr TO 0 STEP –1
230 DISP “Coefficiente di grado:”;I;
240 INPUT P(I)
250 PRINT “a(“;I;”):”;TAB(20);P(I)
260 NEXT I
270 PRINT
280 PRINT
290 Ko=1/P(Gr)
300 FOR I=0 TO Gr
310 P(I)=P(I)*Ko
320 NEXT I
330 A(Gr)=1
340 Nfq=INT(Gr/2)!Numero fattori quadratici
350 Nfu=Gr-2*Nfq!Numero fattori unitari
360 MAT A=P
370 GOSUB Bairstow
380 IF Err$=”SI” THEN 810
390 PRINT “Fattori quadratici trovati”;
400 PRINT
410 PRINT
420 FOR I=1 TO Nfq
430 PRINT “p^2”;
440 IF H(I,1)<0 THEN 460
450 PRINT “+”;
460 PRINT H(I,1);”p”;
470 IF H(I,0)<0 THEN 490
480 PRINT “+”;
490 PRINT H(I,0)
500 NEXT I
510 IF Nfu=0 THEN 560
520 PRINT “p”;
530 IF Ao<0 THEN 550
540 PRINT “+”
550 PRINT Ao
560 PRINT
570 PRINT
580 PRINT “Radici del polinomio”
590 PRINT
600 PRINT TAB(25);”Parte reale”;TAB(50);”Parte immaginaria”
610 PRINT
620 FOR I=1 TO Nfq
630 Delta=H(I,1)^2-4*H(I,0)
640 I1=2*I
650 I2=I1+1
660 IF Delta<0 THEN 700
670 Q(I1)=-.5*(H(I,1)+SQR(Delta))
680 Q(I2)=-.5*(H(I,1)-SQR(Delta))
690 GOTO 740
700 Q(I1)=-.5*H(I,1)
710 Q(I2)=Q(I1)
720 T(I1)=-.5*SQR(-Delta)
730 T(I2)=-T(I1)
740 PRINT TAB(20);Q(I1);TAB(45);T(I1)
750 PRINT TAB(20);Q(I2);TAB(45);T(I2)
760 NEXT I
770 IF Nfq=0 THEN 800
780 PRINT TAB(20);-Ao;TAB(45);”.00000000”
790 PRINT
800 GOTO 840
810 PRINT “Algoritmo fermatosi, non è possibile fornire radici”
820 PRINT
830 PRINT
840 STOP
850 Bairstow: Ip=Gr
860 FOR I=1 TO Nfq
870 Nit=0
880 Ao=RND!L’istruzione RND genera un numero casuale
890 A1=RND!Idem come sopra
900 DISP “I:”;I;”Nit:”;Nit
910 Ip1=Ip-1
920 Ip2=Ip-2
930 B(Ip)=0
940 B(Ip1)=0
950 FOR J=Ip2 TO 0 STEP –1
960 J1=J+1
970 J2=J+2
980 B(J)=A(J2)-A1*B(J1)-Ao*B(J2)
990 NEXT J
1000 B1= A(1)-A1*B(0)-Ao*B(1)
1010 B2=A(0)-A1*B1-Ao*B(0)!Abbiamo b1 e b2
1020 IF B1=0 AND B2=0 THEN 1730
1030 C(Ip)=0
1040 C(Ip1)=0
1050 FOR J=Ip2 TO 0 STEP –1
1060 J1=J+1
1070 J2=J+2
1080 C(J)=-B(J1)-A1*C(J1)-Ao*C(J2)
1090 NEXT J
1100 C1=-B(0)-A1*C(0)-Ao*C(1)
1110 C2=-B1-A1*C1-Ao*C(0)
1120 Delta= C1^2-C(0)*C2
1130 Md=SQR(2*C1^2+C(0)^2+C2^2)
1140 IF ABS(Delta/Md)<1.E-8 THEN 1320
1150 Da1=(C(0)*B2-C1*B2)/Delta
1160 Dao=(C2*B1-C1*B2)/Delta
1170 A1=A1+Da1
1180 Ao=Ao+Dao
1190 IF A1=0 OR Ao=0 THEN 1250
1200 Amax=MAX(A1,Ao)
1210 IF ABS(Da1/Amax)<1.E-12 AND ABS(Dao/Amax)<1.E-12 THEN 1250
1220 GOTO 910
1230 IF ABS(Da1)<1.E-12 AND ABS(Dao)<1.E-12 THEN 1250
1240 GOTO 910
1250 H(I,1)=A1
1260 H(I,0)=Ao
1270 Ip=Ip2
1280 FOR J=0 TO Ip2
1290 A(J)=B(J)
1300 NEXT J
1310 GOTO 1350
1320 IF Nit=1000 THEN 1400
1330 Nit=Nit+1
1340 GOTO 880
1350 NEXT I
1360 IF Nfu=0 THEN 1400
1370 Ao=A(0)
1380 GOTO 1400
1390 Err$=”SI”
1400 RETURN
1410 STOP
1420 END
Accedi a tutti gli appunti
Tutor AI: studia meglio e in meno tempo