Télécharger la plaquette Cast3M


fr en

Télécharger supdep.eso

Retour à la liste

Numérotation des lignes : 

supdep

  1. C SUPDEP    SOURCE    PV090527  24/10/23    21:15:08     12044          
  2.       SUBROUTINE SUPDEP
  3. c=================================================================
  4. c cette procedure est appelée par SUPER
  5. c A partir du champ de deplacements interface et du champ d efforts
  6. c sur la sous structure, elle determine le champ de déplacement
  7. c sur toute la sous-structure
  8. c=================================================================
  9.       IMPLICIT INTEGER(I-N)
  10.       IMPLICIT REAL*8(A-H,O-Z)
  11. -INC SMSUPER
  12. -INC SMELEME
  13.  
  14. -INC PPARAM
  15. -INC CCOPTIO
  16. -INC SMMATRI
  17. -INC SMCHPOI
  18. -INC SMVECTD
  19. -INC SMRIGID
  20. -INC CCREEL
  21. c ????????????????????????????????????????????????????????????????
  22.       SEGMENT/BID/(BIDON(IIMAX+10)*D)
  23. c ????????????????????????????????????????????????????????????????
  24. c de la meme facon que dans MONDES , on va charger la maximum de
  25. c lignes en memoire pendant que l'on va gerer les efforts
  26. c IPLIG : première ligne en mémoire
  27. c IDLIG : dernière ligne en mémoire
  28. c----------------------------------------------------------------
  29.         PARAMETER (LPCRAY=10000)
  30.         segment ITRAA(NENS)
  31.         integer OOOVAL
  32.         logical leffort
  33.         character*4 mcle(1)
  34.         data mcle/'NOER'/
  35.  
  36.         dnorma=0.d0
  37.         dnormb=0.d0
  38. c creation d'une pile MRU
  39.         call OOOMRU(1)
  40.         call lirobj('SUPERELE',MSUPER,1,IRETOU)
  41.         if (ierr.ne.0) RETURN
  42. c lecture du déplacement interface
  43.         call lirobj('CHPOINT ',MDEPI,1,IRETOU)
  44.         if (ierr.ne.0) RETURN
  45. c lecture du  champ d'efforts si il existe
  46. c  27/11/2009  lecture obligatoire car risque d'erreur utilisateur
  47. c              due  a une mauvaise comprehension de l'operation
  48.         call lirobj('CHPOINT ',MFORC,1,IRETOU)
  49.         if (ierr.ne.0) RETURN
  50.         if (iretou .eq. 0) then
  51.             leffort = .false.
  52.         else
  53.             leffort = .true.
  54.         end if
  55.         noer=0
  56.         call lirmot(mcle,1,noer,0)
  57.  
  58.         segact MSUPER
  59. **      write(6,*) MRIGTO,MSUPEL,MSURAI,MBLOQU,MSUMAS,MCROUT
  60.         if (MCROUT.EQ.0) then
  61. *  pas de super element, on rend les LX du  champ d'entree
  62.           CALL EXCOPP(MDEPI,'LX  ',0,MDEPSO,'LX  ',0,1)
  63.           NAT=2
  64.           NSOUPO=0
  65.           call ecrobj('CHPOINT',MDEPSO)
  66.           return
  67.         endif
  68. c
  69. c recherche du nombre d'inconnues maitres
  70.         MRIGID = MSURAI
  71.         segact MRIGID
  72.         xMATRI = IRIGEL(4,1)
  73.         lagdua=msuper.islag
  74. *       write (6,*) ' msurai lagdua dans supdep ', mrigid,lagdua
  75.         segdes MRIGID
  76.         segact xMATRI
  77. *        XMATRI = IMATTT(1)
  78. *        segdes IMATRI
  79. *        segact XMATRI
  80.         nligra = RE(/1)
  81.         segdes XMATRI
  82. c
  83.         MCROUX = MCROUT
  84.         if(mcroux.eq.0) call erreur(1123)
  85.         if (ierr.ne.0) return
  86. * creation du meleme  des noeuds maitres
  87.        MRIGID=MRIGTO
  88.        SEGACT MRIGID
  89.  
  90.         meleme=msupel
  91.  
  92.         MRIGID = MSURAI
  93. c sauvegarde du deplacement interface
  94.         call reduir(MDEPI,MELEME,MDEPint)
  95.         call adchpo(mdepi,mdepint,mdext,1.d0,-1.d0)
  96. *       call ecchpo(mdepi,0)
  97. *       call ecchpo(mdepint,0)
  98. *        write (6,*) ' mdext '
  99. *       call ecchpo(mdext,0)
  100. c transformation du champ de déplacement en vecteur
  101. *       call ecmail(lagdua,0)
  102. **      call dbbch(mdepint,lagdua)
  103.         call chv1(MCROUX,MDEPInt,MVECTX,1)
  104. c normalisation
  105.         MMATRI = MCROUT
  106.         segact MMATRI
  107.         MDNOR = IDNORM
  108.         segact MDNOR
  109.         INC = DNOR(/1)
  110. c inbine : intger nombre inconnues esclaves
  111.  
  112.         inbine = INC - Nligra
  113. *       write (6,*) ' inbine inc nligra ',inbine,inc,nligra
  114.         MVECTD = MVECTX
  115.         segact MVECTD*mod
  116.         dnormb= 0.D0
  117.         do 1 ii1 = inbine+1,INC
  118.            VECTBB(ii1) = VECTBB(ii1) / DNOR(ii1)
  119. **         dnormb= max(dnormb,abs(VECTBB(ii1)))
  120.  1      continue
  121. **       dnormb = dnormb * xzprec*xzprec
  122.          dnormb = xpetit / xzprec
  123.  
  124.         MVECT1 = MVECTD
  125.         MILIGN = IILIGN
  126.         segact MILIGN
  127.         NNOE = ILIGN(/1)
  128.  
  129. c ????? idem que dans MONDES ??????
  130.         LTRK = MAX(1+LPCRAY,OOOVAL(1,4))
  131.         IIMAX = (((IJMAX+LTRK)/LTRK)+1)*LTRK
  132. c
  133. ***********************************************************************
  134. ***********************************************************************
  135. c pb à résoudre : L1t X +  L2t DEPI = D(-1) L1(-1) F
  136. c
  137. *  commencons par traiter le terme en effort
  138. c transformation du CHPOINT en vecteur
  139.         if (leffort) then
  140.           MCHPOI = MFORC
  141.           MRIGID = MRIGTO
  142.           segdes MRIGID
  143.           call copie2(mchpoi,mchpo1)
  144.           SEGACT MCHPO1
  145.         DO 432 I=1,mchpo1.IPCHP(/1)
  146.            MSOUPO=mchpo1.IPCHP(I)
  147.            SEGACT MSOUPO*MOD
  148.            IPT4=IGEOC
  149.            SEGINI,ipt5=ipt4
  150.            SEGDES ipt4
  151.            IGEOC=ipt5
  152.   432     CONTINUE
  153.           call dbbch(mchpo1,lagdua)
  154. *         call ecchpo(mchpo1 ,0)
  155.           call chv2(MCROUX,MCHPO1,MVECTX,1)
  156.           MVECTD = MVECTX
  157. c attention CHV2 desactive tout
  158.           segact MMATRI
  159.           segact MILIGN
  160.           segact MVECTD*MOD
  161.           MDIAG = IDIAG
  162.           segact MDIAG
  163. c
  164. c il faut normaliser ce vecteur
  165. c on va de plus recuperer l'indice du premier terme non nul
  166.           dnorma = 0.D0
  167.           inbi=vectbb(/1)
  168. *         write (6,*) ' inbine ',inbine
  169. *         write (6,*) ' vectbb-1 ',(vectbb(i),i=1,vectbb(/1))
  170. *         write (6,*) ' dnor     ',(dnor  (i),i=1,vectbb(/1))
  171.  
  172.  
  173.           do 2 ii1=1,inbine
  174.               VECTBB(ii1) = VECTBB(ii1) * DNOR(ii1)
  175. **            dnorma = max(dnorma,abs(VECTBB(ii1)))
  176.  2        continue
  177. *         write (6,*) ' vectbb-2 ',(vectbb(i),i=1,vectbb(/1))
  178. **        dnorma = dnorma * xzprec*xzprec
  179.           dnorma = xpetit/xzprec
  180. c
  181. c recherche du premier terme non nul
  182.           do 3 ii1=1,inbine
  183.              if (abs(VECTBB(ii1)).GE.dnorma) GOTO 4
  184.  3        continue
  185. *          write(6,*) ' attention vecteur force nulle'
  186.  4        continue
  187.           ipremf = ii1
  188.           if(dnorma.eq.0.d0) ipremf=max(1,inbine)
  189.  
  190. * effectuons maintenant le produit (L1)-1 EFFORT
  191. *-------------------------------------------------
  192. c ????????????
  193.           segini BID
  194.           ivalma = 0
  195. c ????????????
  196. c d'ou le bloc associé
  197.           INPREM = IPNO(ipremf)
  198. c et le bloc de la dernière ligne
  199.           INDERn=-1
  200.           if(inbine.ne.0)INDERN = IPNO(inbine)
  201.  
  202. c
  203.           IPLIG = INPREM
  204.           do 10 ii0=INPREM,INDERN
  205.              LIGN = ILIGN(ii0)
  206.              segact /ERR=20/ LIGN
  207.              IPRELL = IPREL
  208.              IVALMA = IVALMA+VAL(/1)
  209. *pvpv            if (IVALMA.GT.NGMAXY) GOTO 20
  210.              NA = IMMM(/1)
  211.              call supde2(IPPVV(1),IVPO(1),VAL(1),VECTBB(1),
  212.      >         na,inumli,inbine,ipremf,iprel,dnorma)
  213.  10       continue
  214.           segsup BID
  215.           IDLIG = INDERN
  216.           GOTO 60
  217. c on a eu des problèmes de memoire
  218.  20       continue
  219.           SEGSUP BID
  220.           IDLIG = ii0 - 1
  221.           itemp1 = ii0
  222. c c'est reparti
  223.           do 30 ii0=itemp1,INDERN
  224.             LIGN = ILIGN(ii0)
  225.             segact /ERR=22/ LIGN
  226.             NA = IMMM(/1)
  227.             IPRELL = IPREL
  228.             call supde2(IPPVV(1),IVPO(1),VAL(1),VECTBB(1),
  229.      >         na,inumli,inbine,ipremf,iprel,dnorma)
  230.             segdes LIGN*(NOMOD,MRU)
  231.             GOTO 30
  232. c il y a encore des problèmes
  233.  22         CONTINUE
  234. c il va falloir faire de la place ==> desactivation d'une ligne
  235. c reste-t-il des lignes à desactiver ?
  236.             if (IDLIG .lt. IPLIG) then
  237.               CALL ERREUR(5)
  238.             else
  239.               LIGN=ILIGN(IDLIG)
  240.               SEGDES LIGN*(NOMOD,MRU)
  241.               IDLIG=IDLIG-1
  242.             end if
  243.  30       continue
  244. c
  245. c muliplions maintenant par D-1
  246. c---------------------------------
  247.  60       continue
  248.           do 70 i=1,inbine
  249.              VECTBB(i) = VECTBB(i) * DIAG(i)
  250.  70       continue
  251.           segdes MDIAG
  252.         else
  253. c il n y a pas d'efforts
  254.            segini MVECTD
  255. c il faut preciser qu il n y a aucune ligne chargée
  256.            IDLIG = 0
  257.            IPLIG = NNOE
  258.         end if
  259. c
  260. c
  261. ***********************************************************************
  262. ***********************************************************************
  263. c Il faut maintenant effectuer le calcul L1t X1 + L2t X2 = Nouveau F
  264. c il ne s'agit que d'une remontée
  265. c
  266. c on va en profiter pour compter les mouvements d'ensemble
  267.         iaa = 0
  268.         segini itraa
  269. c
  270. c inumli : numero de la ligne en cours
  271.         inumli = INC
  272.         do 120 ii0=NNOE,1,-1
  273.            LIGN = ILIGN(ii0)
  274.  122       continue
  275.            if ((ii0.gt.IDLIG).or.(ii0.lt.IPLIG)) then
  276.              segact /ERR=124/ LIGN
  277.            else
  278.              IDLIG = IDLIG - 1
  279.            end if
  280.            NA = IMMM(/1)
  281.          IFIB=IVPO(/1)
  282.            call supde1(IPPVV(1),IVPO(1),VAL(1),VECTBB(1),
  283.      &      MVECT1.VECTBB(1),na,inumli,inbine,iprel,ifib,dnormb)
  284.            do ibb=1,NA
  285.             if (IMMM(ibb) .ne. 0) then
  286.              iaa = iaa+1
  287.              itraa(iaa)=iprel+ibb-1
  288.             end if
  289.            end do
  290.            segdes LIGN*(NOMOD,MRU)
  291.            GOTO 120
  292.  124       CONTINUE
  293. c encore des problèmes mémoires
  294.            if  (IDLIG .lt. IPLIG) then
  295.              CALL ERREUR(5)
  296.            else
  297.              LIG1=ILIGN(IDLIG)
  298.              SEGDES LIG1*(NOMOD,MRU)
  299.              IDLIG=IDLIG-1
  300.                  GOTO 122
  301.            end if
  302.  120    continue
  303.  
  304.         segsup MVECTD
  305.         MVECTD = MVECT1
  306. c
  307. *******************************************************************************
  308. *******************************************************************************
  309. c gestion des déplacements de corps rigide
  310. c meme chose que dans MONDES
  311. c                     ------
  312. c
  313.       if (nens .ne. 0) then
  314.         MINCPO = IINCPO
  315.         MIMIK = IIMIK
  316.         segact MINCPO,MIMIK
  317.         ipt2 = IGEOMA
  318.         segact ipt2
  319.         NNOE = INCPO(/2)
  320.         IINC1=INCPO(/1)
  321. c
  322.         XMA = xpetit
  323.         XMAL = xpetit
  324.         inan=0
  325.         do kk=1,INC
  326.           if (NOER.EQ.0.OR.(noer.eq.1.and.abs(vectbb(kk)).lt.xgrand))
  327.      &           goto 500
  328.           inan = inan + 1
  329.           vectbb(kk)=0.D0
  330.  500   continue
  331.          if (ittr(kk).eq.0) then
  332.            XMA=MAX(XMA,ABS(VECTBB(kk)))
  333.          else
  334.            XMAL=MAX(XMAL,ABS(VECTBB(kk)))
  335.          endif
  336.         end do
  337.         XMA = XMA * 1.d-10
  338.         XMAL = XMAL * 1.d-10
  339.         xmal = max(xma*1d-2,xmal)
  340. *       write (6,*) ' supdep cma cmal ',xma,xmal
  341. c boucle sur les mouvements d'ensemble
  342.         do 200 ia=1,NENS
  343.           i1=itraa(ia)
  344.           j=IPNO(i1)
  345.           do 210 k=1,iinc1
  346.             if(INCPO(K,J).eq.i1) goto 220
  347.  210      continue
  348.           call erreur(5)
  349.           return
  350.  220      continue
  351.           if ((ittr(i1).eq.0).and.(ABS(VECTBB(i1)).le.XMA)) GOTO 250
  352.           if ((ittr(i1).ne.0).and.(ABS(VECTBB(i1)).le.XMAL)) GOTO 250
  353.           MOTERR(1:4)=IMIK(k)
  354.           INTERR(1)=ipt2.NUM(1,J)
  355.           write (6,*) ' vectbb) ',vectbb(i1)
  356. **        CALL ERREUR(149)
  357. **        RETURN
  358.           call soucis(149)
  359.  251      continue
  360.  250      continue
  361.           jjk = ipt2.NUM(1,J)
  362.           IF(ITTR(I1).EQ.0) WRITE(IOIMP,280) JJK,IMIK(K)
  363.  
  364.           IF (IIMPI.NE.0 .AND. ITTR(I1).NE.0)
  365.      &                WRITE(IOIMP,290) JJK,IMIK(K)
  366.  
  367.  200    continue
  368.   280 FORMAT(' INDETERMINATION DETECTEE AU NOEUD ',I6,' INCONNUE ',
  369.      *   A4,/,' INDETERMINATION LEVEE PAR LA MISE A ZERO DE ',
  370.      * 'LA SUSDITE')
  371.   290 FORMAT(' INDETERMINATION ENTRE CONDITIONS IMPOSEES DETECTEE ',
  372.      * 'AU NOEUD ',I6,' INCONNUE ',A4,/,' INDETERMINATION LEVEE ',
  373.      * 'PAR LA SUPPRESSION DE LA CONDITION REDONDANTE ')
  374.  
  375.         segdes MINCPO,MIMIK
  376.           segdes ipt2
  377.       end if
  378. c
  379. *******************************************************************************
  380. *******************************************************************************
  381. * un petit coup de normalisation
  382.         do 300 ii1=1,INC
  383.             VECTBB(ii1) = VECTBB(ii1) * DNOR(ii1)
  384.  300    continue
  385.  
  386.         segdes MDNOR
  387.         segdes MVECTD
  388.         MVECTX = MVECTD
  389.         MMATRX = MMATRI
  390.         MRIGTX = MRIGTO
  391. *        write(6,*) ' mrigto mmatri mvectx ', mrigto,mmatri,mvectx
  392. c il faut transformer ce vecteur en chpoint
  393.         call VCH1(MMATRX,MVECTX,ISOLU,MRIGTX)
  394.         if (lagdua.ne.0) call dbbcf(isolu,lagdua)
  395. **      write (6,*) ' mdext '
  396. **      call ecchpo(mdext,0)
  397. **      write (6,*) ' isolu '
  398. **      call ecchpo(isolu,0)
  399.         call adchpo(isolu,mdext,mchpoi,1.d0,1.d0)
  400. *       mchpoi=isolu
  401.         segact MCHPOI*mod
  402.         JATTRI(1)=1
  403.         segdes MCHPOI
  404.  
  405.         segdes MMATRI
  406. c il faut maintenant rajouter les déplacements imposés éliminés (jrcond)
  407.         MRIGID = MRIGTO
  408. *        segact MRIGID
  409.         segdes MRIGID
  410. c
  411.         segdes MSUPER
  412. c il faut ajuster le type du champ
  413.  
  414. c on rajoute les multiplicateurs interfaces
  415.         call ecrobj('CHPOINT ',mchpoi)
  416.  
  417. c désactivation de la pile MRU
  418.         call OOOMRU(0)
  419.         return
  420.         end
  421.  
  422.  
  423.  
  424.  
  425.  
  426.  
  427.  
  428.  
  429.  
  430.  
  431.  
  432.  
  433.  
  434.  
  435.  
  436.  
  437.  
  438.  
  439.  
  440.  
  441.  
  442.  

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