Télécharger la plaquette Cast3M


fr en

Télécharger ytmx.eso

Retour à la liste

Numérotation des lignes : 

ytmx

  1. C YTMX      SOURCE    CB215821  20/11/25    13:44:36     10792          
  2.       SUBROUTINE YTMX(IRE1,IRE2,IRE3,VA)
  3.       IMPLICIT INTEGER(I-N)
  4.       IMPLICIT REAL*8(A-H,O-Z)
  5. C
  6. C  ****  MULTIPLICATION D'UNE MATRICE(IRE3) PAR UN CHAMPPOINT (IRE1) A
  7. C  ****  GAUCHE ET PAR UN CHAMPPOINY (IRE2) A DROITE.
  8. C  ****  VA=    IRE1 *IRE3 *IRE2
  9. C  ****  POUR EFFECTUER L'OPERATION ON ELIMINE LES COMPOSANTES LX
  10. C  ****  DU CHPOINT ET DE LA RIGIDITE. ON TESTE QUE LES AUTRES INCONNUES
  11. C  ****  DU CHPOINT SONT INCLUSES DANS CELLES DE L OBJET RIGIDITE
  12. C  ON SUPPOSE QUE :
  13. C   1. YT ET X SONT DE MEME TYPE
  14. C   2. LA MATRICE EST CARREE
  15. C   3. LA MATRICE POSSEDE DES CORRESPONDANCES SUR LES INCONNUES
  16. C      (C'EST A DIRE QUE LA IEME LIGNE EST LA DUALE DE LA IEME COLONNE)
  17. C
  18. C   BP , avril 2010 : on supprime l hypothese 3.
  19. C  (pour permettre l utilisation de matrices crees par imped par ex.)
  20. C
  21.  
  22. -INC PPARAM
  23. -INC CCOPTIO
  24. -INC CCREEL
  25. -INC SMELEME
  26. -INC SMCHPOI
  27. -INC SMRIGID
  28. -INC SMCOORD
  29. -INC CCHAMP
  30. c
  31.       CHARACTER*(LOCOMP) NI,INCOM
  32.       SEGMENT,IHAR(0)
  33.       SEGMENT,SIINC
  34.         CHARACTER*(LOCOMP) IINC(0)
  35.       ENDSEGMENT
  36.       SEGMENT ICPR(nbpts)
  37.       SEGMENT/ITRAV/(CC(NLIGMA)*D,DD(NLIGMA)*D,VAA(NNIN,ITES)*D,
  38.      *VBB(NNIN,ITES)*D)
  39.       SEGMENT IPOS(NLIGMA)
  40.       SEGMENT IPOS2(NLIGMA)
  41.       SEGMENT SIINCO
  42.         CHARACTER*(LOCOMP) IINCO(NNIN)
  43.       ENDSEGMENT
  44.  
  45. C ITES = NONBRE DE NOEUD DU CHPOINT
  46. C NLIGMA = TAILLE MAX D'UNE LIGNE DE MATRICE DE RIGIDITE ELEMENTAIRE
  47.  
  48. C **** INITIALISATION DU RESULTAT
  49.       VA=0.D0
  50. C
  51. C
  52. C  **** ON RETIRE DES CHPOINTS LES MULT. DE LAGRANGE S'IL Y EN A.
  53. C
  54. C
  55. C  **** ON CREE LES TABLEAUX :
  56. c
  57. C  **** ICPR(I)=J VEUT DIRE QUE LE NOEUD I A LE NUMERO LOCAL J.
  58. C
  59.       SEGINI ICPR
  60.       KMAX=nbpts
  61.       DO 6 K=1,KMAX
  62.       ICPR(K)=0
  63.     6 CONTINUE
  64.       IK=0
  65.       IPA=0
  66. C  ** ON COMMENCE PAR RECENSER LES NOEUDS DU 1er CHPOINT
  67.       MCHPOI=IRE1
  68.   50  CONTINUE
  69.       IPA=IPA+1
  70.       SEGACT,MCHPOI
  71.       NSOUPO=IPCHP(/1)
  72.       DO 1 ISOU=1,NSOUPO
  73.       MSOUPO=IPCHP(ISOU)
  74.       SEGACT,MSOUPO
  75.       MELEME=IGEOC
  76.       SEGACT,MELEME
  77.       N2=NUM(/2)
  78.       DO 2 I2=1,N2
  79.       K=NUM(1,I2)
  80. *     on ajoute le noeud K a ICPR(K) si pas deja vu
  81.       IF(ICPR(K).NE.0)  GO TO 2
  82.       IK=IK+1
  83.       ICPR(K)=IK
  84.    2  CONTINUE
  85.       SEGDES,MELEME
  86.       SEGDES,MSOUPO
  87.    1  CONTINUE
  88.       SEGDES,MCHPOI
  89. C  ** ON CONTINUE AVEC LES NOEUDS DU 2nd CHPOINT
  90.       MCHPOI=IRE2
  91.       IF(IPA.EQ.1)  GO TO 50
  92.  
  93.       ITES=IK
  94.       NLIGMA=0
  95.  
  96. C  **** REMPLISSAGE DE IINC et IHAR
  97. C       = couple(inconnue primale + harmonique) de la matrice MRIGID
  98.       SEGINI,SIINC
  99.       SEGINI,IHAR
  100.       MRIGID=IRE3
  101.       SEGACT,MRIGID
  102.       NRIGE=IRIGEL(/1)
  103.       NRIGEL=IRIGEL(/2)
  104.       DESCR=IRIGEL(3,1)
  105.       SEGACT,DESCR
  106. *     Initialisation de la 1ere valeur
  107.       IINC(**)=LISINC(1)
  108.       IHAR(**)=IRIGEL(5,1)
  109.       ININC=1
  110. *     boucle sur les rigidites elementaires
  111.       DO 3 IRI=1,NRIGEL
  112.         MELEME=IRIGEL(1,IRI)
  113.         SEGACT,MELEME
  114.         DESCR=IRIGEL(3,IRI)
  115.         NOHA=IRIGEL(5,IRI)
  116.         SEGACT,DESCR
  117.         NLIGRE=LISINC(/2)
  118.         IF(NLIGRE.GT.NLIGMA) NLIGMA=NLIGRE
  119.         DO 7 I2=1,NLIGRE
  120.           DO 8 I1=1,ININC
  121.             IF(LISINC(I2).NE.IINC(I1)) GO TO 8
  122.             IF(NOHA.EQ.IHAR(I1)) GO TO 7
  123.    8      CONTINUE
  124. *         on ajoute le couple INC+HAR si pas deja vu
  125.           IINC(**)=LISINC(I2)
  126.           IHAR(**)=NOHA
  127.           ININC=ININC+1
  128.    7    CONTINUE
  129.         SEGDES,DESCR
  130.   11    SEGDES MELEME
  131.    3  CONTINUE
  132.       SEGDES,MRIGID
  133. C
  134. C  ****  ON INITIALISE LE SEGMENT ITRAV
  135. C
  136.       NNIN=ININC
  137.       SEGINI SIINCO,IPOS,IPOS2
  138.       DO  I=1,NNIN
  139.         IINCO(I)=IINC(I)
  140.       enddo 
  141.       SEGINI ITRAV
  142. C      CALL ZERO(CC,NLIGMA,1)
  143. C      CALL ZERO(DD,NLIGMA,1)
  144. C      CALL ZERO(VAA,NNIN,ITES)
  145. C      CALL ZERO(VBB,NNIN,ITES)
  146. C
  147. C  ****  ON INITIALISE VAA (= 1er chpoint) et VBB (= 2nd chpoint)
  148. C
  149. C  ** LE 1er CHPOINT
  150.       MCHPOI=IRE1
  151.       IPA=0
  152.   51  CONTINUE
  153.       IPA=IPA+1
  154.       SEGACT,MCHPOI
  155.       NSOUPO=IPCHP(/1)
  156. c --- boucle sur les zones ---
  157.       DO 15 ISOU=1,NSOUPO
  158.         MSOUPO=IPCHP(ISOU)
  159.         SEGACT,MSOUPO
  160.         MELEME=IGEOC
  161.         SEGACT,MELEME
  162.         MPOVAL=IPOVAL
  163.         SEGACT,MPOVAL
  164.         N2=VPOCHA(/1)
  165.         NC=VPOCHA(/2)
  166. c    -- boucle sur les composantes --
  167.         DO 16 INC=1,NC
  168.           INCOM=NOCOMP(INC)
  169.           NOHA=NOHARM(INC)
  170.           DO 17 IH=1,NNIN
  171.             IF(INCOM.NE.IINCO(IH)) GO TO 17
  172.             IF(NOHA.EQ.IHAR(IH)) GO TO 18
  173.    17     CONTINUE
  174.           GO TO 16
  175. c         on a trouvé le bon couple inconnue primale+harmonique : IH
  176.    18     CONTINUE
  177.           IF(IPA.EQ.1)  GO TO 190
  178. c         On recupere dans VBB(IH,IK)
  179. c         la valeur du 2nd chpoint pour l inconnue IH et le point IK
  180.           DO 191 I2=1,N2
  181.           IK=ICPR(NUM(1,I2))
  182.           VBB(IH,IK)=VPOCHA(I2,INC)
  183. c          VBBB=VPOCHA(I2,INC)
  184.  1111     FORMAT(1X,I5,1X,I5,1X,I5,1X,1PE12.5)
  185.   191     CONTINUE
  186.           GO TO 16
  187.   190     CONTINUE
  188. c         On recupere dans VAA(IH,IK)
  189. c         la valeur du 1er chpoint pour l inconnue IH et le point IK
  190.           DO 19 I2=1,N2
  191.           IK=ICPR(NUM(1,I2))
  192.           VAA(IH,IK)=VPOCHA(I2,INC)
  193. c          VAAA=VPOCHA(I2,INC)
  194.    19     CONTINUE
  195.    16   CONTINUE
  196.         SEGDES,MSOUPO
  197.         SEGDES,MPOVAL
  198.         SEGDES,MELEME
  199.    15 CONTINUE
  200.       SEGDES,MCHPOI
  201. c     on recommence pour le 2nd chpoint
  202.       MCHPOI=IRE2
  203.       IF(IPA.EQ.1) GO TO 51
  204. C
  205. C  ****  BOUCLE 20 SUR LES OBJETS RIGIDITES ELEMENTAIRES
  206. C
  207.       SEGACT,MRIGID
  208.       DO 20 IRI=1,NRIGEL
  209.         IANTI=0
  210.         IF (NRIGE.GE.7) THEN
  211.            IANTI=IRIGEL(7,IRI)
  212.         ENDIF
  213.         MELEME=IRIGEL(1,IRI)
  214.         SEGACT,MELEME
  215.         DESCR=IRIGEL(3,IRI)
  216.         SEGACT,DESCR
  217.         LISI=LISINC(/2)
  218. C
  219. C    ** ON VERIFIE QUE:
  220. C      -LA MATRICE EST CARREE
  221.         LISD=LISDUA(/2)
  222.         IF ( LISI.NE.LISD) THEN
  223.           CALL ERREUR(21)
  224.           RETURN
  225.         ENDIF
  226. C      -NOELED ET NOELEP SONT LES MEMES
  227.         DO ITEFR=1,LISI
  228.           IF( NOELED(ITEFR).NE.NOELEP(ITEFR) ) THEN
  229.             CALL ERREUR(21)
  230.             RETURN
  231.           ENDIF
  232.         ENDDO
  233. C
  234. C    ** ON REMPLIT IPOS(I)=J QUI DIT QUE LA IEME INCONNUE PRIMALE
  235. C       DE LA MATRICE ELEMENTAIRE EST LA JEME DE IINC
  236.         NOHA=IRIGEL(5,IRI)
  237.         DO 21 IN=1,LISI
  238.             NI=LISINC(IN)
  239.             DO 22 IJ=1,ININC
  240.                IF(NI.NE.IINC(IJ))  GO TO 22
  241.                IF(NOHA.EQ.IHAR(IJ))  GO TO 23
  242.  22         CONTINUE
  243.  23         CONTINUE
  244.             IPOS(IN)=IJ
  245.  21     CONTINUE
  246. C
  247. C    ** ON ETABLIT LA CORRESPONDANCE INCONNUES PRIMALES ET DUALES
  248. C       (important si hypothèse 3 non vérifiée)
  249. C    ** ON REMPLIT IPOS2(I)=J QUI DIT QUE LA IEME INCONNUE DUALE
  250. C       DE LA MATRICE ELEMENTAIRE EST "NATURELLEMENT" ASSOCIEE A LA
  251. C       JEME INCONNUE PRIMALE DE IINC
  252.         if(IIMPI.ge.5) write(6,*) 'Pour la rigidite elementaire ',IRI
  253.         DO IN=1,LISI
  254.            NI=LISDUA(IN)
  255.            if(IIMPI.ge.5)
  256.      &     write(6,*) 'l inconnue primale ',LISINC(IN),
  257.      &                ' produit le dual ',NI
  258.            do idu=1,LNOMDU
  259.               if(NOMDU(idu).eq.NI) goto 25
  260.            enddo
  261.            write(6,*) 'ERREUR : NOM D INCONNUE DUALE INCONNUE',NI
  262.            CALL ERREUR(21)
  263.            return
  264.  25        continue
  265. c          on a trouve le numero du dual -> on en deduit le primal
  266. C          "naturellement" associé pour le produit scalaire
  267. c          il faut le chercher dans le chpoint VBB cad dans IINC
  268.            NI=NOMDD(idu)
  269.            DO 26 IJ=1,ININC
  270.               IF(NI.NE.IINC(IJ))  GO TO 26
  271.               IF(NOHA.EQ.IHAR(IJ))  GO TO 27
  272.  26        CONTINUE
  273.            write(6,*) 'ERREUR : NOM D INCONNUE PRIMALE INCONNUE',NI
  274.            CALL ERREUR(21)
  275.            return
  276.  27        CONTINUE
  277.            IPOS2(IN)=IJ
  278.         ENDDO
  279.  
  280.  
  281. C
  282. C  ****  BOUCLE 30 SUR LES PETITES MATRICES D'UNE RIGIDITE ELEMENTAIRE
  283.          N1=NUM(/1)
  284.          N2=NUM(/2)
  285.          xMATRI=IRIGEL(4,IRI)
  286.          SEGACT,xMATRI
  287.          COER=COERIG(IRI)
  288.          DO 30 I2=1,N2
  289. C
  290. C       **  AVANT D'EFFECTUER LE PRODUIT ON VERIFIE QU'IL EST A FAIRE
  291.             DO 31 I1=1,N1
  292.                IF(ICPR(NUM(I1,I2)).NE.0) GO TO 32
  293.  31         CONTINUE
  294.             GO TO 30
  295.  32         CONTINUE
  296. C
  297. C       **  FABRICATION D'UN VECTEUR ISSU DU 1er CHPOINT
  298.             DO 33 JN=1,LISI
  299.                CC(JN)=0.D0
  300.                J2=ICPR(NUM(NOELEP(JN),I2))
  301.                IF(J2.EQ.0) GO TO 33
  302.                J1=IPOS(JN)
  303.                CC(JN)=VAA(J1,J2)
  304.  33         CONTINUE
  305.             if(IIMPI.ge.5)
  306.      &      write(6,*) 'CC=',(CC(iou),iou=1,LISI)
  307. C       **  FABRICATION D'UN VECTEUR ISSU DU 2nd CHPOINT
  308.             DO 330 IN=1,LISI
  309.                DD(IN)=0.D0
  310.                J2=ICPR(NUM(NOELEP(IN),I2))
  311.                IF(J2.EQ.0)  GO TO 330
  312.                J1=IPOS2(IN)
  313.                DD(IN)=VBB(J1,J2)
  314.  330        CONTINUE
  315.             if(IIMPI.ge.5)
  316.      &      write(6,*) 'DD=',(DD(iou),iou=1,LISI)
  317. C
  318. C  **** BOUCLE 35 SUR LES LIGNES D'UNE MATRICE ELEMENTAIRE
  319. *            XMATRI=IMATTT(I2)
  320. *            SEGACT,XMATRI
  321.             DO 35 IN=1,LISI
  322. *               IF (ABS(DD(IN)).GT.1.D-10) THEN
  323.                IF (ABS(DD(IN)).GT.XPETIT) THEN
  324.                      VB=0.D0
  325. C                 ** BOUCLE 38 SUR LES COLONNES D'UNE MATRICE ELEMENTAIRE
  326.                      DO 38 JN=1,LISI
  327.                            VB=VB+CC(JN)*RE(IN,JN,i2)
  328.  38                  CONTINUE
  329.                      VA=VA+ DD(IN)*VB*COER
  330.                ENDIF
  331.  35         CONTINUE
  332.             if(IIMPI.ge.5)
  333.      &      write(6,*) 'VA=',VA
  334. *            SEGDES,XMATRI
  335.  30      CONTINUE
  336.          SEGDES,xMATRI
  337.          SEGDES,DESCR
  338.  24      SEGDES MELEME
  339.  20   CONTINUE
  340.       SEGDES,MRIGID
  341.       SEGSUP,ITRAV
  342.       SEGSUP,SIINC
  343.       SEGSUP,IHAR
  344.       SEGSUP ICPR,SIINCO,IPOS
  345.  5000 CONTINUE
  346. *
  347.       END
  348.  
  349.  
  350.  
  351.  
  352.  
  353.  
  354.  

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