Télécharger la plaquette Cast3M


fr en

Télécharger versy2.eso

Retour à la liste

Numérotation des lignes : 

versy2

  1. C VERSY2    SOURCE    GOUNAND   24/11/12    21:15:10     12076          
  2. *
  3. ************************************************************************
  4. * verification symetrie de la matrice
  5. *
  6. * IANTI = 0  : matrice symetrique       ->  Aij =  Aji  pour i<j
  7. * IANTI = 1  : matrice antisymetrique   ->  Aij = -Aji  (et  Aii=0)
  8. * IANTI = 2  : matrice quelconque       ->  pas de verification
  9. *
  10. ************************************************************************
  11.       SUBROUTINE versy2(re,ire1,ire2,ire3,IANTI,kerre)
  12.  
  13.       IMPLICIT REAL*8 (a-h,o-z)
  14. -INC PPARAM
  15. -INC CCOPTIO
  16. -INC CCREEL
  17.       REAL*8 re(ire1,ire2,ire3)
  18. *     facteur de blocage pour optimiser le cache
  19.       PARAMETER (NINC=32)
  20.       PARAMETER (LIMTAIL=65532)
  21.  
  22.       kerre=0
  23.       mdim=ire1
  24. *     parametres
  25.       xzref=(xpetit/xzprec)
  26.       xzpref=1d1*xzprec
  27.  
  28. *     boucle sur les termes ir,jr (limite sur ir pour optimiser cache)
  29.  
  30. *   - Cas symetrique
  31.       IF(IANTI.eq.0.and.ire1.lt.limtail) THEN
  32.         do iel=1,ire3
  33.           do i=1,mdim,ninc
  34.             do jr=i+1,mdim
  35.                 redj=abs(re(jr,jr,iel))
  36.               do ir=i,min(jr-1,i+ninc-1)
  37.                 re1=re(ir,jr,iel)
  38.                 re2=re(jr,ir,iel)
  39.                 redi=abs(re(ir,ir,iel))
  40.                 xn=abs(re1-re2)
  41.                 xd=max(redi,redj,abs(re1),abs(re2))+xzref
  42.                 if (xn/xd.gt.xzpref) then
  43. c                   MOTERR=' '
  44. c                   MOTERR(1:10)='SYMETRIQUE'
  45.                   MOTERR='SYMETRIQUE'
  46.                   reaerr(1)=re1
  47.                   reaerr(2)=re2
  48.                   reaerr(3)=xn
  49.                   kerre=1044
  50.                   return
  51.                 endif
  52.                 if (xn/xd.gt.xzref) then
  53.                    rem=(re1+re2)/2.d0
  54.                    re(ir,jr,iel)=rem
  55.                    re(jr,ir,iel)=rem
  56.                 endif
  57.               enddo
  58.             enddo
  59.           enddo
  60.         enddo
  61.  
  62. *   - Cas anti-symetrique (on change les signes)
  63.       ELSEIF(IANTI.eq.1) THEN
  64.         do iel=1,ire3
  65.           do i=1,mdim,ninc
  66.             do jr=1,mdim
  67.               do ir=i,min(jr-1,i+ninc-1)
  68.                 re1=re(ir,jr,iel)
  69.                 re2=re(jr,ir,iel)
  70.                 xn=abs(re1+re2)
  71.                 xd=max(abs(re1),abs(re2))+xzref
  72.                 if (xn.gt.xzpref*xd) then
  73. c                   MOTERR=' '
  74. c                   MOTERR(1:15)='ANTI-SYMETRIQUE'
  75.                   MOTERR='ANTI-SYMETRIQUE'
  76.                   reaerr(1)=re1
  77.                   reaerr(2)=re2
  78.                   reaerr(3)=xn
  79.                   kerre=1044
  80.                   return
  81.                 endif
  82.               enddo
  83.             enddo
  84.           enddo
  85. c         verif diagonale nulle
  86.           do ir=1,mdim
  87.                 re1=re(ir,ir,iel)
  88.                 if (abs(re1).gt.xzref) then
  89. c                   MOTERR=' '
  90. c                   MOTERR(1:15)='ANTI-SYMETRIQUE'
  91.                   MOTERR='ANTI-SYMETRIQUE'
  92.                   reaerr(1)=re1
  93.                   reaerr(2)=0.D0
  94.                   reaerr(3)=abs(re1)
  95.                   kerre=1044
  96.                   return
  97.                 endif
  98.           enddo
  99.        enddo
  100. * Matrice non verifiee
  101.       ELSE
  102.          kerre=-1
  103.       ENDIF
  104.  
  105.       END
  106.  
  107.  

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