Appunti Metodi numerici - Parte 3

%% integrazione nel tempo con metodo di Eulero

x(1)=x0; y(1)=y0; Teta(1)=Teta0;

hold on

fige=figure;imshow('noccino.jpg'); set(gca,'YDir','normal')

axis equal

set(fige,'KeyPressFcn',@keypress)

for campo1=percorso

t=0:dt:4

x(int1,:) =

y(int1,:) =

Teta(t+1)=Teta(t)+dt+Vt*Tau(phikl);

carpo5 (x(rt1),y(rt4),37)

carro7ogx=([x(t+1)+100*cos(Teta(t+1)),y(t+1)+100*sin(Teta(t1)),1])

plot(...)

Title([' t=',num2str2(t(t+1)),' ph=',num2str2(180/pi*phi1)])

drawnow;pause(.01)

end

function :

end

EQ. DIFF. ORDINARIE AI VALORI DI CONTORNO

17/11/2022

Fino ad ora abbiamo studiato:

dy/dt = -a.y

y(t=0)=yo

y(t=s.2)=28

Queste ode descrivono una geometria di curve

la cond. iniziale individua una posizione all'interno di una famiglia

Oppure abbiamo visto:

dxm/dt ... →descrivono lo stato del sistema

dym/dt ...

le cond. iniziali sono:

x(0)=x₀

y(0)=v₀

Se specifico le cond. iniziali per 2 istanti diversi del tempo onde t₁ e t₂?

x(0)=x₀, x(5,2)=v₀

Non mi permette di identificare lo stato del sistema dal quale partire

Ragioniamo sul seguente problema:

Immaginiamo di caricare la Trave e vogliamo studiare le deformazioni

Sotto il carico verticale la Trave si inflette. Ma non abbiamo una risposta statica perchè la trave inizia ad oscillare

Quindi il sistema fisico che stiamo studiando dipende sia dalla spazio che dal tempo. Preso preso le oscillazioni si smorzano arrivando ad una situazione statica.

Quindi:

M(x,t) → def. verticale della trave

PDE → problema della derivata parziale

Spazio del sistema all'istante tf → Ecco calcolius

Se c'è una funzione posso partire da T_b e calcolarmi il valore esatto della q_a, con l'intersezione delle rette trovate.

I problemi ai valori di contorno ODE-BV si possono risolvere con 3 metodi diversi.

1. Metodo a tentativi (non ci sono esercizi dell'esame)
2. Metodi diretti/impliciti
3. Metodi iterativi

Metodi Diretti

Abbiamo sviluppato il seguente sistema:

dq/dx = 2h/z(T_∞-T)

dT/dx = -q/k

• Sistema di 2 eq. in 2 incognite (T, q), ODE lineari.

Per tali metodi dobbiamo usare la versione del II ordine, bisogna scrivere di grado.

d²T/dx² = -1/k dq/dx

d²T/dx² = -1/k * 2h/z (T_∞-T)

d²T/dx² + (2h/kz)(T_∞-T) = 0

• -> Eq. lineare del II ordine

avere:

d²T/dx² - h'T = -h'T_∞

-> decodificare questa eq. differenziale

k(1,1)=1; B(1)=Ta;

k(Nx,Nx)=1; B(Nx)=Tb;

alfa=-2*hp*dx^2;

for ix=2:Nx-1

k(ix,ix)=alfa;

k(ix,ix-1)=1;

k(ix,ix+1)=1;

B(ix)=-dx^2*hp*TiuB;

end

%% Risoluzione per T

T=k\B;

(oppure T=inv(k)*B;)

end

Allora:

5 · T₄ + αT₅ + T₄ = -h' Δx² T_∞

Raccolgo i coefficienti delle incognite

5 · (2)T₄ + αT₅ = -h' Δx² T_∞

Il sistema è fatto nel seguente modo:

1. αT₁ + 2T₂ = -h' Δx² T_∞
2. T₁ + αT₂ + T₃ = -h' Δx² T_∞
3. T₂ + αT₃ + T₄ = -h' Δx² T_∞
4. T₃ + αT₄ + T₅ = -h' Δx² T_∞
5. 2T₄ + αT₅ = -h' Δx² T_∞

Posso riscrivere il sistema in forma matriciale:

X = T₁T₂T₃T₄T₅→ K X = B

K α20001α10001α10001α10002αX T₁T₂T₃T₄T₅= [-h' Δx² T_∞::::]

Quindi:

Condizioni al contorno: T₀ e T₆ sono calcolate in base alle BC Neumann

for ix = 2: Nx - 1

T(ix) = 1/α (-Told(ix-1) - Told(ix+1) - hp*dx^2*Tiuf)

* end

Told=T;

end

Test di convergenza(non viene chiesto spesso)Sì continuo nelle iterazioniwhile differenza > ε

NO mi fermo

MATLAB

ODE-BV-04

function [x, T] = Termico (l, dx, hp, Tiuf, Ta, Tb, mie)

x = 0 : dx : l; x = x'; Nx = length(x);

Told = zeros(Nx, 1); T = zeros(Nx, 1);

Told(1) = Ta; Told(Nx) = Tb;

T(1) = Ta; T(Nx) = Tb;

alfa = -2 - dx^2 * hp;

for T = 1; mie;

for ix = 2: Nx -1

T(ix) = 1/α * (-Told(ix-1) - Told(ix+1) - hp*dx^2*Tiuf);

end

Told = T;

plot(x, T);

pause(0.1)

end

Ci sono dei metodi che velocizzano queste operazioni?

Sì, sono il metodo di Gauss-Seidel e il metodo SOR (successive over relaxation)