Télécharger la plaquette Cast3M


fr en

Télécharger rela2.eso

Retour à la liste

Numérotation des lignes : 

rela2

  1. C RELA2     SOURCE    PV        22/04/07    21:15:01     11331          
  2.       SUBROUTINE RELA2(LEMOT)
  3. ************************************************************************
  4. *   cet operateur impose des mouvements de corps rigide
  5. *   ou l accrochage de deux maillages
  6. *
  7. *   syntaxe :
  8. *          rig1 = rela   |   cori   depl    (rota)      |
  9. *                        |                              |  geo1 (geo2)
  10. *                        | ense ux uy uz ut rx ry rz rt |
  11. *
  12. *          cori  = mot cle
  13. *          ense  = mot cle
  14. *          depl  = pour imposer que la distance entre les points
  15. *                  reste constante
  16. *          rota  = pour imposer que tous les points ont les memes
  17. *                  rotations
  18. *          ux    = pour imposer que la valeur sur la composante ux
  19. *                  est la meme pour tous les points
  20. *          geo1  = objet de type meleme ou point
  21. *          geo2  = objet de type meleme
  22. *          rig1  = objet de type rigidite
  23. *
  24. *   accrochage syntaxe :
  25. *          rig1 = RELA  geo1 ACCR  geo2  ; ( geo2 massif )
  26. *
  27. ************************************************************************
  28.  
  29.       IMPLICIT INTEGER(I-N)
  30.       IMPLICIT REAL*8 (A-H,O-Z)
  31.  
  32. -INC PPARAM
  33. -INC CCOPTIO
  34. -INC CCGEOME
  35. -INC SMELEME
  36. -INC SMCOORD
  37. -INC SMRIGID
  38. -INC CCHAMP
  39. c
  40.       SEGMENT MOTDDV
  41.       CHARACTER*4 MOTDDL(0)
  42.       ENDSEGMENT
  43.       SEGMENT NOMINV
  44.       CHARACTER*4 NOMINC(0)
  45.       ENDSEGMENT
  46.       SEGMENT NUMAIT(NBMAIT)
  47.       SEGMENT NUESCL(NBESCL)
  48. c
  49.       CHARACTER*4 ITCORI(3),LEMOT
  50.       CHARACTER*4 MOACCR(2)
  51.       CHARACTER*4 MODEPL(6),MODEDU(6)
  52.       CHARACTER*8 CTYP
  53.       DIMENSION COEF(3)
  54. c
  55.       DATA ITCORI(1)/'DEPL'/
  56.       DATA ITCORI(2)/'ROTA'/
  57.       DATA ITCORI(3)/'NOVE'/
  58.       DATA MODEPL(1)/'UX  '/
  59.       DATA MODEPL(2)/'UY  '/
  60.       DATA MODEPL(3)/'UZ  '/
  61.       DATA MODEPL(4)/'UR  '/
  62.       DATA MODEPL(5)/'UZ  '/
  63.       DATA MODEPL(6)/'UT  '/
  64.       DATA MODEDU(1)/'FX  '/
  65.       DATA MODEDU(2)/'FY  '/
  66.       DATA MODEDU(3)/'FZ  '/
  67.       DATA MODEDU(4)/'FR  '/
  68.       DATA MODEDU(5)/'FZ  '/
  69.       DATA MODEDU(6)/'FT  '/
  70.       DATA MOACCR(1)/'FORT'/
  71.       DATA MOACCR(2)/'FAIB'/
  72. c
  73.       DATA ZERO,UNO   /0.D0,1.D0/
  74.  
  75.  
  76. * on teste la compatibilité des données ******************************
  77. * en fourier le mouvement de corps rigide n'a pas toujours de sens
  78. cbp       IF ( LEMOT.NE.'ENSE') THEN
  79.       IF ( LEMOT.EQ.'CORI') THEN
  80. * cette remarque ne vaut que pour l'option cori
  81.          IF (((IDIM.EQ.3).AND.(IFOMOD.NE.2)).OR.
  82.      +        ((IDIM.EQ.2).AND. (IFOMOD.EQ.2)).OR.
  83.      +        ((IFOUR.EQ.1).AND.( ABS(NIFOUR).GT.1))) THEN
  84.             CALL ERREUR(21)
  85.             RETURN
  86.          ENDIF
  87.       ENDIF
  88. c
  89. c on cherche a lire les directions ***********************************
  90.       SEGINI MOTDDV,NOMINV
  91.       IDEPL=0
  92.       IROTA=0
  93.       IENSE=0
  94.       JVER=1
  95.  
  96. c
  97.       IF(IFOUR.NE.-2.AND.IFOUR.NE.-1.AND.IFOUR.NE.-3)GOTO 101
  98. c
  99. c     deformations planes ou contraintes planes ou def. planes generalisées
  100.       LDEPL=2
  101.       IADEPL=0
  102.       IAROTA=2
  103.       GOTO 107
  104.  
  105.   101 CONTINUE
  106.       IF(IFOUR.NE.1) GOTO 103
  107. c
  108. c     fourier
  109.       LDEPL=3
  110.       IADEPL=3
  111.       IAROTA=2
  112.       GOTO 107
  113.  
  114.   103 CONTINUE
  115.       IF(IFOUR.NE.2) GOTO 107
  116. c
  117. c     tridim
  118.       LDEPL=3
  119.       IADEPL=0
  120.       IAROTA=0
  121.  
  122.   107 CONTINUE
  123.  
  124. C option ACCRO ou GLIS *********************************************
  125.       IF(LEMOT.EQ.'ACCR'.OR.LEMOT.EQ.'GLIS') THEN
  126.       igli=0
  127.       iaccr=0
  128.       if (LEMOT.EQ.'GLIS' ) then
  129.          igli=1
  130.       else
  131. *       ACCRO : lecture facultative du mot clé FORT ou FAIBLE
  132.          CALL LIRMOT(MOACCR,2,iaccr,0)
  133.       endif
  134.       if (iaccr.eq.0) then
  135. c       RELA GLIS ou RELA ACCRO => ACCRO
  136.          CALL ACCRO (igli)
  137.       elseif (iaccr.eq.1) then
  138. c       RELA ACCRO 'FORT' => ACCRO3
  139.          CALL ACCRO3
  140.       elseif (iaccr.eq.2) then
  141. c       RELA ACCRO 'FAIBL' => ACCRO2
  142.          CALL ACCRO2
  143.       else
  144. c       Option indisponible
  145.          call erreur(19)
  146.       endif
  147.       RETURN
  148.       ENDIF
  149.  
  150. c option ENSE ******************************************************
  151.       IF(LEMOT.EQ.'ENSE') THEN
  152.  
  153.       IENSE = 1
  154.   445 CONTINUE
  155.       CALL LIRMOT(NOMDD,LNOMDD,IMOT,0)
  156.       IF(IMOT.EQ.0) GO TO 446
  157.       MOTDDL(**)=NOMDD(IMOT)
  158.       MOTDDL(**)=NOMDU(IMOT)
  159.       GO TO 445
  160.   446 CONTINUE
  161.       IBDDL=MOTDDL(/2)
  162.  
  163. c option CORI ******************************************************
  164.       ELSE IF(LEMOT.EQ.'CORI') THEN
  165.  
  166. cccccccc   a t on une directive depl
  167. c
  168.       IOBL=1
  169.  4480 CONTINUE
  170.       CALL LIRMOT(ITCORI,3,IMOT,IOBL)
  171.       IF(IERR.NE.0) RETURN
  172.       IF(IMOT.EQ.0) GOTO 448
  173.       IOBL=0
  174.       GO TO (44801,44802,44803),IMOT
  175. 44801 continue
  176. c
  177. cccccccccccc on a trouve le mot deplacement
  178. c
  179.       IDEPL=1
  180.       DO 4481 IA=1,LDEPL
  181.          NOMINC(**)=MODEPL(IADEPL+IA)
  182.          NOMINC(**)=MODEDU(IADEPL+IA)
  183.  4481 CONTINUE
  184.       GO TO 4480
  185. 44802 continue
  186. c
  187. cccccccc      on a trouve le mot rotation
  188. c
  189.       IROTA=1
  190.       GO TO 4480
  191. 44803 continue
  192.       JVER=0
  193.       GO TO 4480
  194.   448 CONTINUE
  195.  
  196.       ENDIF
  197. c endif des options ENSE et CORI ***********************************
  198.  
  199. c
  200. c on cherche a lire le(s) maillage(s) ******************************
  201. c
  202. * cas de mot clé rota ou fourier 0 ou axis 0
  203.       IF ((IDEPL.EQ.1).AND.((IROTA.EQ.1).OR.
  204.      +     (NIFOUR.EQ.0 .AND. IFOMOD.EQ.1).OR.(IFOMOD.EQ.0))) THEN
  205.          CALL QUETYP(CTYP,1,IRETOU)
  206.          IF (IRETOU.EQ.0) GOTO 555
  207.          IF(CTYP.EQ.'MAILLAGE') THEN
  208.             CALL LIROBJ('MAILLAGE',KOBJET,1,IRETOU)
  209.             IF (IRETOU.EQ.0) GOTO 555
  210.             LP=0
  211. * on regarde si on a donné juste un point
  212.          ELSE IF (CTYP.EQ.'POINT   ') THEN
  213.             CALL LIROBJ('POINT   ',KOBJET,1,IRETOU)
  214.             IF (IRETOU.EQ.0) GOTO 555
  215.             LP=1
  216.          ELSE
  217.             GOTO 555
  218.          ENDIF
  219.          CALL LIROBJ('MAILLAGE',KOBJE2,0,IRETOU)
  220.          IF((IRETOU.EQ.0).AND.(LP.EQ.1)) THEN
  221.             CALL ERREUR(26)
  222.             GOTO 555
  223.          ENDIF
  224.          CALL RELA3 (KOBJET,KOBJE2,LP,IROTA)
  225.          GOTO 555
  226.       ELSE
  227. * cas sans mot clé rota (sauf  fourier 0 et axis 0)
  228.          CALL LIROBJ('MAILLAGE',KOBJET,1,IRETOU)
  229.          IF(IRETOU.EQ .0) GOTO 555
  230.       ENDIF
  231.  
  232.       MELEME=KOBJET
  233.       SEGACT MELEME
  234.       IF(ITYPEL .NE.1) CALL CHANGE(MELEME,1)
  235.       NBPOIN=NUM(/2)
  236. c
  237. c  dans le cas depl on ecrit que la distance entre les points reste
  238. c  constante ,1 matrice avec
  239. c  en 2d , 3 noeuds maitres et 2(n-3)+3 relations
  240. c  en 3d , 4 noeuds maitres et 3(n-4)+6 relations
  241. c   on prend comme noeuds maitres possibles
  242. c       les 3 premiers noeuds non alignes en 2d
  243. c       les 4 premiers noeuds non coplanaires en 3d
  244. c
  245.       IF (IDEPL.EQ.1) THEN
  246.          IF (IDIM.EQ.2) NBMAIT=3
  247.          IF (IDIM.EQ.3) NBMAIT=4
  248. c   erreur nb de pts du maillage inferieur au nb de noeuds maitres
  249.          IF (NBPOIN.LT.NBMAIT.AND.JVER.EQ.1) THEN
  250.             CALL ERREUR(364)
  251.             GOTO 998
  252.          ENDIF
  253.          IF(JVER.EQ.1) THEN
  254.             NBESCL=NBPOIN-NBMAIT
  255.             SEGINI NUMAIT,NUESCL
  256. cc   cas 2d determination des 3 noeuds maitres
  257.             IF (IDIM.EQ.2) THEN
  258.                IP1=NUM(1,1)
  259.                IP2=NUM(1,2)
  260.                IREFP1=(IP1-1)*(IDIM+1)
  261.                IREFP2=(IP2-1)*(IDIM+1)
  262.                X12=XCOOR(IREFP2+1)-XCOOR(IREFP1+1)
  263.                Y12=XCOOR(IREFP2+2)-XCOOR(IREFP1+2)
  264.                P12N=SQRT(X12**2+Y12**2)
  265.                DO 300 I=3,NBPOIN
  266.                   IP3=NUM(1,I)
  267.                   IREFP3=(IP3-1)*(IDIM+1)
  268.                   X13=XCOOR(IREFP3+1)-XCOOR(IREFP1+1)
  269.                   Y13=XCOOR(IREFP3+2)-XCOOR(IREFP1+2)
  270.                   P13N=SQRT(X13**2+Y13**2)
  271.                   TEST=X12*X13+Y12*Y13
  272.                   IF(ABS(TEST).GT.0.9997D0*(P12N*P13N)) GOTO 300
  273.                   GOTO 399
  274.  300           continue
  275.  310           continue
  276. cccc   erreur  tous les points sont alignes ou quasiment
  277.                CALL ERREUR(362)
  278.                GOTO 999
  279. cc   cas 3d determination des 4 noeuds maitres
  280.             ELSE IF (IDIM.EQ.3) THEN
  281.                IP1=NUM(1,1)
  282.                IP2=NUM(1,2)
  283.                IREFP1=(IP1-1)*(IDIM+1)
  284.                IREFP2=(IP2-1)*(IDIM+1)
  285.                X12=XCOOR(IREFP2+1)-XCOOR(IREFP1+1)
  286.                Y12=XCOOR(IREFP2+2)-XCOOR(IREFP1+2)
  287.                Z12=XCOOR(IREFP2+3)-XCOOR(IREFP1+3)
  288.                P12N=SQRT(X12**2+Y12**2+Z12**2)
  289.                DO 350 I=3,NBPOIN
  290.                   IP3=NUM(1,I)
  291.                   IREFP3=(IP3-1)*(IDIM+1)
  292.                   X13=XCOOR(IREFP3+1)-XCOOR(IREFP1+1)
  293.                   Y13=XCOOR(IREFP3+2)-XCOOR(IREFP1+2)
  294.                   Z13=XCOOR(IREFP3+3)-XCOOR(IREFP1+3)
  295.                   P13N=SQRT(X13**2+Y13**2+Z13**2)
  296.                   TEST=X12*X13+Y12*Y13+Z12*Z13
  297.                   IF(ABS(TEST).GT.0.9997D0*(P12N*P13N)) GOTO 350
  298.                   GOTO 360
  299.  350           continue
  300.                GOTO 380
  301.  360           continue
  302.                XN=Y12*Z13-Z12*Y13
  303.                YN=Z12*X13-X12*Z13
  304.                ZN=X12*Y13-Y12*X13
  305.                PN=XN**2+YN**2+ZN**2
  306.                DO 370 I=3,NBPOIN
  307.                   IP4=NUM(1,I)
  308.                   IREFP4=(IP4-1)*(IDIM+1)
  309.                   X14=XCOOR(IREFP4+1)-XCOOR(IREFP1+1)
  310.                   Y14=XCOOR(IREFP4+2)-XCOOR(IREFP1+2)
  311.                   Z14=XCOOR(IREFP4+3)-XCOOR(IREFP1+3)
  312.                   P14N=SQRT(X14**2+Y14**2+Z14**2)
  313.                   TEST=XN*X14+YN*Y14+ZN*Z14
  314.                   IF(ABS(TEST).LT.0.025D0*(PN*P14N)) GOTO 370
  315.                   GOTO 399
  316.  370           continue
  317.  380           continue
  318. cccc   erreur  tous les points sont coplanaires ou quasiment
  319.                CALL ERREUR(363)
  320.                GOTO 999
  321.             ENDIF
  322.  399        continue
  323.             NUMAIT(1)=IP1
  324.             NUMAIT(2)=IP2
  325.             NUMAIT(3)=IP3
  326.             IF (IDIM.EQ.3) NUMAIT(4)=IP4
  327.             I1=0
  328.             DO 397 I=3,NBPOIN
  329.                DO 398 J=1,NBMAIT
  330.                   IF(NUM(1,I).EQ.NUMAIT(J)) GOTO 397
  331.  398           continue
  332.                I1=I1+1
  333.                NUESCL(I1)=NUM(1,I)
  334.  397        continue
  335.          ELSEIF(JVER.EQ.0) THEN
  336.             NBESCL=NBPOIN-NBMAIT
  337.             IF(NBESCL.GT.0) THEN
  338.                SEGINI NUMAIT,NUESCL
  339.             ELSE
  340.                NBMAIT=NBPOIN
  341.                NBESCL=0
  342.                SEGINI NUMAIT
  343.             ENDIF
  344.             DO 396 I=1,NBMAIT
  345.                NUMAIT(I)=NUM(1,I)
  346.  396        CONTINUE
  347.             IF(NBESCL.GT.0) THEN
  348.                DO 395 I=1,NBESCL
  349.                   NUESCL(I)=NUM(1,I+NBMAIT)
  350.  395           CONTINUE
  351.             ENDIF
  352.          ENDIF
  353.       ENDIF
  354. c
  355. c  dans le cas ense ,on ecrit que la valeur sur un ddl
  356. c  donne est la meme pour tous les points ,ibddl/2 matrices avec
  357. c  (n-1) relations
  358. c
  359.       IF(IENSE.EQ.1) NNMAT=IBDDL/2
  360.       IF(IDEPL.EQ.1) THEN
  361.          NNMAT=1
  362. * en fourier 1 on a en plus que ur=-ut a tout point et une relation d'analogie
  363. * sur uz entre deux noeuds (p.ex. 2 noeuds maitres)
  364.          IF(ABS(NIFOUR).EQ.1) NNMAT=3
  365.       ENDIF
  366.       IF(NNMAT.EQ.0) GO TO 555
  367.       NRIGE=7
  368.       NRIGEL=NNMAT
  369.       SEGINI MRIGID
  370.       ICHOLE=0
  371.       IMGEO1=0
  372.       IMGEO2=0
  373.       ISUPEQ=0
  374.       IFORIG=IFOUR
  375.       MTYMAT='RIGIDITE'
  376.       KRIGI=MRIGID
  377. c
  378. cccc on rajoute les noeuds associes aux multiplicateurs au maillage
  379. c
  380.       segact mcoord*mod
  381. ccccccccc   nombre de noeuds
  382.       NBNO=nbpts
  383.       NBNOI=NBNO
  384.       IF(IENSE.EQ.1) NBPTS=NBNO+NNMAT*(NBPOIN-1)
  385.       IF(IDEPL.EQ.1) THEN
  386.          IF(IDIM.EQ.2) NBPTS=NBNO+(2*(NBPOIN-3)+3)
  387.          IF(IDIM.EQ.3) NBPTS=NBNO+(3*(NBPOIN-4)+6)
  388.          IF(NBPOIN.EQ.2) NBPTS=NBNO +1
  389.          IF(ABS(NIFOUR).EQ.1)NBPTS = NBPTS + NBPOIN + 1
  390.       ENDIF
  391.       SEGADJ MCOORD
  392. c
  393. cccc  boucle sur les rigidites elementaires
  394. c   si idepl=1 et irota=1 la 1ere rigidite elementaire est associee a
  395. c   l option depl et les suivantes a l option rota
  396. c
  397.  
  398.       IF (IENSE.EQ.1) NBRELA=NBPOIN-1
  399.       I7=0
  400.       DO 7 IAA=1,NNMAT
  401.          IF (IDEPL.EQ.1) THEN
  402.             IF (IAA.EQ.1) THEN
  403.                IF(IDIM.EQ.2) NBRELA=2*(NBPOIN-3)+3
  404.                IF(IDIM.EQ.3) NBRELA=3*(NBPOIN-4)+6
  405.                IF(NBPOIN.EQ.2) NBRELA=1
  406.             ELSE IF (IAA.eq.2) THEN
  407. * fourier 1 ur=-ut
  408.                NBRELA=NBPOIN
  409. * fourier 1 relation d'analogie sur uz
  410.             ELSE IF (IAA.eq.3) THEN
  411.                nbrela = 1
  412.             ENDIF
  413.          ENDIF
  414. *      print*,'idepl=',idepl
  415. *      print*,'iaa=',iaa
  416. cccccccc on cree le noeud nbno+1 et nbno+2
  417. cccccccc on les met a l origine
  418.          DO 5 IA=1,NBRELA
  419.             XCOOR(NBNOI*(IDIM+1)+1)=0.D0
  420.             XCOOR(NBNOI*(IDIM+1)+2)=0.D0
  421.             IF (IDIM.EQ.3) XCOOR(NBNOI*(IDIM+1)+3)=0.D0
  422.             XCOOR(NBNOI*(IDIM+1)+(IDIM+1))=0.D0
  423.             NBNOI=NBNOI+1
  424.  5       CONTINUE
  425. c
  426. c     on initialise le segment meleme associe aux blocages
  427. c
  428.          NBSOUS=0
  429.          NBREF=0
  430.          NBNN=3
  431.          NBELEM=NBRELA
  432.          SEGINI IPT1
  433.          IRIGEL(1,IAA)=IPT1
  434.          IPT1.ITYPEL=22
  435.          DO 400 I4=1,NBRELA
  436.             IPT1.ICOLOR(I4)=IDCOUL
  437. c punti associati ai moltiplicatori
  438.             IPT1.NUM(1,I4)=NBNO+I7+1
  439.             I7=I7+1
  440.  400     continue
  441. c caricamento nodi pseudo-elementi
  442. * option ensemble
  443.          IF (IENSE.EQ.1) THEN
  444.             DO 420 I9=1,NBRELA
  445.                IPT1.NUM(2,I9)=NUM(1,I9)
  446.                IPT1.NUM(3,I9)=NUM(1,I9+1)
  447.  420        continue
  448. * fourier 1 ur = -ut
  449.          ELSE IF (IDEPL.EQ.1.AND.IAA.EQ.2) THEN
  450.             DO 421 I9=1,NBRELA
  451.                IPT1.NUM(2,I9)=NUM(1,I9)
  452.                IPT1.NUM(3,I9)=NUM(1,I9)
  453.  421        continue
  454.          ELSE IF (IDEPL.EQ.1.AND.IAA.EQ.1) THEN
  455.             I9=0
  456.             IF(NBPOIN.GT.NBMAIT) THEN
  457. c pseudo-elements pour les relations des autres noeuds avec
  458. c         les noeuds maitres
  459.                DO 430 IC=1, NBESCL
  460.                   ICOREC = 0
  461.                   DO 431 IB=1,NBMAIT-1
  462.                      I9=I9+1
  463.                      IPT1.NUM(3,I9)=NUESCL(IC)
  464.                      IPT1.NUM(2,I9)=NUMAIT(IB)
  465. * on teste si le noeud esclave est colineaire avec les 2 maitres
  466.                      IF (IB.EQ.2) THEN
  467.                         IREFES = (NUESCL(IC)-1)*(IDIM+1)
  468.                         XX13=XCOOR(IREFES+1)-XCOOR(IREFP1+1)
  469.                         YY13=XCOOR(IREFES+2)-XCOOR(IREFP1+2)
  470.                         PP13N=SQRT(XX13**2+YY13**2)
  471.                         TTEST=X12*XX13+Y12*YY13
  472.                         IF (IDIM.EQ.3) THEN
  473.                            ZZ13 = XCOOR(IREFES+3)-XCOOR(IREFP1+3)
  474.                            PP13N=PP13N+ZZ13**2
  475.                            TTEST=TTEST+Z12*ZZ13
  476.                         ENDIF
  477.                         IF(ABS(TTEST).GT.0.9997D0*(P12N*PP13N)) THEN
  478.                            IPT1.NUM(2,I9)=NUMAIT(IB+1)
  479.                            ICOREC = 1
  480.                         ENDIF
  481.                      ELSE IF (IB.EQ.3) THEN
  482.                         IF (ICOREC.EQ.1) THEN
  483.                            IPT1.NUM(2,I9)=NUMAIT(IB+1)
  484. * on teste si le noeud esclave est coplanaire avec les 3 maitres
  485.                         ELSE
  486.                            TTEST=XN*XX13+YN*YY13+ZN*ZZ13
  487.                            IF(ABS(TTEST).LT.0.025D0*(P12N
  488.      $                          *P13N))IPT1.NUM(2,I9)=NUMAIT(IB+1)
  489.                         ENDIF
  490.                      ENDIF
  491.  431              continue
  492.  430           continue
  493.             ENDIF
  494. c pseudo-elements pour les relations entre les noeuds maitres
  495.             DO 432 IB=2,NBMAIT
  496.                I9=I9+1
  497.                IPT1.NUM(2,I9)=NUMAIT(1)
  498.                IPT1.NUM(3,I9)=NUMAIT(IB)
  499.  432        continue
  500.             IF(NBMAIT.GE.3) THEN
  501.                DO 433 IB=3,NBMAIT
  502.                   I9=I9+1
  503.                   IPT1.NUM(2,I9)=NUMAIT(2)
  504.                   IPT1.NUM(3,I9)=NUMAIT(IB)
  505.  433           continue
  506.             ENDIF
  507.             IF (NBMAIT.GE.4) THEN
  508.                DO 434 IB=4,NBMAIT
  509.                   I9=I9+1
  510.                   IPT1.NUM(2,I9)=NUMAIT(3)
  511.                   IPT1.NUM(3,I9)=NUMAIT(IB)
  512.  434           continue
  513.             ENDIF
  514.          ENDIF
  515. * fourier 1, relation d'analogie sur uz
  516. * on cherche les 2 noeuds maitres avec xmax et xmin
  517. * pour etre sur de ne pas avoir choisi 2 noeuds ayan la
  518. * meme coordonnee x
  519.          if ((idepl.eq.1).and.(iaa.eq.3)) then
  520.             xx = xcoor(irefp1+1)
  521.             x1 = xx
  522.             ipp1 = ip1
  523.             if (xcoor(irefp2+1).gt.xx) then
  524.                x1 =xcoor(irefp2+1)
  525.                ipp1 = ip2
  526.             else if (xcoor(irefp3+1).gt.xx) then
  527.                x1 =xcoor(irefp3+1)
  528.                ipp1 = ip3
  529.             endif
  530.             x2 = xx
  531.             ipp2 =ip1
  532.             if (xcoor(irefp2+1).lt.xx) then
  533.                x2 =xcoor(irefp2+1)
  534.                ipp2 = ip2
  535.             else if (xcoor(irefp3+1).lt.xx) then
  536.                x2 =xcoor(irefp3+1)
  537.                ipp2 = ip3
  538.             endif
  539.             ipt1.num(2,1)=ipp2
  540.             ipt1.num(3,1)=ipp1
  541.          endif
  542. *  correction on place les pts mult au barycentre des noeuds associes
  543. *   a chaque relation
  544.          DO 425 IA=1,NBRELA
  545.             IREF1=(IPT1.NUM(1,IA)-1)*(IDIM+1)
  546.             IREF2=(IPT1.NUM(2,IA)-1)*(IDIM+1)
  547.             IREF3=(IPT1.NUM(3,IA)-1)*(IDIM+1)
  548.             XCOOR(IREF1+1)=(XCOOR(IREF2+1)+XCOOR(IREF3+1))/2
  549.             XCOOR(IREF1+2)=(XCOOR(IREF2+2)+XCOOR(IREF3+2))/2
  550.             IF(IDIM.EQ.3) XCOOR(IREF1+3)=(XCOOR(IREF2+3)+XCOOR(IREF3+3))
  551.      $           /2
  552.  425     CONTINUE
  553. cccccccc   on vient de crer les segment meleme resultat
  554. cccccccc   on va creer les raideurs
  555.          IRIGEL(2,IAA)=0
  556.          IRIGEL(5,IAA)=NIFOUR
  557.          IRIGEL(6,IAA)=0
  558.          IRIGEL(7,IAA)=0
  559.          IF(IENSE.EQ.1.OR.IAA.GT.1) NLIGRE=3
  560.          IF(IDEPL.EQ.1.AND.IAA.EQ.1) NLIGRE=2*IDIM+1
  561.          NLIGRP=NLIGRE
  562.          NLIGRD=NLIGRE
  563.          SEGINI DESCR
  564.          IRIGEL(3,IAA)=DESCR
  565. c         on remplit le tableau des descripteurs de rig
  566.  
  567.          LISINC(1)='LX'
  568.          LISDUA(1)='FLX'
  569.          NOELEP(1)=1
  570.          NOELED(1)=1
  571.          IF(IENSE.EQ.1) THEN
  572.             JAA=2*(IAA-1)
  573.             LISINC(2)=MOTDDL(JAA+1)
  574.             LISINC(3)=LISINC(2)
  575.             LISDUA(2)=MOTDDL(JAA+2)
  576.             LISDUA(3)=LISDUA(2)
  577.             NOELEP(2)=2
  578.             NOELEP(3)=3
  579.             NOELED(2)=2
  580.             NOELED(3)=3
  581.          ELSE IF (IDEPL.EQ.1.AND.IAA.EQ.1) THEN
  582.             DO 250 I=1,IDIM
  583.                JAA=2*(I-1)
  584.                J=1+2*(I-1)
  585.                LISINC(J+1)=NOMINC(JAA+1)
  586.                LISINC(J+2)=LISINC(J+1)
  587.                LISDUA(J+1)=NOMINC(JAA+2)
  588.                LISDUA(J+2)=LISDUA(J+1)
  589.                NOELEP(J+1)=2
  590.                NOELEP(J+2)=3
  591.                NOELED(J+1)=2
  592.                NOELED(J+2)=3
  593.  250        continue
  594.  
  595. * fourier 1 ur = -ut
  596.          ELSE IF (IDEPL.EQ.1.AND.IAA.EQ.2) THEN
  597.             LISINC(2)=MODEPL(4)
  598.             LISINC(3)=MODEPL(6)
  599.             LISDUA(2)=MODEDU(4)
  600.             LISDUA(3)=MODEDU(6)
  601.             NOELEP(2)=2
  602.             NOELEP(3)=3
  603.             NOELED(2)=2
  604.             NOELED(3)=3
  605.          ELSE IF (IDEPL.EQ.1.AND.IAA.EQ.3) THEN
  606.             LISINC(2)=MODEPL(3)
  607.             LISINC(3)=MODEPL(3)
  608.             LISDUA(2)=MODEDU(3)
  609.             LISDUA(3)=MODEDU(3)
  610.             NOELEP(2)=2
  611.             NOELEP(3)=3
  612.             NOELED(2)=2
  613.             NOELED(3)=3
  614.          ENDIF
  615.          SEGDES DESCR
  616. c     
  617.          NELRIG=NBRELA
  618.          SEGINI xMATRI
  619.          IRIGEL(4,IAA)=xMATRI
  620.          COERIG(IAA)=1.D0
  621. c
  622. c  remplissage de  RE
  623.  
  624.          IF (IENSE.EQ.1.OR.IAA.GT.1) THEN
  625.             do iou=1,nelrig
  626. *      SEGINI XMATRI
  627. *      IXMATR=XMATRI
  628.                RE(1,1,iou)=  0.D0
  629.                RE(2,1,iou)= 1.D0
  630.                IF (IENSE.EQ.1) THEN
  631.                   RE(3,1,iou)= -1.D0
  632. * ur = - ut
  633.                ELSE IF (IAA.EQ.2) THEN
  634.                   RE(3,1,iou)= 1.D0
  635. * analogie sur uz en fourier 1
  636.                ELSE
  637.                   RE(3,1,iou)= -x2/x1
  638.                ENDIF
  639.                DO 777 I=1,1
  640.                   DO 666 J=2,NLIGRE
  641.                      RE(I,J,iou)=RE(J,I,iou)
  642.  666              CONTINUE
  643.  777           CONTINUE
  644.             enddo
  645. *      SEGDES XMATRI
  646. *      DO 3 IA=1,NELRIG
  647. *      IMATTT(IA)=IXMATR
  648. *    3 CONTINUE
  649.          ELSE IF (IDEPL.EQ.1.AND.IAA.EQ.1) THEN
  650.             DO 740 IA=1,NELRIG
  651.                cmax=0.D0
  652.                DO 741 I=1,IDIM
  653.                   IREF2=(IPT1.NUM(2,IA)-1)*(IDIM+1)
  654.                   IREF3=(IPT1.NUM(3,IA)-1)*(IDIM+1)
  655.                   COEF(I)=XCOOR(IREF2+I)-XCOOR(IREF3+I)
  656.                   cmax=max(cmax,abs(coef(i)))
  657.  741           continue
  658. *      SEGINI XMATRI
  659. *      IMATTT(IA)=XMATRI
  660.                RE(1,1,ia)= 0.D0
  661. *  normalisation par le terme max de la relation
  662.                DO 760 I1=2,NLIGRE
  663.                   I4=INT(I1/2)
  664.                   COEF1=COEF(I4)
  665.                   IF((2*I4).LT.I1) COEF1=-COEF(I4)
  666.                   RE(I1,1,ia)=COEF1/cmax
  667.                   RE(1,I1,ia)=COEF1/cmax
  668.  760           continue
  669. *     SEGDES XMATRI
  670.  740        continue
  671.          ENDIF
  672.          SEGDES xMATRI
  673. c     
  674.          SEGDES IPT1
  675.     7 CONTINUE
  676.  
  677. c
  678. *  nettoyage eventuel des relations
  679.       call relasi(mrigid)
  680.       SEGDES MRIGID
  681.       CALL ECROBJ('RIGIDITE',MRIGID)
  682.  999  continue
  683.       IF (IDEPL.EQ.1) THEN
  684.          SEGSUP NUMAIT
  685.          IF(NBESCL.GT.0) SEGSUP NUESCL
  686.       ENDIF
  687.  998  continue
  688.       SEGDES MELEME
  689.  555  CONTINUE
  690.       SEGSUP MOTDDV,NOMINV
  691.       RETURN
  692.       END
  693.  
  694.  
  695.  
  696.  
  697.  
  698.  
  699.  
  700.  
  701.  
  702.  
  703.  
  704.  
  705.  
  706.  
  707.  
  708.  
  709.  

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