Umberto Marzo - Tesi di Laurea

Introduzione | Cap.1 | Cap.2 | Cap.3 | Cap.4 | Cap.5 | Appendice | Bibliografia
APPENDICE AL CAPITOLO 1
Riportiamo in queste pagina il programma utilizzato (Cap.1) per il calcolo delle derivate ai nodi di una matrice rettangolare, e per la relativa analisi degli errori. E' da notare che per il caso studiato era possibile la costruzione di un algoritmo piu' efficiente, tenendo conto dell' alta simmetria di distribuzione dei punti e delle proprieta' del prodotto vettoriale. Di conseguenza i tempi di calcolo riportati non vanno intesi in senso assoluto ma rapportati fra loro.
'    Data una matrice rettangolare,modello di una superficie,viene
'    costruita la matrice delle derivate prime rispetto a x (o ad y)
'    della superficie calcolate ai nodi della griglia.
'
'    Significato delle variabili utilizzate:
'
'    grid()     La matrice che campiona la superficie.
'    zx()       La matrice che conterra' le derivate calcolate dall'algoritmo.
'    zxder()    La matrice che contiene le derivate esatte ed e' usata per
'               l'analisi degli errori.
'
'    Funzioni utilizzate:
'
'    PrVetX     Sono le function che calcolano le componenti (rispett.x,y,z)
'    PrVetY     del prodotto vettoriale dei vettori indicati nella loro
'    PrVetZ     chiamata (call).
'    VetSumX    Sono le function che calcolano le componenti del vettore
'    VetSumY    somma dei vettori indicati nella call.
'    VetSumZ
'    alfa       Sono le funzioni che calcolano i coseni direttori del
'    beta       vettore risultante.
'    gamma
'    funz       La funzione utilizzata per costruire il modello della
'               superficie (la matrice grid).Modificando la definizione e'
'               possibile campionare funzioni diverse.
'    funzder    La funzione utilizzata per la costruzione della matrice
'               delle derivate esatte.
'               ATTENZIONE se si desidera campionare una nuova funzione
'               introdurre la sua definizione nelle function funz e la
'               definizione della sua derivata nella function funzder.
'               Se non c'e' corrispondenza tra la funzione e la sua derivata
'               i risultati saranno privi di significato.
'    norma1     La function che calcola la norma della matrice passata
'               come argomento nella call.


DECLARE FUNCTION norma1! (zx!(), ndim AS INTEGER, mdim AS INTEGER)
DECLARE FUNCTION funzder! (x!, y!)
DECLARE FUNCTION alfa! (vx!, vy!, vz!)
DECLARE FUNCTION beta! (vx!, vy!, vz!)
DECLARE FUNCTION gamma! (vx!, vy!, vz!)
DECLARE FUNCTION PrVetX! (x1!, y1!, z1!, x2!, y2!, z2!)
DECLARE FUNCTION PrVetY! (x1!, y1!, z1!, x2!, y2!, z2!)
DECLARE FUNCTION PrVetZ! (x1!, y1!, z1!, x2!, y2!, z2!)
DECLARE FUNCTION VetSumX! (vx!(), n%)
DECLARE FUNCTION VetSumY! (vy!(), n%)
DECLARE FUNCTION VetSumZ! (vz!(), n%)
DECLARE FUNCTION funz! (x!, y!)


DIM grid(20, 20), zx(20, 20), zxder(20, 20)
DIM VettorX(10), VettorY(10), VettorZ(10)
DIM InnerProdX(30), InnerProdY(30), InnerProdZ(30)
DIM xdim1 AS INTEGER, ydim1 AS INTEGER, xdim AS INTEGER, ydim
AS INTEGER
DIM nc AS INTEGER

CLS
PRINT "dominio in cui campionare la funzione:"
INPUT "intervallo delle ascisse (a,b) "; xx, xx1
INPUT "intervallo delle ordinate (c,d)"; yy, yy1
IF xx >= xx1 OR yy >= yy1 THEN RUN
INPUT "num. punti sulle x"; xdim
INPUT "num.punti sulle y"; ydim
INPUT "Valore per nc (2<=nc<=8) "; nc
IF nc < 2 OR nc > 8 THEN
    PRINT "nc="; nc; " e' fuori del range permesso.Verra'assunto nc=2"
    nc = 2
END IF
hx = (xx1 - xx) / (xdim - 1)
hy = (yy1 - yy) / (ydim - 1)

'  Costruisce una reticolazione del dominio in base ai dati forniti.
'  I valori ai nodi della griglia sono forniti dalla funzione funz
'  I valori per la matrice delle derivate esatte dalla funzder
ydim1 = ydim - 1
xdim1 = xdim - 1
FOR j% = 0 TO ydim1
    y = yy + hy * j%
    FOR i% = 0 TO xdim1
       x = xx + hx * i%
       grid(j%, i%) = funz(x, y)
       zxder(i%, j%) = funzder(x, y)
    NEXT i%
NEXT j%

Tempoinizio = TIMER

FOR j0% = 0 TO ydim1
            
   FOR i0% = 0 TO xdim1


'     a seconda del valore di nc mi sposto sulla riga data opportuna

      RESTORE: a% = 0: DO UNTIL a% = nc: READ a%: LOOP
'    costruisco i vettori PoPi per i=1....nc  le componenti di
'    questi vettori sono negli array VettorX,VettorY,VettorZ

      FOR i% = 1 TO nc
         READ jj%, ii%
         IF (j0% + jj%) < 0 THEN jj% = 1
         IF (j0% + jj%) > ydim1 THEN jj% = -1
         IF (i0% + ii%) < 0 THEN ii% = 1
         IF (i0% + ii%) > xdim1 THEN ii% = -1
         VettorX(i%) = hx * ii%
         VettorY(i%) = hy * jj%
         VettorZ(i%) = grid(j0% + jj%, i0% + ii%) - grid(j0%,i0%)
      NEXT i%

'    quindi faccio tutti i possibili prodotti vettoriali e salvo le
'    loro componenti negli array InnerProdX,InnerProdY,InnerProdZ

      k% = 0
      FOR i% = 1 TO nc
         FOR j% = i% + 1 TO nc
            k% = k% + 1
            InnerProdZ(k%) = PrVetZ(VettorX(i%), VettorY(i%),VettorZ(i%),VettorX(j%),VettorY(j%), VettorZ(j%))
            IF InnerProdZ(k%) < 0 THEN
               InnerProdX(k%) = PrVetX(VettorX(j%), VettorY(j%),VettorZ(j%),VettorX(i%), VettorY(i%), VettorZ(i%))
               InnerProdY(k%) = PrVetY(VettorX(j%), VettorY(j%),VettorZ(j%),VettorX(i%), VettorY(i%), VettorZ(i%))
               InnerProdZ(k%) = PrVetZ(VettorX(j%), VettorY(j%),VettorZ(j%),VettorX(i%), VettorY(i%), VettorZ(i%))
            ELSE
               InnerProdX(k%) = PrVetX(VettorX(i%), VettorY(i%),VettorZ(i%),VettorX(j%), VettorY(j%), VettorZ(j%))
               InnerProdY(k%) = PrVetY(VettorX(i%), VettorY(i%),VettorZ(i%),VettorX(j%), VettorY(j%), VettorZ(j%))
            END IF
         NEXT j%
      NEXT i%
 
'    quindi si calcola il vettore risultante dei prodotti vettoriali

      RisulX = VetSumX(InnerProdX(), k%)
      RisulY = VetSumY(InnerProdY(), k%)
      RisulZ = VetSumZ(InnerProdZ(), k%)

'    alfa beta e gamma sono i coseni direttori della direzione
'    della risultante

      alfa1 = alfa(RisulX, RisulY, RisulZ)
      beta1 = beta(RisulX, RisulY, RisulZ)
      gamma1 = gamma(RisulX, RisulY, RisulZ)

'    le derivate parziali prime della nostra superficie sono quelle del
'    piano passante per Po e ortogonale al vettore risultante

      zx = -alfa1 / gamma1
      zx(i0%, j0%) = zx
'     zy = -beta1 / gamma1
   NEXT i0%
NEXT j0%


'    calcolo del tempo impiegato

    tempotrascorso = TIMER - Tempoinizio
    minuti = INT(tempotrascorso / 60)
    secondi = INT(tempotrascorso)
    centesimi = INT((tempotrascorso - secondi) * 100)


'   output dei risultati ottenuti


PRINT "      ";
FOR i0% = 0 TO xdim1
    PRINT USING " ###.##"; xx + i0% * hx;
NEXT i0%
PRINT
FOR j0% = 0 TO ydim1
    PRINT USING " ###.##"; yy + j0% * hy;
    FOR i0% = 0 TO xdim1
        PRINT USING " ###.##"; zx(i0%, j0%);
'       PRINT USING " ###.##"; -beta1 / gamma1
    NEXT i0%
    PRINT
NEXT j0%

PRINT
PRINT "Tempo impiegato per la costruzione della matrice:";
minuti; "'"; secondi; "''"; centesimi


'  qui parte il calcolo degli errori


ErrMax = 0: ErrMin = 100
FOR j0% = 1 TO ydim1 - 1
FOR i0% = 1 TO xdim1 - 1
'     PRINT zx(i0%, j0%); zxder(j0%, i0%)
     Errore = ABS(zx(i0%, j0%) - zxder(i0%, j0%))
     IF Errore > ErrMax THEN ErrMax = Errore
     IF Errore < ErrMin THEN ErrMin = Errore
NEXT i0%
NEXT j0%


PRINT "Errore max:"; ErrMax
PRINT "Errore min:"; ErrMin
PRINT "Norma matrice delle derivate:"; norma1(zx(),ydim1 - 1, xdim1 - 1)
PRINT "Norma della matrice delle derivate esatte:";norma1(zxder(), ydim1 - 1, xdim1 - 1)


'    questi sono gli indici dei punti che hanno proiezione piu' vicina a Po


DATA 2,-1,0,0,-1
DATA 3,-1,0,1,0,0,-1
DATA 4,-1,0,1,0,0,-1,0,1
DATA 5,-1,0,1,0,0,-1,0,1,-1,-1
DATA 6,-1,0,1,0,0,-1,0,1,-1,-1,1,1
DATA 7,-1,0,1,0,0,-1,0,1,-1,-1,1,1,1,-1
DATA 8,-1,0,1,0,0,-1,0,1,-1,-1,1,1,1,-1,-1,1


FUNCTION alfa (vx, vy, vz)
alfa = vx / (SQR(vx * vx + vy * vy + vz * vz))
END FUNCTION


FUNCTION beta (vx, vy, vz)
beta = vy / (SQR(vx * vx + vy * vy + vz * vz))
END FUNCTION


' questa e' la definizione della funzione utilizzata per la
' costruzione della griglia
'
FUNCTION funz (x, y)
'funz = 10   ' this is a plane
'funz = x + y
'funz = 3 / 2 * (COS(3 / 5 * (y - 1)) + 5 / 4) / (1 + ((x - 4)
/ 3) ^ 2)
'funz = ABS(x + y)
'funz = SIN(x ^ 2) + COS(y ^ 2) + SIN(x)
funz = SIN(x * y) + COS(x)
' ex = 2.71 ^ (y - x): ex1 = 2.71 ^ (-(y - x)): funz = (ex - ex1)
/ (ex + ex1) + 1
'funz = 8 / 3.14 * EXP(-8 * (x * x + y * y))
'funz = 9 * (3 / 4 * EXP(-(x - 3) ^ 2 - (y - 3) ^ 2) / 4) + EXP(-(x
/ 7) ^ 2 - (y / 10) - 1 / 5    * EXP(-(x - 5) ^ 2 - (y - 8) ^
2) + 1 / 2 * EXP((-(x - 8) ^ 2 - (y -     4) ^ 2) / 4))
'funz = x * x * y
END FUNCTION


FUNCTION funzder (x, y)
funzder = y * COS(x * y) - SIN(x)
END FUNCTION


FUNCTION gamma (vx, vy, vz)
gamma = vz / (SQR(vx * vx + vy * vy + vz * vz))
END FUNCTION


FUNCTION norma1 (zx(), ndim AS INTEGER, mdim AS INTEGER)
FOR j% = 1 TO ndim
FOR i% = 1 TO mdim
s = s + ABS(zx(i%, j%))
'PRINT zx(i%, j%)
NEXT i%
IF s > nor THEN nor = s
s = 0
NEXT j%
norma1 = nor
END FUNCTION


FUNCTION PrVetX (x1, y1, z1, x2, y2, z2)
PrVetX = y1 * z2 - z1 * y2
END FUNCTION


FUNCTION PrVetY (x1, y1, z1, x2, y2, z2)
PrVetY = z1 * x2 - x1 * z2
END FUNCTION


FUNCTION PrVetZ (x1, y1, z1, x2, y2, z2)
PrVetZ = x1 * y2 - y1 * x2
END FUNCTION


FUNCTION VetSumX (vx(), n%)
FOR i% = 1 TO n%
sum = sum + vx(i%)
NEXT i%
VetSumX = sum
END FUNCTION


FUNCTION VetSumY (vy(), n%)
FOR i% = 1 TO n%
sum = sum + vy(i%)
NEXT i%
VetSumY = sum
END FUNCTION


FUNCTION VetSumZ (vz(), n%)
FOR i% = 1 TO n%
sum = sum + vz(i%)
NEXT i%
VetSumZ = sum
END FUNCTION
Indietro | Home | Bibliografia
Introduzione | Cap.1 | Cap.2 | Cap.3 | Cap.4 | Cap.5 | Appendice | Bibliografia