Extraits
[...] December 2009 input . - array of coefficients for matrix A n - dimension output . - inverse matrix of A comments . [...]
[...] Notamment on n'a pas de coefficients pour les degrés impairs dans le développement de Taylor du cosinus. 3.3 Voir SUBROUTINE sol_gen 3.4 Voir FUNCTION erreurQuadLinGen Le graphe suivant montre l'évolution de l'angle pour la solution générale (en bleu) et pour la solution harmonique (en rouge) : On constate que la période du mouvement est un peu plus élevée pour la solution non linéaire. L'erreur quadratique entre les deux solutions et la suivante : 10.5123 Annexe : code source types.f90 Name : types.f90 Author : Version : 1.0 Copyright : Description : Module contenant les types derives MODULE types IMPLICIT NONE Definition du type derive "pendule" TYPE pendule REAL l longueur du pendule REAL g = 9.81 constante de pesanteur REAL theta0 angle initial END TYPE pendule END MODULE types procedures.f90 Name : procedures.f90 Author : Version : 1.0 Copyright : Description : Module contenant les procedures (fonctions et subroutines) MODULE procedures USE types IMPLICIT NONE DOUBLE PRECISION, PARAMETER PI = 3.141592653589793 CONTAINS FONCTIONS Fonction calculant l'integrale par la formule des trapezes composites DOUBLE PRECISION FUNCTION integraleTrapeze(t, ft, idxMin, idxMax) DOUBLE PRECISION, DIMENSION(:), INTENT(IN) ft INTEGER, INTENT(IN) idxMin, idxMax INTEGER i variables locales DOUBLE PRECISION dt pas integraleTrapeze = 0.0 dt = - calcul de l'integrale par la methode des trapezes integraleTrapeze = dt/2.0 * (ft(idxMin) + ft(idxMax)) DO i = idxMin+1, idxMax-1 integraleTrapeze = integraleTrapeze + dt * ft(i) END DO END FUNCTION integraleTrapeze Fonction calculant l'erreur quadratique et l'affichant DOUBLE PRECISION FUNCTION erreurQuadratique(t, theta, polyCoef) DOUBLE PRECISION, DIMENSION(:), INTENT(IN) theta DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN) polyCoef INTEGER N dernier chiffre du numero d'etudiant INTEGER Nt variables locales DOUBLE PRECISION erreurQuadratiqueCarree = 0.0, temp = 0.0 Nt = SIZE(t) calcul de l'erreur quadratique au carre DO i = Nt temp = 0.0 DO j = N+6 temp = temp + * polyCoef(j,1) END DO erreurQuadratiqueCarree = erreurQuadratiqueCarree + (temp - theta(i))**2 END DO calcul de l'erreur quadratique erreurQuadratique = SQRT(erreurQuadratiqueCarree) affichage de l'erreur quadratique PRINT erreurQuadratique END FUNCTION erreurQuadratique Fonction permettant de calculer numeriquement la periode dans le cas non-lineaire DOUBLE PRECISION FUNCTION periode(pendule1, dtheta) TYPE (pendule), INTENT(IN) pendule1 un pendule DOUBLE PRECISION, INTENT(IN) dtheta pas INTEGER i variable locale INTEGER Nt nombre de pas DOUBLE PRECISION integrale DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE theta, ftheta INTEGER ok1 = ok2 = 1 variables utilisees pour l'allocation des deux tableaux Nt = FLOOR(pendule1%theta0 / dtheta) allocation des tableaux ALLOCATE(theta(1:Nt), STAT = ok1); IF (ok1 0 ) STOP ALLOCATE(ftheta(1:Nt), STAT = ok2); IF (ok2 0 ) STOP DO i = Nt theta(i) = dtheta * ftheta(i) = 1.0/SQRT(COS(theta(i)) - COS(pendule1%theta0)) END DO integrale = integraleTrapeze(theta, ftheta Nt) periode = 4 * SQRT(pendule1%l/(2.0*pendule1%g)) * integrale END FUNCTION periode Fonction calculant l'erreur quadratique entre solution lineaire et generale et l'affichant DOUBLE PRECISION FUNCTION erreurQuadLinGen(theta_lin, theta) DOUBLE PRECISION, DIMENSION(:), INTENT(IN) theta_lin, theta INTEGER Nt variables locales DOUBLE PRECISION erreurQuadratiqueCarree = 0.0 Nt = SIZE(theta_lin) calcul de l'erreur quadratique au carre DO i = Nt erreurQuadratiqueCarree = erreurQuadratiqueCarree + (theta_lin(i) - theta(i))**2 END DO calcul de l'erreur quadratique erreurQuadLinGen = SQRT(erreurQuadratiqueCarree) affichage de l'erreur quadratique PRINT F8.4)", "Erreur quadratique entre solution lineaire et generale erreurQuadLinGen END FUNCTION erreurQuadLinGen SUBROUTINES Sous-programme permettant de calculer sin(nt) SUBROUTINE sin_nt(ti, tf, nbval, theta) DOUBLE PRECISION, INTENT(IN) ti, tf extremites de l'intervalle INTEGER, INTENT(IN) nbval nombre de valeurs DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE, INTENT(OUT) theta INTEGER N dernier chiffre du numero d'etudiant INTEGER dernier chiffre du numero d'etudiant DOUBLE PRECISION dt pas INTEGER i variable locale INTEGER ok1 = ok2 = 1 variables utilisees pour l'allocation des deux tableaux allocation des tableaux ALLOCATE(t(1:nbval), STAT = ok1); IF (ok1 0 ) STOP ALLOCATE(theta(1:nbval), STAT = ok2); IF (ok2 0 ) STOP calcul du pas dt = (tf - ti)/(nbval - calcul des tableaux t et theta DO i = nbval = ti + theta(i) = END DO END SUBROUTINE sin_nt Calcule la solution de la theorie linearisee SUBROUTINE sol_lin (pendule1, ti, tf, theta_lin) IMPLICIT NONE TYPE (pendule), INTENT(IN) pendule1 DOUBLE PRECISION, INTENT(IN) ti, tf extremites de l'intervalle INTEGER, INTENT(IN) N longueur du vecteur t DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE, INTENT(OUT) theta_lin INTEGER i variable locale DOUBLE PRECISION dt pas DOUBLE PRECISION w0 pulsation INTEGER ok1 = ok2 = 1 variables utilisees pour l'allocation de tableau dt = (tf - - w0 = SQRT(pendule1%g / pendule1%l) allocation des tableaux ALLOCATE(t(1:N), STAT = ok1); IF (ok1 0 ) STOP ALLOCATE(theta_lin(1:N), STAT = ok2); IF (ok2 0 ) STOP DO i = N = theta_lin(i) = pendule1%theta0 * COS(w0*t(i)) END DO END SUBROUTINE sol_lin Calcule la solution de la theorie generale (non linearisee) SUBROUTINE sol_gen (pendule1, theta) IMPLICIT NONE TYPE (pendule), INTENT(IN) pendule1 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE, INTENT(OUT) theta INTEGER, INTENT(IN) N longueur du vecteur t INTEGER i variable locale DOUBLE PRECISION dt pas DOUBLE PRECISION w0 pulsation DOUBLE PRECISION per periode INTEGER ok1 = ok2 = 1 variables utilisees pour l'allocation de tableau per = periode(pendule1, 1.0e-5*PI/180.0) dt = per/(N - w0 = SQRT(pendule1%g / pendule1%l) allocation des tableaux ALLOCATE(t(1:N), STAT = ok1); IF (ok1 0 ) STOP ALLOCATE(theta(1:N), STAT = ok2); IF (ok2 0 ) STOP theta(1) = pendule1%theta0 theta(2) = pendule1%theta0 DO i = N = END DO DO i = N-1 theta(i+1) = 2*theta(i) - theta(i-1) - dt*dt*w0**2*SIN(theta(i)) END DO END SUBROUTINE sol_gen Calcule la derivee premiere en utilisant les differences finies progressives SUBROUTINE derivee ft, dft) DOUBLE PRECISION, DIMENSION(:), INTENT(IN) ft DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE, INTENT(OUT) dft INTEGER N longueur du vecteur t INTEGER i variable locale DOUBLE PRECISION dt pas INTEGER ok = 1 variables utilisees pour l'allocation de tableau N = SIZE(t) dt = - allocation du tableau derivee ALLOCATE(dft(1:N-1), STAT = IF (ok 0 ) STOP DO i = N-1 dft(i) = - / dt; END DO END SUBROUTINE derivee Calcule la derivee seconde en utilisant les differences finies centrees SUBROUTINE deriveeSeconde ft, d2ft) DOUBLE PRECISION, DIMENSION(:), INTENT(IN) ft DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE, INTENT(OUT) d2ft INTEGER N longueur du vecteur t INTEGER i variable locale DOUBLE PRECISION dt pas INTEGER ok = 1 variables utilisees pour l'allocation de tableau N = SIZE(t) dt = - allocation du tableau derivee seconde ALLOCATE(d2ft(1:N-2), STAT = IF (ok 0 ) STOP calcul derivee DO i = N-1 d2ft(i-1) = - 2*ft(i) + / END DO END SUBROUTINE deriveeSeconde Determine les coefficients du polynome au sens des moindres carres SUBROUTINE polynomeMoindreCarres theta, polyCoef) IMPLICIT NONE DOUBLE PRECISION, DIMENSION(:), INTENT(IN) theta DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE, INTENT(OUT) polyCoef INTEGER N dernier chiffre du numero d'etudiant INTEGER Nt variables locales INTEGER ok = 1 variables utilisees pour l'allocation de tableau DOUBLE PRECISION, DIMENSION(SIZE(t),N+6) X DOUBLE PRECISION, DIMENSION(N+6,SIZE(t)) Xt DOUBLE PRECISION, DIMENSION(SIZE(t),1) Y DOUBLE PRECISION, DIMENSION(N+6,N+6) XtX DOUBLE PRECISION, DIMENSION(N+6,N+6) invXtX Nt = SIZE(t) DO i = Nt DO j = N+6 = END DO = theta(i) END DO Xt = TRANSPOSE(X) XtX = MATMUL(Xt,X) CALL inverseMatrice(XtX,invXtX,N+6) allocation du tableau de coeff ALLOCATE(polyCoef(1:N+6,1), STAT = IF (ok 0 ) STOP calcul derivee seconde polyCoef = MATMUL(invXtX,MATMUL(Xt,Y)) END SUBROUTINE polynomeMoindreCarres SUBROUTINE inverseMatrice(a,c,n) Inverse matrix Method: Based on Doolittle LU factorization for Ax=b Alex G. [...]
[...] I.4° Voir SUBROUTINE polynomeMoindreCarres et FUNCTION erreurQuadratique Le polynôme de degré 7 au sens des moindres carrés pour approcher le tableau est : L'erreur quadratique obtenue est d'environ . PARTIE II. APPLICATIONS AU PENDULE II.1° Pour simplifier on peut choisir un pendule de masse et de longueur de telle sorte que . [...]
[...] Différents pas ont été testés, de deg à deg, les périodes suivantes ont été obtenues Pas : 1.000E-01 deg ; Periode : 5.496s Pas : 1.000E-02 deg ; Periode : 6.043s Pas : 1.000E-03 deg ; Periode : 6.215s Pas : 1.000E-04 deg ; Periode : 6.277s Pas : 1.000E-05 deg ; Periode : 6.288s La période attendue était de s donc on voit bien que plus le pas est petit, plus la période calculée se rapproche de l'attendue. On note aussi que la période pour la théorie non-linéaire est légèrement plus élevée que pour la théorie linéaire. [...]
[...] Méthodes numériques pour la mécanique - Fortran Note : dans tout le projet les graphiques ont été généré à l'aide du logiciel Octave (version opensource de Matlab). Le projet a cependant été réalisé en langage Fortran 90. PARTIE I. [...]
