Esercizio analisi numerica
Esercizio analisi numerica
Ciao a tutti, c'è qualcuno cos'ì gentile da spiegarmi (come lo spiegherebbe a sua nonna ) la soluzione di questo esercizio svolto con Mathematica?
Ringrazio anticipatamente tutti coloro che mi risponderanno.
Vi riporto il testo e la soluzione dell'esercizio.
Risolvere il sistema Ax=b utizizzando il metodo di richardson.
Stabilire empiricamente il miglior valore di "alfa" per la convergenza.
A{aij}ij = 1.N aij=1 i diverso da j
aij= n+1 x= 1 ossia [1,1,1...1]
bi=2n per ogni i.
SOLUZIONE
e = 0.1;
α = 2;
n = 5;
b = Table[1, {i, 1, n}];
A = Table[1 + n*KroneckerDelta[i, j], {i, 1, n}, {j, 1, n}];
x0 = Table[1.5, {i, 1, n}];
r0 = b - A.x0;
x = x0;
r = r0;
i = 0;
While[Norm[r] ≥ ϵ*Norm,
{
i++;
x = x + α*r;
r = b - A.x;
}
]
Print["Numero di iterazioni:"];
Print
Print["Soluzione raggiunta:"]
MatrixForm[x]
Ringrazio anticipatamente tutti coloro che mi risponderanno.
Vi riporto il testo e la soluzione dell'esercizio.
Risolvere il sistema Ax=b utizizzando il metodo di richardson.
Stabilire empiricamente il miglior valore di "alfa" per la convergenza.
A{aij}ij = 1.N aij=1 i diverso da j
aij= n+1 x= 1 ossia [1,1,1...1]
bi=2n per ogni i.
SOLUZIONE
e = 0.1;
α = 2;
n = 5;
b = Table[1, {i, 1, n}];
A = Table[1 + n*KroneckerDelta[i, j], {i, 1, n}, {j, 1, n}];
x0 = Table[1.5, {i, 1, n}];
r0 = b - A.x0;
x = x0;
r = r0;
i = 0;
While[Norm[r] ≥ ϵ*Norm,
{
i++;
x = x + α*r;
r = b - A.x;
}
]
Print["Numero di iterazioni:"];
Print["Soluzione raggiunta:"]
MatrixForm[x]
Ciao. Dunque, devo confessare la mia ignoranza sul metodo di Richiardson. Conosco però un po' di Mathematica, e cerco di fare l'esegesi del tuo algoritmino.Wuque ha scritto:A{aij}ij = 1.N aij=1 i diverso da j
aij= n+1 x= 1 ossia [1,1,1...1]
bi=2n per ogni i.
Intanto, così come lo hai enunciato, il quesito non è affatto chiaro.
Suppongo che A sia la matrice che vale 1 vuori dalla diagonale e n+1 sulla diagonale principale.
b invece è un vettore con tutte le entrate uguali a 2n.
Codice: Seleziona tutto
e = 0.1;
a = 2;
n = 5;
a non ho idea di cosa sia. Potrebbe essere un "magic number" tipico del metodo? Ho il forte sospetto che qui ci sia puzza d'errore. Tornerò su questo punto poi.
Codice: Seleziona tutto
b = Table[1, {i, 1, n}];
A = Table[1 + n*KroneckerDelta[i, j], {i, 1, n}, {j, 1, n}];
?b = Table[2n, {i, 1, n}];
Codice: Seleziona tutto
x0 = Table[1.5, {i, 1, n}];
r0 = b - A.x0;
Codice: Seleziona tutto
x = x0;
r = r0;
i = 0;
Codice: Seleziona tutto
While[Norm[r] >= ?*Norm[b],
{
i++;
x = x + a*r;
r = b - A.x;
}
]
Ho però dei forti dubbi sulla convergenza di questa roba (e quindi ho il forte sospetto che l'algoritmo giri senza mai uscire).
Ti spiego perché: l'algoritmo calcola la successione
$ $\mathbf x_{i+1} = \mathbf x_i + a \mathbf r_i $
$ $\mathbf r_{i+1} = \mathbf b - A \mathbf x_{i+1} $
L'equazione di r continua a essere quella del resto e mi torna. Non capisco però la riga che aggiorna x. x viene modificato con il resto corrente, moltiplicato per un multiplo scalare [il famoso a che mi puzza].
Ora, ci sono scelte di che possono in qualche modo portarti a qualche sorta di convergenza, ma non certo un valore arbitrario di a.
Nel tuo specifico esempio, nota che hai preso un fortunatissimo valore del punto iniziale x0 e del termine noto b, dato che sono entrambi autovettori del medesimo autospazio di A. L'autovalore è n+2, quindi il tuo problema si risolve banalmente con soluzione $ \frac{1}{n+2}\mathbf b $. Se invece supponi di far girare l'algoritmo, tu hai:
$ $\mathbf x_i = \alpha_i \mathbf b $
$ $\mathbf r_i = (1 - \alpha_i(n+2)) \mathbf b $
$ $\alpha_{i+1} = \alpha_i + 2 ( 1 - \alpha_i(n+2) $.
Affinché il tutto converga, ci attenderemmo che $ \alpha_i $ converga a $ 1/(n+2) $.
Invece i valori saltellano da una parte all'altra, allontanandosi con velocità esponenziale dal punto di equilibrio. Evidentemente, c'è qualcosa che non torna.
[i:2epswnx1]già ambasciatore ufficiale di RM in Londra[/i:2epswnx1]
- - - - -
"Well, master, we're in a fix and no mistake."
- - - - -
"Well, master, we're in a fix and no mistake."
Il codice l'ho scritto io, tutte le varibili inserite sono intrinseche al metodo di Richardson. In effetti il charset di mathematica non viene correttamente riprodotto nel post.
I parametri $ \alpha, \epsilon $ vanno aggiustati in modo che il metodo converga (cosa che avevo lasciato da sistemare a Wuque).
Il vettore b ho sbagliato a digitarlo... :p
I parametri $ \alpha, \epsilon $ vanno aggiustati in modo che il metodo converga (cosa che avevo lasciato da sistemare a Wuque).
Il vettore b ho sbagliato a digitarlo... :p
Breve riassunto sul metodo di Richardson:
è un metodo iterativo più "ignorante" dei comuni Jacobi e Gauss-Seidel. L'idea è fare per l'appunto quella iterazione,
$ x_{i+1}=x_i+a (b-Ax_i) $
cioè
$ x_{i+1}=(I-aA)x_i+ab $
partendo da un punto x_0 a caso.
Per "cose note" sui metodi iterativi, il resto converge a zero (e quindi x converge alla soluzione) per ogni scelta di x_0 sse la matrice I-aA ha raggio spettrale minore di 1. Se A ha tutte le entrate positive, come nel nostro caso, basta prendere un a abbastanza piccolo e questa condizione è verificata.
ciao,
è un metodo iterativo più "ignorante" dei comuni Jacobi e Gauss-Seidel. L'idea è fare per l'appunto quella iterazione,
$ x_{i+1}=x_i+a (b-Ax_i) $
cioè
$ x_{i+1}=(I-aA)x_i+ab $
partendo da un punto x_0 a caso.
Per "cose note" sui metodi iterativi, il resto converge a zero (e quindi x converge alla soluzione) per ogni scelta di x_0 sse la matrice I-aA ha raggio spettrale minore di 1. Se A ha tutte le entrate positive, come nel nostro caso, basta prendere un a abbastanza piccolo e questa condizione è verificata.
ciao,
--federico
[tex]\frac1{\sqrt2}\bigl(\left|\text{loves me}\right\rangle+\left|\text{loves me not}\right\rangle\bigr)[/tex]
[tex]\frac1{\sqrt2}\bigl(\left|\text{loves me}\right\rangle+\left|\text{loves me not}\right\rangle\bigr)[/tex]
Catraga se solo fossi in grado lo farei...
Cmq se ho ben capito c'è un errore, qual'è la "soluzione" scritta nel modo corretto?
Me la spiegate? Io sto leggendo manuali, e il mio libro di analisi numerica, ma... più lo leggo più mi sento stupido ;_; Inoltre come fate a scrivere ai nel modo giusto??
Cercate di farmi capire vi prego!
P.S
Per questioni di lavoro risponderò ai futuri post (sperando che ci siano ) dal venerdi alla domenica.
Grazie a tutti coloro che mi daranno una mano.
Cmq se ho ben capito c'è un errore, qual'è la "soluzione" scritta nel modo corretto?
Me la spiegate? Io sto leggendo manuali, e il mio libro di analisi numerica, ma... più lo leggo più mi sento stupido ;_; Inoltre come fate a scrivere ai nel modo giusto??
Cercate di farmi capire vi prego!
P.S
Per questioni di lavoro risponderò ai futuri post (sperando che ci siano ) dal venerdi alla domenica.
Grazie a tutti coloro che mi daranno una mano.
Re: Esercizio analisi numerica
Non ho capito esattamente di cosa hai bisogno ma ci provo...
Esegue l'iterazione fintantoché la norma del "resto" (cioè la discrepanza tra Ax e b) non è abbastanza "piccola". Qui credo che si usi il valore "e" che hai inizializzato all'inizio, al posto del carattere sgorbio. Nota che visto che b potrebbe contenere entrate molto grandi, il valore soglia viene moltiplicato per la norma di b.
Questo è il metodo di Richardson. Se non l'hai mai visto, dovresti riguardarlo su un buon libro di analisi numerica.
Scrive il valore di i (numero di iterazioni) e quello di x (la stima corrente della soluzione)
spero che ora sia ok?
ciao,
Inizializza qualche variabile: n è la dimensione della matrice, e credo che sia una specie di "soglia epsilon" anche se non lo ritrovo dopo, a è il parametro che compare nel metodo di RichardsonWuque ha scritto: e = 0.1;
α = 2;
n = 5;
Inizializza matrici e vettori A,b e x_0 come dalle ipotesi del problema. KroneckerDelta è un nome altisonante per la funzione che ritorna 1 se i suoi due argomenti i e j sono uguali, zero altrimenti.b = Table[1, {i, 1, n}];
A = Table[1 + n*KroneckerDelta[i, j], {i, 1, n}, {j, 1, n}];
x0 = Table[1.5, {i, 1, n}];
Inizializza un altro po' di roba. x e r sono la stima di x e del "resto" b-Ax che verranno updatate durante l'algoritmo. r_0 è pressoché inutile, perché la prima e la terza riga qui sopra potevano essere condensate in una. i è un contatore che conta il numero di iterazioni fatte (difatti a ogni iterazione qui sotto c'è l'istruzione i++, cioè i=i+1.r0 = b - A.x0;
x = x0;
r = r0;
i = 0;
While[Norm[r] ≥ ϵ*Norm,
Esegue l'iterazione fintantoché la norma del "resto" (cioè la discrepanza tra Ax e b) non è abbastanza "piccola". Qui credo che si usi il valore "e" che hai inizializzato all'inizio, al posto del carattere sgorbio. Nota che visto che b potrebbe contenere entrate molto grandi, il valore soglia viene moltiplicato per la norma di b.
{
i++;
x = x + α*r;
r = b - A.x;
}
}
Questo è il metodo di Richardson. Se non l'hai mai visto, dovresti riguardarlo su un buon libro di analisi numerica.
Print["Numero di iterazioni:"];
Print["Soluzione raggiunta:"]
MatrixForm[x]
Scrive il valore di i (numero di iterazioni) e quello di x (la stima corrente della soluzione)
spero che ora sia ok?
ciao,
--federico
[tex]\frac1{\sqrt2}\bigl(\left|\text{loves me}\right\rangle+\left|\text{loves me not}\right\rangle\bigr)[/tex]
[tex]\frac1{\sqrt2}\bigl(\left|\text{loves me}\right\rangle+\left|\text{loves me not}\right\rangle\bigr)[/tex]
Grazie 1000 dell'aiuto fph, qui sotto ti riporo le domande di cosa non mi è ancora chiaro.
Mi potresti spiegare la sintassi di tale scrittura b = Table[1, {i, 1, n}]?
KroneckerDelta è il nome della funzione, ok, ma come fa a ritornare 1 se i e j dal testo dell'esercizio vengono dichiarati diversi?
Unendo la prima e la terza cosa avrei dovuto scrivere? r = r0 = b-A.x0 ?
Cosa si intende per discrepanza tra Ax e b?
Cos'è la norma di b?
MatrixForm[x] è il comando di mathematica che fa apparire la matrice x?
Infine quando catraga scirsse:
Cosa intendeva?
Grazie ancora per l'aiuto
Le matrici sono del tipo [1,1,1..] ?Citazione:
b = Table[1, {i, 1, n}];
A = Table[1 + n*KroneckerDelta[i, j], {i, 1, n}, {j, 1, n}];
x0 = Table[1.5, {i, 1, n}];
Inizializza matrici e vettori A,b e x_0 come dalle ipotesi del problema. KroneckerDelta è un nome altisonante per la funzione che ritorna 1 se i suoi due argomenti i e j sono uguali, zero altrimenti.
Mi potresti spiegare la sintassi di tale scrittura b = Table[1, {i, 1, n}]?
KroneckerDelta è il nome della funzione, ok, ma come fa a ritornare 1 se i e j dal testo dell'esercizio vengono dichiarati diversi?
Cosa si intende per stima e resto?Citazione:
r0 = b - A.x0;
x = x0;
r = r0;
i = 0;
Inizializza un altro po' di roba. x e r sono la stima di x e del "resto" b-Ax che verranno updatate durante l'algoritmo. r_0 è pressoché inutile, perché la prima e la terza riga qui sopra potevano essere condensate in una. i è un contatore che conta il numero di iterazioni fatte (difatti a ogni iterazione qui sotto c'è l'istruzione i++, cioè i=i+1.
Unendo la prima e la terza cosa avrei dovuto scrivere? r = r0 = b-A.x0 ?
Citazione:
While[Norm[r] ≥ ϵ*Norm,
Esegue l'iterazione fintantoché la norma del "resto" (cioè la discrepanza tra Ax e b) non è abbastanza "piccola". Qui credo che si usi il valore "e" che hai inizializzato all'inizio, al posto del carattere sgorbio. Nota che visto che b potrebbe contenere entrate molto grandi, il valore soglia viene moltiplicato per la norma di b.
Cosa si intende per discrepanza tra Ax e b?
Cos'è la norma di b?
Citazione:
Print["Numero di iterazioni:"];
Print["Soluzione raggiunta:"]
MatrixForm[x]
Scrive il valore di i (numero di iterazioni) e quello di x (la stima corrente della soluzione)
MatrixForm[x] è il comando di mathematica che fa apparire la matrice x?
Infine quando catraga scirsse:
I parametri vanno aggiustati in modo che il metodo converga.
Il vettore b ho sbagliato a digitarlo... :p
Cosa intendeva?
Grazie ancora per l'aiuto
Non conosco mathematica, ma mi sembra di capire che "Table[f(i), {i, 1, n}]" restituisce un vettore lungo n, che all'entrata i-esima (per i=1..n) contiene il valore f(i).Wuque ha scritto: Le matrici sono del tipo [1,1,1..] ?
Mi potresti spiegare la sintassi di tale scrittura b = Table[1, {i, 1, n}]?
KroneckerDelta è il nome della funzione, ok, ma come fa a ritornare 1 se i e j dal testo dell'esercizio vengono dichiarati diversi?
Quindi b=[1 ... 1 ], A=E+nI (dove I è l'identità e E è la matrice con tutti uni) e x0=[1.5 1.5 ... 1.5].
Nel caso della dichiarazione di A, i e j variano a seconda delle entrate della tabella (=matrice) nel modo in cui sei abituato, quindi a volte KroneckerDelta restituisce 1 e a volte 0.
Basta r=b-A.x0, r0 poi non lo usi più. Stima=vuoi approssimare la soluzione, x è la tua approssimazione corrente. Il resto è b-A.x, che dovrebbe essere uguale a 0 se x è davvero una soluzione.Cosa si intende per stima e resto?
Unendo la prima e la terza cosa avrei dovuto scrivere? r = r0 = b-A.x0 ?
La norma del vettore resto b-Ax. Questa non è terminologia standard, è la prima parola che mi è venuta.Cosa si intende per discrepanza tra Ax e b?
ehm... dovresti veramente ripassare qualcosa mi sa.Cos'è la norma di b?
Esistono diverse "norme" per un vettore, quella che si usa di solito è la cosiddetta norma-2, cioè
$ \left\Vert v\right\Vert_2=\left(\sum_{i=1}^n v_i^2\right)^{1/2} $
Boh, non conosco mathematica ma ad occhio direi di sì. Guarda sull'help di mathematica magari.MatrixForm[x] è il comando di mathematica che fa apparire la matrice x?
(Domanda: ma il codice che stiamo commentando l'hai scritto tu o lui?)Infine quando catraga scirsse:Cosa intendeva?I parametri vanno aggiustati in modo che il metodo converga.
Il vettore b ho sbagliato a digitarlo... :p
Che nel testo dell'esercizio b=[2n 2n ... 2n], mentre nel programma si usa un altro vettore, b=[1 1 ... 1]. Probabilmente una distrazione. (tra l'altro nota che basta dividere tutto per 2n e si riconduce una soluzione all'altra).
--federico
[tex]\frac1{\sqrt2}\bigl(\left|\text{loves me}\right\rangle+\left|\text{loves me not}\right\rangle\bigr)[/tex]
[tex]\frac1{\sqrt2}\bigl(\left|\text{loves me}\right\rangle+\left|\text{loves me not}\right\rangle\bigr)[/tex]
Il codice lo ha scritto catraga, la formula A=E+nI da dove l'hai presa?è quella generica? "E" Nel nostro caso sarebbe b giusto? mentre l'identità qual'è?
Se ho ben capito per rimediare alla diamenticanza di catraga dobrei dichiarare b in questo modo b= Table[2n, {i, 1, n}] giusto?
P.S
Qualcuno mi sa spiegare il perchè scrivendo questo codice su mathematica e cliccando shif+invio non succede nulla??
Fph Grazie molte del tuo aiuto mi sei moooolto Utile!
Se ho ben capito per rimediare alla diamenticanza di catraga dobrei dichiarare b in questo modo b= Table[2n, {i, 1, n}] giusto?
P.S
Qualcuno mi sa spiegare il perchè scrivendo questo codice su mathematica e cliccando shif+invio non succede nulla??
Fph Grazie molte del tuo aiuto mi sei moooolto Utile!
E è la matrice con tutte le entrate uguali a 1, E_{ij}=1 per tutti gli i e j. I è la matrice identità (1 sulla diagonale, 0 altrove). Quindi E+nI, se ci pensi bene, è proprio uguale alla matrice A che serve a te (guarda il testo dell'esercizio...).Wuque ha scritto:la formula A=E+nI da dove l'hai presa?è quella generica? "E" Nel nostro caso sarebbe b giusto? mentre l'identità qual'è?
Il codice "corretto" per inizializzare b è quello giusto.
ciao,
--federico
[tex]\frac1{\sqrt2}\bigl(\left|\text{loves me}\right\rangle+\left|\text{loves me not}\right\rangle\bigr)[/tex]
[tex]\frac1{\sqrt2}\bigl(\left|\text{loves me}\right\rangle+\left|\text{loves me not}\right\rangle\bigr)[/tex]
OK, un'ultima domanda, ma perchè su Mathematica 5.2 caricando questo esercizio e cliccando shift+invio non succede nulla?? non mi dovrebbe comparire il numero delle iterazioni e la matrice?? So che fph non fa uso di mathematica, quindi se qualcun'altro può rispondermi gli sarei grato.
Può dipendere dalla versione non ancora registrata?
Grazie
Può dipendere dalla versione non ancora registrata?
Grazie
Mi hanno detto che il motivo per cui questo progetto non converge dipende da Alpha che non e' settato correttmante per la risposta, e devo fare dei tentativi, che rientrano anche tra i "to do" dell'esercizio.
Io ho provato con 1, 0.5 e 0.1, ma non mi restituisce nulla, con 0.1 mi da errore. Dov'è che sbaglio? Cosa si intende per "rientra nei "to do" dell'esercizio? Quale ragionamento dovrei seguire per capire in che modo devo settare correttamente alfa?
Io ho provato con 1, 0.5 e 0.1, ma non mi restituisce nulla, con 0.1 mi da errore. Dov'è che sbaglio? Cosa si intende per "rientra nei "to do" dell'esercizio? Quale ragionamento dovrei seguire per capire in che modo devo settare correttamente alfa?