Télécharger la plaquette Cast3M


fr en

Télécharger rela1.eso

Retour à la liste

Numérotation des lignes : 

rela1

  1. C RELA1     SOURCE    PV090527  24/02/19    21:15:02     11841          
  2.       SUBROUTINE RELA1
  3. C
  4. C     RELAZIONE LINEARE TRA DDL
  5. C
  6. C     ICONR  NUMERO OGGETTI MELEME DELLA RELAZIONE
  7. C     ITYESR(NBSOUS) LISTA TIPI ELEMENTI CONTENUTI IN MELEME-N
  8. C     MUNESR(NBSOUS)  NUMERO ELEMENTI IN OGNI SOTTOSTRUTTURA DI MELEME-N
  9. C     LNODSR(NNR) LISTA NODI MELEME-N
  10. C     IPOR1(N)  PUNTATORE SU MWREL1
  11. C     IPOR2(N)  PUNTATORE SU MWREL2
  12. C     IPOR3(N)  PUNTATORE SU MWREL3
  13. C     COEFR(N)  COEFFICIENTI RELAZIONE
  14. C     INCREL(N) IDENTIFICATORE NOME DDL
  15. C
  16.       IMPLICIT INTEGER(I-N)
  17.       IMPLICIT REAL*8(A-H,O-Z)
  18.  
  19. -INC SMELEME
  20. -INC SMCOORD
  21. -INC SMRIGID
  22. -INC TMTRAV
  23. -INC SMCHPOI
  24.  
  25. -INC PPARAM
  26. -INC CCOPTIO
  27. -INC CCREEL
  28. -INC CCGEOME
  29. -INC CCHAMP
  30.  
  31.       SEGMENT /MWGGM1/(IPORE1(0))
  32.       SEGMENT /MWGGM2/(IPORE2(0))
  33.       SEGMENT /MWGGM3/(IPORE3(0))
  34.       SEGMENT /MWGGM4/(INCREL(0))
  35.       SEGMENT /MWGGM5/(COEFR(0)*D)
  36.       SEGMENT /MWREL1/(ITYESR(0)),IEL11.MWREL1,
  37.      +     IEL21.MWREL1,iel01.MWREL1
  38.       SEGMENT /MWREL2/(MUNESR(0)),IEL12.MWREL2,
  39.      +     IEL22.MWREL2,iel02.MWREL2
  40.       SEGMENT /MWREL3/(LNODSR(0)),IEL13.MWREL3,
  41.      +     IEL23.MWREL3
  42.       CHARACTER*4 MOTPV(3)
  43.       CHARACTER*4 MOCORI(7)
  44.       CHARACTER*4 MOTPM(2)
  45.       CHARACTER*4 MOTDDL(3)
  46.       CHARACTER*4 MOTBLO(3)
  47.       CHARACTER*4 MODEPL(8)
  48.       CHARACTER*4 MOROTA(5)
  49.       CHARACTER*4 MODUAL(1)
  50.       CHARACTER*8 MILMOT
  51.       DIMENSION XNOR(3)
  52.       DATA MOTPV /'MINI','MAXI','FROT'/
  53.       DATA MOTPM /'+   ','-   '/
  54.       DATA MOCORI/'CORI','ENSE','ACCR','GLIS','BARY','MILI','TUYA'/
  55. C
  56.       DATA MOTBLO/'DEPL','ROTA','DIRE'/
  57.       DATA MODEPL/'UX  ','UY  ','UZ  ','UR  ','UZ  ','UT  ',
  58.      &     'ALFA','BETA'/
  59.       DATA MOROTA/'RX  ','RY  ','RZ  ','RT  ','RS  '/
  60.       DATA MODUAL/'DUAL'/
  61. c
  62.       iel01=0
  63.       iel02=0
  64. C
  65. C     EST-CE UNE RELATION DE CORPS RIGIDE ?
  66. C
  67.  
  68.       CALL LIRMOT(MOCORI,7,ICORI,0)
  69. *      write(ioimp,*) 'rela1, icori=',icori
  70.       IF (ICORI.NE.0) THEN
  71. c     BARY
  72.          if (icori.eq.5) then
  73.             call relaba
  74. c     MILI
  75.          ELSE IF (ICORI.EQ.6) THEN
  76.             CALL RELAMI
  77. c     TUYA
  78.          elseif(icori.eq.7) then
  79.             call reltuy
  80. c     CORI, ENSE, ACCRO, GLIS
  81.          else
  82.             CALL RELA2(MOCORI(ICORI))
  83.          endif
  84.          RETURN
  85.       ENDIF
  86. C
  87. C     EST-CE UNE CONDITION UNILATERALE ?
  88. C
  89.       NILATE=0
  90.       CALL LIRMOT (MOTPV,3,IPO,0)
  91.       IF(IPO.EQ.1) NILATE=-1
  92.       IF(IPO.EQ.2) NILATE=+1
  93.       IF(IPO.EQ.3) NILATE=+2
  94. C
  95. C     INITIALISATIONS
  96. C
  97.       COEX=1.D0
  98. C
  99. C     Deformations planes ou contraintes planes ou defo. plane gene :
  100.       IF (IFOUR.EQ.-1.OR.IFOUR.EQ.-2.OR.IFOUR.EQ.-3) THEN
  101.          LDEPL=2
  102.          IADEPL=0
  103.          LROTA=1
  104.          IAROTA=2
  105. C     Axisymetrique :
  106.       ELSE IF (IFOUR.EQ.0) THEN
  107.          LDEPL=2
  108.          IADEPL=3
  109.          LROTA=1
  110.          IAROTA=3
  111. C     Fourier :
  112.       ELSE IF (IFOUR.EQ.1) THEN
  113.          LDEPL=3
  114.          IADEPL=3
  115.          LROTA=2
  116.          IAROTA=2
  117. C     Tridimensionnel :
  118.       ELSE IF (IFOUR.EQ.2) THEN
  119.          LDEPL=3
  120.          IADEPL=0
  121.          LROTA=3
  122.          IAROTA=0
  123. C     Massif 1D (IDIM=1) :
  124.       ELSE IF (IFOUR.GE.3.AND.IFOUR.LE.15) THEN
  125.          LDEPL=1
  126.          IADEPL=0
  127.          IF (IFOUR.GE.12) IADEPL=3
  128.          LROTA=0
  129.          IAROTA=0
  130. C     Autres cas :
  131.       ELSE
  132.          LDEPL=0
  133.          IADEPL=0
  134.          LROTA=0
  135.          IAROTA=0
  136.       ENDIF
  137. C
  138.       SEGINI MWGGM1
  139.       SEGINI MWGGM2
  140.       SEGINI MWGGM3
  141.       SEGINI MWGGM4
  142.       SEGINI MWGGM5
  143. *
  144. *     LECTURE EVENTUELLE D'UN MODELLE
  145. *
  146.  
  147.       CALL LIROBJ('MMODEL',IPOMOD,0,IREMOD)
  148.       IF(IREMOD.NE.0) THEN
  149.          CALL RELMOD(IPOMOD)
  150.          RETURN
  151.       ENDIF
  152. *
  153. *     LECTURE EVENTUELLE D'UN CHPOINT
  154. *
  155.       CALL LIROBJ('CHPOINT',IPOCHP,0,IRECHP)
  156.       IF(IRECHP.NE.0) THEN
  157.         CALL ACTOBJ('CHPOINT',IPOCHP,1)
  158. *  lecture eventuelle mot cle DUAL
  159.          IRECHD=0
  160.          CALL LIRMOT (MODUAL,1,IRECHD,0)
  161.          IF(IRECHD.EQ.1) GOTO 500
  162. *
  163. *     ON MET LE CHPOINT SOUS FORME DE TABLEAU DE TRAVAIL
  164. *
  165.          CALL TRACHP(IPOCHP,MTRAV)
  166.          SEGACT MTRAV
  167.          DO 1990 I=1,IBIN(/1)
  168.             DO 19901 J=1,IBIN(/2)
  169.                IF(IBIN(I,J).EQ.0) GO TO 19901
  170.                CALL PLACE(NOMDD,LNOMDD,IPO,INCO(I))
  171.                if (ipo.eq.0) then
  172.                  CALL PLACE(NOMDU,LNOMDU,IPO,INCO(I))
  173.                  ipo = -ipo
  174.                endif
  175.                IF(IPO.EQ.0) THEN
  176.                   moterr(1:4)=inco(i)
  177.                   SEGSUP MTRAV
  178.                   call erreur(108)
  179.                   GO TO 559
  180.                ENDIF
  181.                SEGINI MWREL1
  182.                SEGINI MWREL2
  183.                SEGINI MWREL3
  184.                IPORE1(**)=MWREL1
  185.                IPORE2(**)=MWREL2
  186.                IPORE3(**)=MWREL3
  187.                COEFR(**)=BB(I,J)
  188.                INCREL(**)=IPO
  189.                ITYESR(**)=1
  190.                MUNESR(**)=1
  191.                LNODSR(**)=IGEO(J)
  192. 19901       CONTINUE
  193.  1990    CONTINUE
  194.          SEGSUP MTRAV
  195.          GO TO 300
  196.       ENDIF
  197.  200  CONTINUE
  198. C
  199. C     LETTURA OPERANDI
  200. C
  201.       SEGINI MWREL1
  202.       SEGINI MWREL2
  203.       SEGINI MWREL3
  204.       IPORE1(**)=MWREL1
  205.       IPORE2(**)=MWREL2
  206.       IPORE3(**)=MWREL3
  207.  1320 CONTINUE
  208.       CALL LIRREE(XR,0,IRETOF)
  209.       IF(IRETOF.EQ.0) THEN
  210.          CALL LIRENT(IR,0,IRETOI)
  211.          XR=1.D0
  212.          IF(IRETOI.NE.0) XR=DBLE(IR)
  213.       ENDIF
  214.       CALL LIRMOT(NOMDD,LNOMDD,IPO,0)
  215.       IF(IPO.NE.0) THEN
  216.          COEFR(**)=XR*COEX
  217.          INCREL(**)=IPO
  218.       ELSE
  219. *
  220. *     ON REGARDE SI IL Y A DEPL OU ROTA SUIVI DE DIRECTION
  221. *
  222. C     En DIMENSION 1, le mot-cle 'DIRE' est interdit.
  223.          IF (IDIM.EQ.1) THEN
  224.             INTERR(1)=IDIM
  225.             MOTERR(1:4)=MOTBLO(3)
  226.             CALL ERREUR(971)
  227.             GOTO 559
  228.          ENDIF
  229.          CALL LIRMOT(MOTBLO,2,IMOT,1)
  230.          IF (IMOT.EQ.0) GOTO 559
  231.          IF (IMOT.EQ.1) THEN
  232.             IBDDL=LDEPL
  233.             DO 4481 IA=1,IBDDL
  234.                MOTDDL(IA)=MODEPL(IADEPL+IA)
  235.  4481       CONTINUE
  236.          ENDIF
  237.          IF (IMOT.EQ.2) THEN
  238.             IBDDL=LROTA
  239.             DO 4482 IA=1,IBDDL
  240.                MOTDDL(IA)=MOROTA(IAROTA+IA)
  241.  4482       CONTINUE
  242.          ENDIF
  243.          CALL LIRMOT(MOTBLO(3),1,IMOT,1)
  244.          IF(IMOT.EQ.0) GO TO 559
  245.          CALL LIROBJ('POINT',KPOINT,1,IRETOU)
  246.          IF(IRETOU.EQ.0) GO TO 559
  247.          YL=0.D0
  248.          XPREC=SQRT(XPETIT/XZPREC)
  249.          segact mcoord
  250.          DO 4484 IA=1,IDIM
  251.             XVAL    =XCOOR((KPOINT-1)*(IDIM+1)+IA)
  252.             XNOR(IA)=XVAL
  253.             IF(ABS(XVAL) .GT. XPREC)THEN
  254.               YL=YL+XVAL**2
  255.             ENDIF
  256.  4484    CONTINUE
  257.          IF (YL .EQ. 0.D0) THEN
  258.             CALL ERREUR(239)
  259.             GOTO 559
  260.          ENDIF
  261.          YL=1.D0/SQRT(YL)
  262.          DO 4485 IA=1,IDIM
  263.             XNOR(IA)=XNOR(IA)*YL
  264.  4485    CONTINUE
  265. *
  266. *     ON LIT LE MAILLAGE
  267. *
  268.          CALL LIROBJ('POINT   ',MILPOI,0,IRETOU)
  269.          IF(IRETOU.NE.0) THEN
  270.             MILMOT='POINT   '
  271.          ELSE
  272.             CALL LIROBJ('MAILLAGE',MILPOI,1,IRETOU)
  273.             IF(IRETOU.EQ.0)  GO TO 559
  274.             MILMOT='MAILLAGE'
  275.          ENDIF
  276. *
  277. *     PUIS ON REMET DES OBJETS DANS LA PILE
  278. *
  279.          DO 4477 IB=1,IBDDL
  280.             IF(IB.GE.2.AND.COEX.EQ. 1.D0) CALL ECRCHA(MOTPM(1))
  281.             IF(IB.GE.2.AND.COEX.EQ.-1.D0) CALL ECRCHA(MOTPM(2))
  282.             CALL ECROBJ(MILMOT,MILPOI)
  283.             CALL ECRCHA(MOTDDL(IB))
  284.             IF(IBDDL.EQ.1) THEN
  285.                XRCOO=XR
  286.             ELSE
  287.                XRCOO=XR*XNOR(IB)
  288.             ENDIF
  289.             CALL ECRREE(XRCOO)
  290.  4477    CONTINUE
  291.          GO TO 1320
  292.       ENDIF
  293. C
  294.       CALL LIROBJ('POINT   ',KPOINT,0,IRETOU)
  295.       IF(IRETOU.NE.0) GO TO 110
  296.       CALL LIROBJ('MAILLAGE',KOBJET,1,IRETOU)
  297.       IF(IRETOU.EQ.0)  GO TO 559
  298.       MELEME=KOBJET
  299.       SEGACT MELEME
  300.       NBSOUS=LISOUS(/1)
  301.       IF(NBSOUS.EQ.0) GO TO 120
  302. C     OBJET COMPLEXE
  303.       NNR=1
  304.       DO 130 IS=1,NBSOUS
  305.          IPT1=LISOUS(IS)
  306.          SEGACT IPT1
  307.          ITYESR(**)=IPT1.ITYPEL
  308.          NBNN      =IPT1.NUM(/1)
  309.          NBELEM    =IPT1.NUM(/2)
  310.          MUNESR(**)=IPT1.NUM(/2)
  311.          IF(NNR.EQ.1)LNODSR(**)=IPT1.NUM(1,1)
  312.          DO 140  I1=1,NBELEM
  313.             DO 1401  I2=1,NBNN
  314.                DO 150  I3=1,NNR
  315.                   IF(LNODSR(I3).EQ.IPT1.NUM(I2,I1)) GO TO 1401
  316.  150           CONTINUE
  317.                NNR=NNR+1
  318.                LNODSR(**)=IPT1.NUM(I2,I1)
  319.  1401       CONTINUE
  320.  140     CONTINUE
  321.          SEGDES IPT1
  322.  130  CONTINUE
  323.       SEGDES MWREL1
  324.       SEGDES MWREL2
  325.       SEGDES MWREL3
  326.       SEGDES MELEME
  327.       GO TO 160
  328. C
  329. C     OBJET SIMPLE
  330. C
  331.  120  CONTINUE
  332.       ITYESR(**)=ITYPEL
  333.       NBNN      =NUM(/1)
  334.       NBELEM    =NUM(/2)
  335.       MUNESR(**)=NUM(/2)
  336.       NNR=0
  337.       DO 170  I1=1,NBELEM
  338.          DO 1701  I2=1,NBNN
  339.             DO 180  I3=1,NNR
  340.                IF(LNODSR(I3).EQ.NUM(I2,I1)) GO TO 1701
  341.  180        CONTINUE
  342.             NNR=NNR+1
  343.             LNODSR(**)=NUM(I2,I1)
  344.  1701    CONTINUE
  345.  170  CONTINUE
  346.       SEGDES MWREL1
  347.       SEGDES MWREL2
  348.       SEGDES MWREL3
  349.       SEGDES MELEME
  350.       GO TO 160
  351.  110  CONTINUE
  352. C
  353. C     OBJET POINT
  354. C
  355.       ITYESR(**)=1
  356.       MUNESR(**)=1
  357.       LNODSR(**)=KPOINT
  358.       SEGDES MWREL1
  359.       SEGDES MWREL2
  360.       SEGDES MWREL3
  361. C     FINE OPERANDO RELAZIONE
  362.  160  CONTINUE
  363.       ICONR=IPORE1(/1)+1
  364. C     LETTURA + O -
  365.       IF(ICONR.EQ.2) CALL LIRMOT(MOTPM,2,IPO,0)
  366.       IF (IPO.EQ.0)  GO TO 300
  367. *     LIRMOT(MOTPM,2,IPO,1) goto 559
  368.       IF(ICONR.GT.2) CALL LIRMOT(MOTPM,2,IPO,0)
  369.       IF (IPO.EQ.0)  GO TO 300
  370.       COEX=1.D0
  371.       IF (IPO.EQ.2) COEX=-1.D0
  372. C     SI CERCA DI LEGGERE UN NUOVO OPERANDO DELLA RELAZIONE
  373.       GO TO  200
  374. C
  375. C     VERIFICA CONGRUENZA OPERANDI
  376. C
  377.  300  CONTINUE
  378. C
  379.       ICONR=IPORE1(/1)
  380.       IEL11=IPORE1(1)
  381.       IEL12=IPORE2(1)
  382.       IEL13=IPORE3(1)
  383. *     on autorise maintenant le point unique qui va etre applique sur
  384. *     chaque relation
  385. *     old      ityes=0
  386. *     old      munes=0
  387.       nsor=0
  388.       nnr=0
  389.       do 305 io=1,iconr
  390. *     write(ioimp,*) 'io=',io
  391.          iel11=ipore1(io)
  392.          iel12=ipore2(io)
  393.          iel13=ipore3(io)
  394.          segact iel11,iel12,iel13
  395.          nsor1=iel11.ityesr(/1)
  396.          nsor=max(nsor,nsor1)
  397.          NNR1=IEL13.LNODSR(/1)
  398.          nnr=max(nnr,nnr1)
  399. *     old         do 315 i1=1,nsor1
  400.          if (io.eq.1) then
  401.             segini,iel01=iel11
  402.             segini,iel02=iel12
  403.          else
  404.             do 315 i1=1,nsor1
  405.                ityes=iel01.ityesr(i1)
  406.                iel01.ityesr(i1)=max(iel11.ityesr(i1),ityes)
  407. *     write(ioimp,*) 'i1=',i1,' ityes=',ityes
  408. *     $              ,' iel11.ityesr(i1)',iel11.ityesr(i1)
  409.                munes=iel02.munesr(i1)
  410.                iel02.munesr(i1)=max(iel12.munesr(i1),munes)
  411. *     write(ioimp,*) 'i1=',i1,' munes=',munes
  412. *     $              ,' iel12.munesr(i1)',iel12.munesr(i1)
  413.  315        continue
  414.          endif
  415.  305  continue
  416.       DO 310 IO=1,ICONR
  417. *     write(ioimp,*) 'io=',io
  418.          IEL21=IPORE1(IO)
  419.          IEL22=IPORE2(IO)
  420.          IEL23=IPORE3(IO)
  421.          NSOR2=IEL21.ITYESR(/1)
  422. *     write(ioimp,*) 'nsor,nsor2=',nsor,nsor2
  423.          IF(NSOR.NE.NSOR2.and.nsor2.ne.1) GO TO 556
  424.          DO 320 I2=1,NSOR2
  425.             ityes=iel01.ityesr(i2)
  426.             munes=iel02.munesr(i2)
  427. *     write(ioimp,*) 'i2=',i2,' ityes=',ityes,'munes=',munes
  428. *     write(ioimp,*) 'iel21.ityesr=',iel21.ityesr(i2)
  429. *     $           ,' IEL22.MUNESR(I2)=',IEL22.MUNESR(I2)
  430.             if (iel21.ityesr(i2).eq.1.and.iel22.munesr(i2).eq.1.and.
  431.      >           nsor2.eq.1) goto 320
  432.             IF(IEL21.ITYESR(I2).NE.ityes) GO TO 556
  433.             IF(IEL22.MUNESR(I2).NE.munes) GO TO 556
  434.  320     CONTINUE
  435.          NNR2=IEL23.LNODSR(/1)
  436. *     write(ioimp,*) 'nnr,nnr2=',nnr,nnr2
  437.          IF(NNR2.NE.NNR.and.nnr2.ne.1) GO TO 556
  438.  310  CONTINUE
  439.       if (iel01.ne.0) SEGSUP,IEL01
  440.       if (iel02.ne.0) SEGSUP,IEL02
  441.       NBRELA=NNR
  442. C     OPERANDI  CONGRUENTI  FINE VERIFICA
  443. C
  444. C     CARICAMENTO MATRICE  RIGIDEZZA
  445. C
  446.       NRIGEL=1
  447.       SEGINI MRIGID
  448.       ICHOLE=0
  449.       IMGEO1=0
  450.       IMGEO2=0
  451.       ISUPEQ=0
  452.       IFORIG=IFOUR
  453.       COERIG(1)=1.D0
  454.       MTYMAT='RIGIDITE'
  455.       KRIGI=MRIGID
  456. C
  457. C     ON INITIALISE LE SEGMENT MELEME ASSOCIE AUX BLOCAGES
  458. C
  459.       SEGACT MCOORD*MOD
  460.       NBPTSO=nbpts
  461. *      write(ioimp,*) 'nbptso=',nbptso
  462.       NBPTS=NBPTSO+NBRELA
  463.       SEGADJ MCOORD
  464.       NBSOUS=0
  465.       NBREF=0
  466.       NBNN=1+ICONR
  467.       NBELEM=NBRELA
  468.       SEGINI IPT2
  469.       IRIGEL(1,1)=IPT2
  470.       IPT2.ITYPEL=22
  471. C     SI CREANO UNO PUNTO PER OGNI RELAZIONE
  472.       DO 400 I4=1,NBRELA
  473.          IPT2.ICOLOR(I4)=IDCOUL
  474.          IPTS=NBPTSO+I4
  475.          JPTS=(IPTS-1)*(IDIM+1)
  476.          DO 410 iidim=1,idim+1
  477.             XCOOR(JPTS+iidim)=0.D0
  478.  410     CONTINUE
  479. C     PUNTI ASSOCIATI AI MOLTIPLICATORI
  480.          IPT2.NUM(1,I4)=IPTS
  481.  400  CONTINUE
  482. C     CARICAMENTO NODI PSEUDO-ELEMENTI
  483.       DO 420 I8=1,ICONR
  484.          MWREL3=IPORE3(I8)
  485.          I10=I8+1
  486.          SEGACT MWREL3
  487.          LNOMAX=LNODSR(/1)
  488.          DO 430 I9=1,NBRELA
  489.             NN=LNODSR(min(I9,lnomax))
  490.             NPN=(NN-1)*(IDIM+1)
  491.             IPT2.NUM(I10,I9)=NN
  492.             NPL1=(IPT2.NUM(1,I9)-1)*(IDIM+1)
  493.             do 4301 iidim=1,idim+1
  494.                XCOOR(NPL1+iidim)=XCOOR(NPL1+iidim)+XCOOR(NPN+iidim)
  495.  4301       continue
  496.  430     CONTINUE
  497.          SEGDES MWREL3
  498.  420  CONTINUE
  499. C
  500. C     COORDINATE DEI BARICENTRI ASSOCIATI ALLE RELAZIONI
  501. C
  502. *      write(ioimp,*) 'iconr,nbrela,idim=',iconr,nbrela,idim
  503.       DO 425 I4=1,NBRELA
  504.          NPL1=(IPT2.NUM(1,I4)-1)*(IDIM+1)
  505.          do 4251 iidim=1,idim+1
  506.             XCOOR(NPL1+iidim)=XCOOR(NPL1+iidim)/ICONR
  507.  4251    continue
  508.  425  CONTINUE
  509.  
  510.       SEGDES IPT2
  511.       IRIGEL(2,1)=0
  512.       IRIGEL(5,1)=NIFOUR
  513.       IRIGEL(6,1)=NILATE
  514.       NLIGRP=ICONR+1
  515.       NLIGRD=NLIGRP
  516.       SEGINI DESCR
  517.       IRIGEL(3,1)=DESCR
  518.       LISINC(1)='LX'
  519.       LISDUA(1)='FLX'
  520.       NOELEP(1)=1
  521.       NOELED(1)=1
  522.       DO 700 I1=1,ICONR
  523.          I2=I1+1
  524.          I3=ABS(INCREL(I1))
  525.          LISINC(I2)=NOMDD(I3)
  526.          LISDUA(I2)=NOMDU(I3)
  527.          NOELEP(I2)=I2
  528.          NOELED(I2)=I2
  529.  700  CONTINUE
  530.       SEGDES DESCR
  531.       NELRIG=NBRELA
  532.       SEGINI xMATRI
  533.       IRIGEL(4,1)=xMATRI
  534. *     LVAL=(NBNN*NBNN+NBNN)/2
  535.       NLIGRP=NBNN
  536.       NLIGRD=NBNN
  537. *     SEGINI XMATRI
  538.       DO 740 I6=1,NELRIG
  539. *     IMATTT(I1)=XMATRI
  540. *     740   CONTINUE
  541.          RE(1,1,i6)= 0.D0
  542.          I3=3
  543.          DO 760 I1=2,NBNN
  544.             I4=I1-1
  545.             I2=1
  546.             RE(I1,I2,i6)=COEFR(I4)
  547.             RE(I2,I1,i6)=COEFR(I4)
  548.  760     CONTINUE
  549.  740  continue
  550.       SEGDES XMATRI
  551. *     SEGDES IMATRI
  552.       call relasi(mrigid)
  553.       SEGDES MRIGID
  554.       CALL ECROBJ('RIGIDITE',KRIGI)
  555.  559  CONTINUE
  556.       ICONR=IPORE1(/1)
  557.       DO 558 I1=1,ICONR
  558.          MWREL1=IPORE1(I1)
  559.          MWREL2=IPORE2(I1)
  560.          MWREL3=IPORE3(I1)
  561.          SEGSUP MWREL1
  562.          SEGSUP MWREL2
  563.          SEGSUP MWREL3
  564.  558  CONTINUE
  565.       SEGSUP MWGGM1
  566.       SEGSUP MWGGM2
  567.       SEGSUP MWGGM3
  568.       SEGSUP MWGGM4
  569.       SEGSUP MWGGM5
  570.       RETURN
  571.  556  CONTINUE
  572. C     SEGDES IEL11,IEL12,IEL13,IEL21,IEL22,IEL23
  573.       CALL ERREUR(324)
  574.       GO TO 559
  575. *
  576. *  on arrive en 500 avec la syntaxe chp1 DUAL chp2
  577. *
  578.  500  continue
  579.       CALL LIROBJ('CHPOINT',IPOCHD,1,IRECHD)
  580.       if (ierr.ne.0) return
  581. * ipochp est le chpoin de la relation
  582. * ipochd est le chpoin du vecteur de controle
  583. * decompte des dimensions
  584.       nbsous=0
  585.       nbref=0
  586.       nbelem=1
  587.       mchpoi=ipochp
  588.       segact mchpoi
  589. *  reserver la place pour le lx
  590.       nligrp=1
  591.       nbnp=0
  592.       do ims=1,ipchp(/1)
  593.        msoupo=ipchp(ims)
  594.        segact msoupo
  595.        mpoval=ipoval
  596.        segact mpoval
  597.        nligrp=nligrp+vpocha(/1)*vpocha(/2)
  598.        meleme=igeoc
  599.        segact meleme
  600.        nbnp=nbnp+num(/2)
  601.       enddo
  602. *     write(6,*) ' nligrp nbelep ',nligrp,nbelep
  603.       mchpoi=ipochd
  604.       segact mchpoi
  605.       nbnd=0
  606.       xnormd = 0.d0
  607. *  reserver la place pour le flx
  608.       nligrd=1
  609.       do ims=1,ipchp(/1)
  610.        msoupo=ipchp(ims)
  611.        segact msoupo
  612.        mpoval=ipoval
  613.        segact mpoval
  614.        nligrd=nligrd+vpocha(/1)*vpocha(/2)
  615.        do j=1,vpocha(/2)
  616.        do i=1,vpocha(/1)
  617.         xnormd=max(xnormd,abs(vpocha(i,j)))
  618.        enddo
  619.        enddo
  620.        if (xnormd.lt.1d-5.or.xnormd.gt.1E5) then
  621.          reaerr(1)=xnormd
  622.          call erreur(1118)
  623.          return
  624.        endif
  625.        meleme=igeoc
  626.        segact meleme
  627.        nbnd=nbnd+num(/2)
  628.       enddo
  629. *     write(6,*) ' nligrd nbeled ',nligrd,nbeled
  630. * creation rigidite
  631.       nligrp=nligrp+nligrd-1
  632.       nligrd=nligrp
  633.       segini descr
  634.       nbnn=nbnp+1
  635. *  maillage primal uniquement
  636.       nbnn=nbnd+1
  637. *  maillage dual   uniquement
  638.       nbnn=nbnp+nbnd+1
  639.       nbref=0
  640.       segini meleme
  641.       itypel=22
  642.       nelrig=1
  643.       segini xmatri
  644.       symre=2
  645.       nrigel=1
  646.       segini mrigid
  647.       ICHOLE=0
  648.       IMGEO1=0
  649.       IMGEO2=0
  650.       ISUPEQ=0
  651.       IFORIG=IFOUR
  652.       COERIG(1)=1.D0
  653.       MTYMAT='RIGIDITE'
  654.       irigel(1,1)=meleme
  655.       irigel(3,1)=descr
  656.       irigel(4,1)=xmatri
  657.       IRIGEL(5,1)=NIFOUR
  658.       irigel(6,1)=nilate
  659.       irigel(7,1)=2
  660. *  remplissage maillage descripteur et valeur
  661. *  le premier noeud sera fait a la fin. C'est le multiplicateur de lagrange
  662.       iel=1
  663.       ire=1
  664.       ire=1
  665.       nbnp=1
  666.       nbnd=1
  667.       mchpoi=ipochp
  668.       do ims=1,ipchp(/1)
  669.        msoupo=ipchp(ims)
  670.        segact msoupo
  671.        ipt1=igeoc
  672.        segact ipt1
  673.        mpoval=ipoval
  674.        segact mpoval
  675.        do i=1,ipt1.num(/2)
  676.         iel=iel+1
  677.         num(iel,1)=ipt1.num(1,i)
  678.         nbnp=nbnp+1
  679.        do j=1,nocomp(/2)
  680.         ire=ire+1
  681.         lisinc(ire)=nocomp(j)
  682.         if (lisinc(ire).eq.'LX  ') call erreur(1125)
  683.         noelep(ire)=iel
  684.         re(1,ire,1)=vpocha(i,j)
  685.         call place(nomdd,lnomdd,ipo,nocomp(j))
  686.         if(ipo.eq.0) then
  687.         moterr(1:4)=nocomp(j)
  688.         moterr(5:11)='PRIMALE'
  689.         call erreur(1117)
  690.         endif
  691.         lisdua(ire)=nomdu(ipo)
  692.         noeled(ire)=iel
  693. *       write(6,*) ' 1 ire inc dua ',lisinc(ire),lisdua(ire)
  694.        enddo
  695.        enddo
  696.       enddo
  697.       mchpoi=ipochd
  698.       do ims=1,ipchp(/1)
  699.        msoupo=ipchp(ims)
  700.        segact msoupo
  701.        ipt1=igeoc
  702.        segact ipt1
  703.        mpoval=ipoval
  704.        segact mpoval
  705.        do i=1,ipt1.num(/2)
  706.         iel=iel+1
  707.         num(iel,1)=ipt1.num(1,i)
  708.         nbnd=nbnd+1
  709.        do j=1,nocomp(/2)
  710.         ire=ire+1
  711.         lisdua(ire)=nocomp(j)
  712.         noeled(ire)=iel
  713.         re(ire,1,1)= -vpocha(i,j)
  714.         call place(nomdu,lnomdu,ipo,nocomp(j))
  715.         if(ipo.eq.0) then
  716.         moterr(1:4)=nocomp(j)
  717.         moterr(5:10)='DUALE'
  718.         call erreur(1117)
  719.         endif
  720.         lisinc(ire)=nomdd(ipo)
  721.         if (lisinc(ire).eq.'LX  ') call erreur(1125)
  722.         noelep(ire)=iel
  723. *       write(6,*) ' 2 ire inc dua ',lisinc(ire),lisdua(ire)
  724.        enddo
  725.        enddo
  726.       enddo
  727. *  multiplicateur de lagrange
  728.       segact MCOORD*MOD
  729.       nbpts=nbpts+1
  730.       segadj mcoord
  731.       xl=0.d0
  732.       yl=0.d0
  733.       zl=0.d0
  734.       dl=0.d0
  735.       do i=2,num(/1)
  736.        ip=num(i,1)-1
  737.        xl=xl+xcoor(ip*(idim+1)+1)
  738.        if (idim.gt.1) yl=yl+xcoor(ip*(idim+1)+2)
  739.        if (idim.gt.2) zl=zl+xcoor(ip*(idim+1)+3)
  740.        dl=dl+xcoor(ip*(idim+1)+idim+1)
  741.       enddo
  742.       nbnnr=num(/1)-1
  743.       xl=xl/nbnnr
  744.       yl=yl/nbnnr
  745.       zl=zl/nbnnr
  746.       dl=dl/nbnnr
  747.       xcoor((nbpts-1)*(idim+1)+1)=xl
  748.       if (idim.gt.1) xcoor((nbpts-1)*(idim+1)+2)=yl
  749.       if (idim.gt.2) xcoor((nbpts-1)*(idim+1)+3)=zl
  750.       xcoor((nbpts-1)*(idim+1)+idim+1)=dl
  751.       lisinc(1)= 'LX'
  752.       lisdua(1)='FLX'
  753.       num(1,1)=nbpts
  754.       noelep(1)=1
  755.       noeled(1)=1
  756.       re(1,1,1)=0.D0
  757. *  un resultat
  758.       call ecrobj('RIGIDITE',mrigid)
  759.       SEGSUP MWGGM1
  760.       SEGSUP MWGGM2
  761.       SEGSUP MWGGM3
  762.       SEGSUP MWGGM4
  763.       SEGSUP MWGGM5
  764.       return
  765.       END
  766.  
  767.  
  768.  
  769.  
  770.  
  771.  
  772.  
  773.  
  774.  
  775.  
  776.  
  777.  
  778.  
  779.  
  780.  
  781.  
  782.  
  783.  
  784.  
  785.  
  786.  
  787.  
  788.  
  789.  
  790.  
  791.  
  792.  
  793.  

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