sistema tridiagonale

[maxithron]

dovrei risolvere un sistema tridiagonale.

Suggerimenti per il linguaggio + consono da utilizzare?http://forum.hwupgrade.it/faccine/23.gif

[cionci]

Se ti va bene: http://www.prato.linux.it/~fcasadei/math/calcolo/html/node38.html

In generale qualsiasi linguaggio che permette di gestire i vettori...ad esempio quello sopra dovrebbe essere Matlab...

[maxithron]

che può andar bene il link :)

anche se con google sto trovando esempi con quasi tutti i linguaggi.....ora:

Se cerco di risolvere con un solo linguaggio, mi pentirò per non aver curiosato anche con gli altri ma......se comincio a curiosare anche con gli altri....mi pentirò perchè non riuscirò a risolvere nemmeno con uno http://forum.hwupgrade.it/faccine/19.gif

[a2000]

Originally posted by "maxithron"

dovrei risolvere un sistema tridiagonale.

Suggerimenti per il linguaggio + consono da utilizzare?http://forum.hwupgrade.it/faccine/23.gif

da volumi finiti 2D o ... colonne a stadi ? :D

è un classico: algoritmo di Thomas.

il Fortran è la morte sua (che vuol dire 7 righe). :cool:

[maxithron]

vramente dovrei utilizzare l'algoritmo di eliminazione di gauss per risolvere:

classico caso Ax = B e la matrice A è tridiagonale.

Dici che il fortran è il miglior metodo?

[a2000]

Nel caso di matrici tridiagonali il metodo più efficiente è quello di Thomas:
due passate una in avanti e una all'indietro e hai risolto.

il linguaggio più adatto è il Fortran poi trovi i listati (si fa per dire: < 10 righe) già bell'e pronti, oltre che in Fortran, anche in Pascal, Basic, C, forse anche Java :p

[maxithron]

Originally posted by "a2000"


il linguaggio più adatto è il Fortran poi trovi i listati (si fa per dire: < 10 righe) già bell'e pronti, oltre che in Fortran, anche in Pascal, Basic, C, forse anche Java :p

Nooooo....non voglio trovare i listati già fatti! Sono un masochista convinto :muro:

Ehm..cmq... sette righe dici?

:p

[a2000]

9 con Subroutine e End Subroutine :D :D

se vuoi te lo passo :)

[maxithron]

Originally posted by "a2000"


se vuoi te lo passo :)

Non c'è sfizio :o

hai detto proprio 9?

no perchè io http://forum.hwupgrade.it/faccine/26.gif

[a2000]

con Thomas sì ....

con Gauss ... 15 :D :D

beh però con la ricerca pivot ... 20 !

http://forum.hwupgrade.it/faccine/44.gif


te lo passo ???? http://forum.hwupgrade.it/faccine/41.gif

[maxithron]

Originally posted by "a2000"



te lo passo ???? http://forum.hwupgrade.it/faccine/41.gif

No......No......No........


Com'è? :D

[a2000]

ma da dove nasce il sistema ?

no, perchè se è a diagonale dominante è garantita la convergenza anche con metodi iterativi: 5 righe :D :D

[a2000]

Originally posted by "maxithron"



Com'è? :D

pezzettino (cioè la metà di andata) :p


Do i = 2, II
DEN = P(i) - W(i) * EE(i - 1)
EE(i) = E(i) / DEN
BB(i) = (B(i) + W(i) * BB(i - 1)) / DEN
End Do

[maxithron]

Grrrrrrrr!!!!! http://forum.hwupgrade.it/faccine/40.gif

5 righe????

[maxithron]

ma non è così?:

program main
real, dimension (10):: a,d,c,b,x
integer ::n = 10

interface
subroutine tri(n,a,d,c,b,x)
integer, intent(in)::n
real, dimension(:), intent(in):: a,d,c,b
real, dimension(:), intent(out)::x
end subroutine tri
end interface
do i=1,n
d(i) = 2.0
a(i) = 0.5
c(i) = 0.5
b(i) = 3.0
end do
b(1) = 2.5
b(n) = 2.5
call tri(n,a,d,c,b,x)
print*, "subroutine tri"
print"(f22.14)",(x(i),i=1,n)
do i=1,n
d(i) = 2.0
c(i) = 0.5
b(i) = 3.0
end do
b(1) = 2.5
b(n) = 2.5
call tri(n,c,d,c,b,x)
print*, "subroutine tri again"
print "(f22.14)",(x(i),i=1,n)
end program main

subroutine tri(n,a,d,c,b,x)
integer, intent(in)::n
real, dimension(:), intent(in)::a,d,c,b
real, dimension(:), intent(out):: x
integer ::i
real :: xmult
do i = 2,n
xmult = a(i-1)/d(i-1)
d(i) = d(i) - xmult*c(i-1)
b(i) = b(i) - xmult*b(i-1)
end do
x(n) = b(n)/d(n)
do i = n-1,1,-1
x(i) = (b(i) - c(i)*x(i+1))/d(i)
end do
end subroutine tri

[a2000]

5 + 4 = 9 :D :D

ah, ma ca++o non mi dai elementi ...

da dove viene sto' sistema o è un esercizio di calcolo numerico

cosa studi ?

[a2000]

Originally posted by "maxithron"

ma non è coì?:
....
subroutine tri(n,a,d,c,b,x)
integer, intent(in)::n
real, dimension(:), intent(in)::a,d,c,b
real, dimension(:), intent(out):: x
integer ::i
real :: xmult
do i = 2,n
xmult = a(i-1)/d(i-1)
d(i) = d(i) - xmult*c(i-1)
b(i) = b(i) - xmult*b(i-1)
end do
x(n) = b(n)/d(n)
do i = n-1,1,-1
x(i) = (b(i) - c(i)*x(i+1))/d(i)
end do
end subroutine tri[/size]

sì bravo, con qualche ridondanza, ma direi che ci siamo ! :p

[maxithron]

E' un esercizio di calcolo numerico.

Cosa studio? :eek:

Quello che mi capita... :p

Scherzi a parte non sono uno studente universitario ma sono divorato da una strana passione schizo-frenetico-psicotica per la matematica e la programmazione :cool:

[a2000]

e figa ?

:D :D :D

[maxithron]

ca@@o ci fanno quelle faccine nel list?

[a2000]

no scherzo, bravo: mathematics is fun !

ci sono anche delle gare di matematica e siti dedicati

ma dopo che vuoi fare, matematica, ingegneria ...

[maxithron]

fermo....finiamola di scrivere in contemporanea :D

per quanto riguarda l'altra materia....visto che in qualche 3d fa si parlava di una cena.....

potrei umilmente offrirmi come PR di fighe!!!
:rolleyes: ;) :) :p :D :pig: http://forum.hwupgrade.it/faccine/44.gif

[maxithron]

Originally posted by "a2000"



ma dopo che vuoi fare, matematica, ingegneria ...

Mi sa che ho passato l'età per iscriivermi all'uni. :cry:

Ho una mia piccola azienda con la quale produco appunto software su richiesta con i linguaggi che (se non ho letto male) nei 3d di questa sezione tu amorevolmente aborri!!!

[a2000]

Originally posted by "maxithron"


...
potrei umilmente offrirmi come PR di fighe!!!
:rolleyes: ;) :) :p :D :pig: http://forum.hwupgrade.it/faccine/44.gif

aaaahhhhhh, mi fai morire :D :D
ma vuoi offrirti come PR di fi++e o come prrrrrr di fi++e :p :p :p

spurcazoun ! :pig: :pig:

[a2000]

Originally posted by "maxithron"



Mi sa che ho passato l'età per iscriivermi all'uni. :cry:
Ho una mia piccola azienda con la quale produco appunto software su richiesta con i linguaggi che (se non ho letto male) nei 3d di questa sezione tu amorevolmente aborri!!!

ah allora fai anche il prrrrr di segretarie :p :D :D

basta che sennò ci cacciano per contenuti osceni ! :pig:

[a2000]

e per quello che ti interessi di sistemi tridiagoanali ! :D :D http://forum.hwupgrade.it/faccine/44.gif

[a2000]

maxithron mi hai rovinato !

[maxithron]

Originally posted by "a2000"

maxithron mi hai rovinato !

io....e perchè? :confused:

[maxithron]

cmq....considera che in FORmula TRANsexuals non sono molto ferrato però devo risolvere sta benedetta tridiagonale.

Davvero con l'algoritmo di Thomas sono 9 righe? Perchè alla fine non le hai + postate!

[maxithron]

Forse sono queste 9 righe?


PROGRAM Testtrd
PARAMETER (nmax=10000)
INTEGER i,n
REAL*8 aval,bval,cval,a(nmax),b(nmax),c(nmax),d(nmax),x(nmax)

2 continue

write(*,*) 'Inserire l''ordine del sistema:n'
read (*,*) n
write(*,*) ' '

aval=1
bval=-4
cval=1

DO 10 i=1,n
a(i)=aval
b(i)=bval
c(i)=cval
IF (i.eq.1) THEN
d(i)=b(i)+c(i)
ELSE IF (i.eq.n) THEN
d(i)=a(i)+b(i)
ELSE
d(i)=a(i)+b(i)+c(i)
END IF
10 continue

call tridiag(n,a,b,c,d,x)

write(*,*) 'SOLUZIONE'
write(*,*) ' '
write(*,*) (x(i),i=1,n)
GOTO 2
END

un prezioso contributo da: http://www.detec.unina.it/meola2/htm/algoritmo_di_thomas.htm

[a2000]

Originally posted by "maxithron"

cmq....considera che in FORmula TRANsexuals non sono molto ferrato però devo risolvere sta benedetta tridiagonale.

Davvero con l'algoritmo di Thomas sono 9 righe? Perchè alla fine non le hai + postate!

FORmula TRANsexuals : ma tu hai deciso di accidermi stasera ! :D

e prima mi hai detto che volevi fare solo onanismo !



Sub Thomas(P, E, W, B, X, II)
ReDim EE(0 To II + 1), BB(0 To II + 1)

For i = 1 To II
DEN = P(i) - W(i) * EE(i - 1)
EE(i) = E(i) / DEN
BB(i) = (B(i) - W(i) * BB(i - 1)) / DEN
Next i

X(II) = BB(II)
For i = II - 1 To 1 Step -1
X(i) = BB(i) - EE(i) * X(i + 1)
Next i

End Sub


non tastato !

:D

[maxithron]

sono ancora sveglio e sto ammattendo....


Sub Thomas(P, E, W, B, X, II)
ReDim EE(0 To II + 1), BB(0 To II + 1)

For i = 1 To II
DEN = P(i) - W(i) * EE(i - 1)
EE(i) = E(i) / DEN
BB(i) = (B(i) - W(i) * BB(i - 1)) / DEN
Next i

X(II) = BB(II)
For i = II - 1 To 1 Step -1
X(i) = BB(i) - EE(i) * X(i + 1)
Next i

End Sub

d'annata....sintesi

[maxithron]

potrebbe andare anche così:

program main
real, dimension (10):: a,d,c,b,x
integer ::n = 10

interface
subroutine tri(n,a,d,c,b,x)
integer, intent(in)::n
real, dimension(:), intent(in):: a,d,c,b
real, dimension(:), intent(out)::x
end subroutine tri
end interface
do i=1,n
d(i) = 2.0
a(i) = 0.5
c(i) = 0.5
b(i) = 3.0
end do
b(1) = 2.5
b(n) = 2.5
call tri(n,a,d,c,b,x)
print*, "subroutine tri"
print"(f22.14)",(x(i),i=1,n)
do i=1,n
d(i) = 2.0
c(i) = 0.5
b(i) = 3.0
end do
b(1) = 2.5
b(n) = 2.5
call tri(n,c,d,c,b,x)
print*, "subroutine tri again"
print "(f22.14)",(x(i),i=1,n)
end program main

subroutine tri(n,a,d,c,b,x)
integer, intent(in)::n
real, dimension(:), intent(in)::a,d,c,b
real, dimension(:), intent(out):: x
integer ::i
real :: xmult
do i = 2,n
xmult = a(i-1)/d(i-1)
d(i) = d(i) - xmult*c(i-1)
b(i) = b(i) - xmult*b(i-1)
end do
x(n) = b(n)/d(n)
do i = n-1,1,-1
x(i) = (b(i) - c(i)*x(i+1))/d(i)
end do
end subroutine tri

[a2000]

sì, sì anche, bravo (verifica i residui però).

comunque una cosa che puoi notare, anche in questo caso, è che ti sono servite più righe di codice per l'input-output dati che per il calcolo.

E questo è, ed è stato, un grande e grave problema.
In generale è assurdo pensare di inputare dati su prompt del codice.
Più sempllice ma noioso, e solo per calcoli di prova, l'input dati da codice.

E' indispensabile utilizzare un ambiente potente e flessibile per l'input-output dati che non può essere un linguaggio tal quale neanche nelle versioni Visual..... (1 riga di codice calcolo, 20 righe di codice ausiliario)

Excel risponde bene alle necessità dell'input-output dati in forma numerica e anche grafica con molti vantaggi per il programmatore e per l'utilizzatore, che si riconducono alla facilità di inputare centinaia di dati senza introdurli da tastiera (ma calcolandoli sul foglio di calcolo o da procedure ausiliarie definite dall'utente).

Qualcosa di meglio si dovrebbe avere per l'input-output grafico 3D.

[maxithron]

Originally posted by "a2000"

sì, sì anche, bravo (verifica i residui però).

comunque una cosa che puoi notare, anche in questo caso, è che ti sono servite più righe di codice per l'input-output dati che per il calcolo.



mi piace vedere lo schermo pieno di robaccia :D

ed inoltre, questi listati mi servono per vincere un concorso come dattilografohttp://forum.hwupgrade.it/faccine/44.gif

cmq... quelle righe che hai postato prima ....con quale "criterio" di sintesi le hai generate?

[a2000]

La subroutine che ho postato è l'implementazione dell'algoritmo di Thomas che si ricava supponendo di potere esprimere la soluzione del sistema tridiagonale di equazioni lineari nella forma:
Xi[/siz]

sostituisci nel sistema e ricavi le espressioni di Pi[/siz]).

Sistemi eptadiagonali scaturiscono dalla discretizzazione a volumi finiti delle equazioni di bilancio di quantità di moto, massa ed energia per il moto dei fluidi.
Un metodo famoso per la discretizzazione e linearizzazione di questo set di equazioni differenziali è quello di Patankar e Spalding.
Il metodo di soluzione più semplice è l'applicazione iterativa dell'algoritmo di Thomas al sistema eptadiagonale risultante.

Se ti piace la matematica e la fisica ti puoi cimentare:
S.V. Patankar - Numerical Heat Transfer and Fluid Flow - McGraw-Hill

alla fine della lettura potrai determinare punto per punto ed istante per istante la velocità, pressione, temperatura, composizione dell'aria nel case del tuo computer o nelle stanze della tua casa o intorno alla tua casa, al tuo quartiere, alla tua città, ... fino all'intero globo terracqueo ! :eek:

(poi c'è qualche problemino :sofico: di unicità e stabilità della soluzione ma tu ... per il momento ... fottitene ! :D )

:)

[maxithron]

Originally posted by "a2000"

alla fine della lettura potrai determinare punto per punto ed istante per istante la velocità, pressione, temperatura, composizione dell'aria nel case del tuo computer o nelle stanze della tua casa o intorno alla tua casa, al tuo quartiere, alla tua città, ... fino all'intero globo terracqueo ! :eek:



e se mi accorgo che la temperatura è troppo delle suddette cose è troppo alta, che dissipatore mi consigli?
:p

Ti ringrazio per il suggerimento del libro.Vedrò di comprarlo martedì.

Grazie!

[a2000]

per il case lo scoperchiamento

per le stanze della tua casa una bella bionda naturale

per il tuo quartiere delle tende a strisce arcobaleno

per la tua città un'astronave di quelle alla indipendence-day

per il globo terracqueo niente: ma lasciali anna' a mori' ammazzati sti' 5 miliardi de stronzi :D
(scherzo ... scherzo ? :confused: )

[maxithron]

Beh....finita la tridiagonale passerò al tensore di Riemann & Christoffel
e al diavolo tutto il resto! :p

[a2000]

allora ti tolgo dal mucchio dei 5 miliardi. :D

rimaniamo io, tu e la bionda naturale ... ti piacciono i panini ? :pig:
o solo molto relativamente ?