Télécharger la plaquette Cast3M


fr en

Télécharger pdire3.eso

Retour à la liste

Numérotation des lignes : 

pdire3

  1. C PDIRE3    SOURCE    CHAT      05/01/13    02:11:26     5004
  2. c---------------------------------------------------------------------
  3. c
  4.       SUBROUTINE PDIRE3 (S, XL, V)
  5. c
  6. c=====================================================================
  7. c                                                                    =
  8. c   This routine computes the eigenvectors.                          =
  9. c                                                                    =
  10. c    Input :     xl  (3)   eigenvalues                               =
  11. c                s   (6)   original matrix                           =
  12. c    Output:     v (3,3)   eigenvectors (normalised)                 =
  13. c                                                                    =
  14. c    Note:   s  = (Sxx, Syy, Szz, Sxy, Sxz, Syz)                     =
  15. c            xl = (S11, S22, S33)                                    =
  16. c                                                                    =
  17. c=====================================================================
  18.       IMPLICIT INTEGER(I-N)
  19.       real*8    s (6), xl (3), v (3,3)
  20. c
  21.       real*8    toler, comp, compi, comp2, a1, a2, a3, dot, xn
  22. c
  23.       parameter (toler = 1.0 d-07)
  24. c
  25.       comp2 = xl (1)**2 + xl (2)**2 + xl (3)**2
  26.       compi = 1.0 d0 / sqrt(comp2)
  27. c
  28.       a1 = abs ((xl(2)-xl(3))*compi)
  29.       a2 = abs ((xl(1)-xl(3))*compi)
  30.       a3 = abs ((xl(1)-xl(2))*compi)
  31. c=====----------------------------------------------------------------
  32. c  case-a: three equal eigenvalues                                   |
  33. c=====----------------------------------------------------------------
  34.       if (a1 .lt. toler .and. a2 .lt. toler)  then
  35.           v (1,1) = 1.0 d0
  36.           v (2,2) = 1.0 d0
  37.           v (3,3) = 1.0 d0
  38.           v (1,2) = 0.0 d0
  39.           v (1,3) = 0.0 d0
  40.           v (2,1) = 0.0 d0
  41.           v (2,3) = 0.0 d0
  42.           v (3,1) = 0.0 d0
  43.           v (3,2) = 0.0 d0
  44.           return
  45.       end if
  46. c=====----------------------------------------------------------------
  47. c  case-b: two equal eigenvalues                                     |
  48. c=====----------------------------------------------------------------
  49.       if (a1 .lt. toler)  then
  50.           call pvecto (s, xl(1), v(1,1), comp2)
  51.           call vecbas (v(1,1), v(1,2), v(1,3))
  52.           return
  53.       end if
  54.       if (a2 .lt. toler)  then
  55.           call pvecto (s, xl(2), v(1,2), comp2)
  56.           call vecbas (v(1,2), v(1,3), v(1,1))
  57.           return
  58.       end if
  59.       if (a3 .lt. toler)  then
  60.           call pvecto (s, xl (3), v(1,3), comp2)
  61.           call vecbas (v(1,3), v(1,1), v(1,2))
  62.           return
  63.       end if
  64. c=====----------------------------------------------------------------
  65. c  case-c: three different eigenvalues                               |
  66. c=====----------------------------------------------------------------
  67.       call pvecto (s, xl(1), v(1,1), comp2)
  68.       call pvecto (s, xl(2), v(1,2), comp2)
  69. c
  70.       dot = v (1,1) * v (1,2) + v (2,1) * v (2,2) + v (3,1) * v (3,2)
  71. c
  72.       v (1,2) = v (1,2) - dot * v (1,1)
  73.       v (2,2) = v (2,2) - dot * v (2,1)
  74.       v (3,2) = v (3,2) - dot * v (3,1)
  75. c
  76.       xn = 1.0 d0 / SQRT (v(1,2)**2+v(2,2)**2+v(3,2)**2)
  77. c
  78.       v (1,2) = v (1,2) * xn
  79.       v (2,2) = v (2,2) * xn
  80.       v (3,2) = v (3,2) * xn
  81.       v (1,3) = v (2,1) * v (3,2) - v (3,1) * v (2,2)
  82.       v (2,3) = v (3,1) * v (1,2) - v (1,1) * v (3,2)
  83.       v (3,3) = v (1,1) * v (2,2) - v (2,1) * v (1,2)
  84. c
  85.       return
  86.       end
  87.  
  88.  
  89.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales