Metodo elementi finiti diffusione-reazione (reazione dominante)
Buongiorno a tutti!
Vorrei scrivere un programma in matlab per il seguente problema (reazione dominante $\sigma >> 1$):
$ { ( -u''(x)+\sigmau(x) =0 \quad\quad 0
In particolare devo fare un'approssimazione ad elementi finiti standard, con il metodo di Galerkin.
Ovvero, devo trovare $u_h \in V_h$ (dipendende dal passo della mesh $h$)
tale che la formulazione debole $a(u_h,v_h)=F(v_h)$ $ (1)$
valga per ogni $v_h \in V_h$ (spazi localmente di polinomi).
Ho letto che se:
$K_(ij)= \int _(0)^(1) varphi'_j varphi'_i dx $,
$M_(ij)= \int _(0)^(1) varphi_j varphi_i dx $,
$f_i=F(varphi_i)$
$i,j=1,2,..N$
allora $(1)$ è equivalente alla forma algebrica
$(K+ \sigma M)U=F$. $(2)$
Per cui risolvere il sistema lineare in funzione di $U$ mi darebbe i coefficienti della mia approssimante.
Però osservo che qualcosa di concettuale mi manca:
ho trovato questo approccio fatto anche nel mio caso quando la funzione $f(x)$(forzante) non c'è, perché ho un'omogenea.
Non dovrei avere a destra di $(2)$ il vettore nullo?
Vorrei scrivere un programma in matlab per il seguente problema (reazione dominante $\sigma >> 1$):
$ { ( -u''(x)+\sigmau(x) =0 \quad\quad 0
In particolare devo fare un'approssimazione ad elementi finiti standard, con il metodo di Galerkin.
Ovvero, devo trovare $u_h \in V_h$ (dipendende dal passo della mesh $h$)
tale che la formulazione debole $a(u_h,v_h)=F(v_h)$ $ (1)$
valga per ogni $v_h \in V_h$ (spazi localmente di polinomi).
Ho letto che se:
$K_(ij)= \int _(0)^(1) varphi'_j varphi'_i dx $,
$M_(ij)= \int _(0)^(1) varphi_j varphi_i dx $,
$f_i=F(varphi_i)$
$i,j=1,2,..N$
allora $(1)$ è equivalente alla forma algebrica
$(K+ \sigma M)U=F$. $(2)$
Per cui risolvere il sistema lineare in funzione di $U$ mi darebbe i coefficienti della mia approssimante.
Però osservo che qualcosa di concettuale mi manca:
ho trovato questo approccio fatto anche nel mio caso quando la funzione $f(x)$(forzante) non c'è, perché ho un'omogenea.
Non dovrei avere a destra di $(2)$ il vettore nullo?
Risposte
Stai dimenticando di imporre le condizioni al bordo.
Si è vero così non ne sto tenendo conto delle informazioni al bordo.
Ma per inserirle, devo per esempio risolvere il sistema lineare, togliendo la prima e l'ultima riga/colonna e inserire i dati al bordo alla soluzione? Nel caso specifico $u_h=[0;u_h;1];$
Altrimenti da qualche parte ho sentito parlare di una cosiddetta funzione di rilevamento $G$, esempio $G(x)=x$ e risolvere con $ hat(u)(x) = u(x)- G(x) $, in modo da avere $hat(u)$ nulla al bordo. Non mi è molto chiaro..
Ma per inserirle, devo per esempio risolvere il sistema lineare, togliendo la prima e l'ultima riga/colonna e inserire i dati al bordo alla soluzione? Nel caso specifico $u_h=[0;u_h;1];$
Altrimenti da qualche parte ho sentito parlare di una cosiddetta funzione di rilevamento $G$, esempio $G(x)=x$ e risolvere con $ hat(u)(x) = u(x)- G(x) $, in modo da avere $hat(u)$ nulla al bordo. Non mi è molto chiaro..
No, devi modificare il sistema lineare, prima di risolverlo, in modo che la soluzione rispetti le condizioni al bordo.
Il rilevamento è un altro modo: se fai quella sostituzione nell'equazione ti salta fuori un termine noto aggiuntivo in cambio delle condizioni al bordo omogenee, ma per problemi concreti non ne vale la pena, meglio modificare il sistema lineare.
Il rilevamento è essenziale nell'analisi dell'equazione, perché ti permette di usare il teorema di Lax-Milgram.
Il rilevamento è un altro modo: se fai quella sostituzione nell'equazione ti salta fuori un termine noto aggiuntivo in cambio delle condizioni al bordo omogenee, ma per problemi concreti non ne vale la pena, meglio modificare il sistema lineare.
Il rilevamento è essenziale nell'analisi dell'equazione, perché ti permette di usare il teorema di Lax-Milgram.
forse l'ultima equazione del sistema deve essere $1$
$(int_(x_N)^(x_(N+1)) (\varphi'_(N+1))^2 + \sigma int_(x_N)^(x_(N+1)) (\varphi_(N+1))^2)u_(N+1)=1 $
Avendo così il vettore $F$ tutti gli elementi nulli, tranne l'ultimo, uguale ad $1$.
Sbaglio?
$(int_(x_N)^(x_(N+1)) (\varphi'_(N+1))^2 + \sigma int_(x_N)^(x_(N+1)) (\varphi_(N+1))^2)u_(N+1)=1 $
Avendo così il vettore $F$ tutti gli elementi nulli, tranne l'ultimo, uguale ad $1$.
Sbaglio?
Sbagli: l'ultima equazione del sistema dev'essere \(u_{N+1} = 1\), che è la condizione al bordo, così il vettore F è, sì, quello che hai scritto tu, ma anche la corrispondente riga della matrice dev'essere modificata.
Accedi a tutti gli appunti
Tutor AI: studia meglio e in meno tempo