Télécharger la plaquette Cast3M


fr en

Télécharger precop.eso

Retour à la liste

Numérotation des lignes : 

precop

  1. C PRECOP    SOURCE    CB215821  24/04/12    21:16:54     11897          
  2.       SUBROUTINE PRECOP(IPMODL,IPCHA1,IPTAB,IPSTRS,IRAN,
  3.      &                  PS1,IPCHC1,IRET)
  4. C======================================================================C
  5. C                                                                      C
  6. C     ENTREES :                                                        C
  7. C                                                                      C
  8. C       IPMODL: POINTEUR SUR UN MMODEL                                 C
  9. C       IPCHA1: Pointeur sur le MCHAML de CARACTERISTIQUES             C
  10. C       PS1   : Force sous vérin
  11. C       IRAN  :  Pointeur sur maillage des points d application
  12. C       IPCHC1: POINTEUR SUR le MCHAML DE PRECONTRAINTE INITIALE
  13. C       IPTAB : pointeur sur la table des parametres frottement ....   C
  14. C                                                                      C
  15. C     SORTIES :                                                        C
  16. C                                                                      C
  17. C       IPSTRS:  MCHAML CONTENANT LES PRECONTAINTES ET LES FORCES      C
  18. C                DU CABLE SUR LE BETON                                 C
  19. C       IRET  :  1   cela c est bien passe
  20. C                0   probleme dans le traitement                       C
  21. C                                                                      C
  22. C======================================================================C
  23.       IMPLICIT INTEGER(I-N)
  24.       IMPLICIT REAL*8(A-H,O-Z)
  25. C
  26.  
  27. -INC PPARAM
  28. -INC CCOPTIO
  29. -INC SMCHAML
  30. -INC SMELEME
  31. -INC SMINTE
  32. -INC SMMODEL
  33. -INC CCHAMP
  34. -INC SMCOORD
  35. -INC SMTABLE
  36. C
  37.       REAL*8 XVALIN,XVALRE
  38.       LOGICAL LOGRE,LOGIN
  39.       CHARACTER*8 TAPIND,TYPOBJ,CHARIN,CHARRE
  40.       CHARACTER*4 MOPAR(6)
  41. *
  42.       DATA MOPAR/'F1  ','F2  ','GANC','RMU0','FPRG','RH10'/
  43.       SEGMENT INFO
  44.          INTEGER INFELL(JG)
  45.       ENDSEGMENT
  46. * stockage  de pointeurs sur des segment sielc  crees dans splitag
  47.       segment siezo
  48.       integer iezon(*)
  49.       endsegment
  50.       segment sielc
  51.       integer ideb(*),ifin(*),nbcz,isens(2,*),idejvu(*)
  52.       endsegment
  53. C
  54.       SEGMENT ALTRAV
  55.          REAL*8 ANG(NAM+1),DANG(NAM),ACUR(NAM+1),DACUR(NAM)
  56.       ENDSEGMENT
  57. C
  58.       SEGMENT NOTYPE
  59.         CHARACTER*16 TYPE(NBTYPE)
  60.       ENDSEGMENT
  61. C
  62.       SEGMENT MPTVAL
  63.          INTEGER IPOS(NS)    ,NSOF(NS)
  64.          INTEGER      IVAL(NCOSOU)
  65.          CHARACTER*16 TYVAL(NCOSOU)
  66.       ENDSEGMENT
  67. C
  68.       SEGMENT WRK3
  69.          REAL*8 X1(3),X2(3),X3(3),RL(3),RS(3)
  70.       ENDSEGMENT
  71.  
  72. C
  73.       CHARACTER*8 CMATE,LNOM
  74.       CHARACTER*(NCONCH) CONM
  75.       PARAMETER ( NINF=3 )
  76.       PARAMETER ( ITEMAX=50 )
  77.       INTEGER INFOS(NINF)
  78.       DIMENSION FGG(ITEMAX),XLAM(ITEMAX)
  79.       logical ivers,lsupco
  80.  
  81.       ivers=.false.
  82.  
  83. cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
  84.  
  85.       if (ipchc1.ne.0) then
  86.       MCHELM  = ipchc1
  87.       segact MCHELM
  88. C   actuellement  un cable n a qu une zone de maillage
  89.       imail = imache(1)
  90.       endif
  91.  
  92.       IRET=0
  93. C     verification on doit avoir des elements soit BARR soit  CERC
  94.       MMODEL = ipmodl
  95.       segact mmodel
  96.       do i = 1,kmodel(/1)
  97.       imodel = kmodel(i)
  98.       segact imodel
  99.       if(ipchc1.ne.0.and.imail.ne.imamod) then
  100.       call erreur(21)
  101.       return
  102.       endif
  103.       meleme = imamod
  104.       segact meleme
  105.       if(itypel.ne.2) then
  106.         call erreur(16)
  107.         return
  108.       endif
  109.       if(ifomod.eq.0.and.itypel.ne.1) then
  110.       call erreur(21)
  111.       return
  112.       endif
  113.       enddo
  114.  
  115. C --- Vérification du support des MCHAML --------------
  116.  
  117.       if(ipchc1.gt.0) then
  118.          call quesup(ipmodl,ipchc1,5,0,iresig,iret2)
  119.          if(iresig.eq.9999) then
  120.             moterr(1:8)='CONTRAIN'
  121.             call erreur(124)
  122.             return
  123.          endif
  124.        endif
  125.  
  126.       call quesup(ipmodl,ipcha1,5,0,irecar,iret2)
  127.       if(irecar.eq.9999) then
  128.          moterr(1:8)='CARACTER'
  129.          call erreur(124)
  130.          return
  131.       endif
  132.       if(irecar.eq.1) then
  133.          ipcara = 0
  134.          iret = 0
  135.          call chasup(ipmodl,ipcha1,ipcara,iret,5)
  136.          if(iret.ne.0) return
  137.       else
  138.          ipcara = ipcha1
  139.       endif
  140. *---------------------------------------------------------
  141.       ipmod0 = ipmodl
  142. *---------------------------------------------------------
  143.       if(ifomod.ne.0.and.ifomod.ne.1) then
  144.       ipcar0 = ipcara
  145.       ipchc0 = ipchc1
  146.       call splitag(ipmodl,ipcara,ipchc1,ipmod2,ipcar2,
  147.      &      ipchc2,siezo)
  148.       ipmodl = ipmod2
  149.       ipcara = ipcar2
  150.       ipchc1 = ipchc2
  151.       else
  152. C  pas elegant mais on duplique  le modele
  153.       MMODEL = ipmodl
  154.       segini , mmode2=MMODEL
  155.       ipmodl = mmode2
  156.       endif
  157. C
  158.       ipmodt = ipmodl
  159.       MMODEL = IPMODL
  160. C
  161.       SEGACT MMODEL
  162.       NSOUS=KMODEL(/1)
  163. C
  164. C ----PREPARATION DU CHAMP DE SORTIE----------------
  165. C
  166.       mmode2 = ipmod0
  167.       segact mmode2
  168.       N1=mmode2.kmodel(/1)
  169.       L1=11
  170.       N3=6
  171.       segini mchel5
  172.       mchel5.titche='CONTRAINTES'
  173.       mchel5.ifoche=ifour
  174.  
  175.       do isous=1,n1
  176.       imodel=mmode2.kmodel(isous)
  177.       segact imodel
  178.       mchel5.imache(isous)=imamod
  179.       mchel5.conche(isous)=conmod
  180.       mele= nefmod
  181. C====
  182.       CALL ELQUOI(MELE,0,5,INFO,IMODEL)
  183.          IF (IERR.NE.0) THEN
  184.             SEGSUP MCHEL5
  185.             IRET=0
  186.             RETURN
  187.       ENDIF
  188.        MINTE= INFELL(11)
  189. C====
  190.          MCHEL5.INFCHE(ISOUS,1)=0
  191.          MCHEL5.INFCHE(ISOUS,2)=0
  192.          MCHEL5.INFCHE(ISOUS,3)=NIFOUR
  193.          MCHEL5.INFCHE(ISOUS,4)=MINTE
  194.          MCHEL5.INFCHE(ISOUS,5)=0
  195.          MCHEL5.INFCHE(ISOUS,6)=5
  196. C
  197.       N2=1
  198.       segini MCHAML
  199.       mchel5.ichaml(isous)=mchaml
  200.       nomche(1)='EFFX'
  201.       typche(1)='REAL*8'
  202.       meleme = imamod
  203.       segact meleme
  204.       N1PTEL=2
  205.       N1EL= NUM(/2)
  206.       N2PTEL=0
  207.       N2EL=0
  208.       segini melval
  209.       ielval(1)=melval
  210.       segsup info
  211.       enddo
  212.  
  213. C
  214. C ----FIN    PREPARATION DU CHAMP DE SORTIE----------------
  215. C
  216.  
  217. C   preparation de la recherche du sens d application des tensions
  218.        if(iran.ne.0)  then
  219.        ipt3=iran
  220.          segact ipt3
  221.          segini , ipt5=ipt3
  222.          ity=ipt5.itypel
  223.            if(ity.ne.1) then
  224.           call change(iran,1)
  225.           ipt5=iran
  226.           segact ipt5
  227.            endif
  228.        else
  229.        ipt5=0
  230.        endif
  231.  
  232. C
  233. C --- CREATION DU MCHELM  DE CALCUL  ( TEMPORAIRE )
  234. C
  235.       N1=NSOUS
  236. C
  237.       L1=11
  238.       N3=6
  239.       SEGINI MCHELM
  240.       TITCHE='CONTRAINTES'
  241.       IFOCHE=IFOUR
  242.       IPSTRS=MCHELM
  243. C
  244. C --- DEBUT DE LA BOUCLE SUR LES DIFFERENTES ZONES
  245. C
  246.       call oooeta(mcoord,ieta,imod)
  247.        if (ieta.ne.1) then
  248.            ieta=1
  249.            segact mcoord
  250.        endif
  251.       DO 500 ISOUS=1,NSOUS
  252. C
  253. C --- INITIALISATION
  254. C
  255.          MOCARA=0
  256.          MOMATR=0
  257.          MOSTRS=0
  258.          IVAMAT=0
  259.          IVACAR=0
  260.          IVASTR=0
  261.          ivasi0=0
  262. c
  263.          ipeffx = 0
  264.          ipyoun = 0
  265.          ipff   = 0
  266.          ipphif = 0
  267.          ipganc = 0
  268.          iprmu0 = 0
  269.          ipfprg = 0
  270.          iprh10 = 0
  271.          ipsect = 0
  272.          IMODEL=KMODEL(ISOUS)
  273.          SEGACT IMODEL
  274.          IPMAIL=IMAMOD
  275.          CONM  =CONMOD
  276.          IMACHE(ISOUS)=IPMAIL
  277.          CONCHE(ISOUS)=CONMOD
  278. C ...... TRAITEMENT DU MODELE ......................
  279.  
  280.          MELE=NEFMOD
  281.          MELEME=IPMAIL
  282. C
  283.          segact meleme
  284.          itpdel=itypel
  285.  
  286.          NFOR=FORMOD(/2)
  287.          NMAT=MATMOD(/2)
  288. C
  289. C    ... on ordonne le nouveau maillage (si on en a besoin) ...
  290.  
  291.          NBELEM = NUM(/2)
  292. C=====================================================
  293.          if(itypel.eq.1) then
  294.  
  295. c    ... dans le cas des POI1 on n'a rien a orienter ...
  296. C=====================================================
  297.          elseif(itypel.eq.2) then
  298.  
  299.            ipt2 = meleme
  300.             NOD1=NUM(1,1)
  301.             NOD2=NUM(2,num(/2))
  302. C           write(6,*)' noeuds extremites ' ,nod1,nod2
  303.          if(iran.eq.0)  then
  304.          ivers=.false.
  305.          else
  306.             do ino=1,ipt5.num(/2)
  307.             inodi=ipt5.num(1,ino)
  308.             if(inodi.eq.nod1) then
  309.             ivers=.false.
  310.             goto  364
  311.             elseif(inodi.eq.nod2) then
  312.             ivers=.true.
  313.             goto  364
  314.             endif
  315.             enddo
  316. C     aucune des extremites n'est  dans la liste
  317. C           on applique au point final
  318.             call erreur(833)
  319. C  attention aux desactivations
  320.             SEGSUP MCHELM
  321.             IRET=0
  322.             RETURN
  323. 364      continue
  324.          endif
  325. C=====================================================
  326. c    ... si ITYPEL différent de 1 et 2 ...
  327.          else
  328.  
  329. c       ... Type d''élément support inconnu !!!!!!!!!!! ............
  330.             SEGSUP MCHELM
  331.             IRET=0
  332.             RETURN
  333.             call erreur(19)
  334.             return
  335.          endif
  336. C=====================================================
  337. C
  338. C ------ NATURE DU MATERIAU
  339. C
  340.          CALL NOMATE(FORMOD,NFOR,MATMOD,NMAT,CMATE,MATE,INAT)
  341.          IF (CMATE.EQ.' ') THEN
  342.           CALL ERREUR(251)
  343.             SEGSUP MCHELM
  344.             IRET=0
  345.             RETURN
  346.           ENDIF
  347. C
  348. C ------ INFORMATION SUR L'ELEMENT FINI
  349. C
  350.          CALL ELQUOI(MELE,0,5,IPINF,IMODEL)
  351.          IF (IERR.NE.0) THEN
  352.             SEGSUP MCHELM
  353.             IRET=0
  354.             RETURN
  355.          ENDIF
  356. C
  357.          INFO=IPINF
  358.          MELE =INFELL(1)
  359.          MFR  =INFELL(13)
  360.          MINTE=INFELL(11)
  361.          NSTRS=INFELL(16)
  362.  
  363.          SEGINI WRK3
  364. C
  365. C ------ CREATION DU TABLEAU INFOS
  366. C
  367.          IPCHE1=0
  368.          IPCHE2=0
  369. c ... ATTENTION !!! A quoi ca sert ????????????? .......
  370.         CALL IDENT(IPMAIL,CONM,IPCHE1,IPCHE2,INFOS,IRTD)
  371.          IF (IRTD.EQ.0) THEN
  372.             SEGSUP MCHELM
  373.             SEGSUP INFO
  374.              IRET=0
  375.             RETURN
  376.          ENDIF
  377.  
  378. c    ... Remplissage du tableau INFCHE du MCHELM ...
  379.  
  380.          INFCHE(ISOUS,1)=0
  381.          INFCHE(ISOUS,2)=0
  382.          INFCHE(ISOUS,3)=NIFOUR
  383.          INFCHE(ISOUS,4)=MINTE
  384.          INFCHE(ISOUS,5)=0
  385.          INFCHE(ISOUS,6)=5
  386. C
  387. C ------ LECTURE de NBPGAU dans
  388. C
  389.          SEGACT MINTE
  390.          NBPGAU=POIGAU(/1)
  391. C
  392. C ------ RECHERCHE DES NOMS DE COMPOSANTES des contraintes ...
  393. C
  394.          if(lnomid(4).ne.0) then
  395.           nomid=lnomid(4)
  396.           segact nomid
  397.           mostrs=nomid
  398.           nstr=lesobl(/2)
  399.           nfac=lesfac(/2)
  400.           lsupco=.false.
  401.         else
  402.          lsupco=.true.
  403.          CALL IDCONT(IMODEL,IFOUR,MOSTRS,NSTR,NFAC)
  404.         endif
  405.          nomid = mostrs
  406.          segact,nomid
  407.          call place(lesobl,nstr,ipeffx,'EFFX')
  408. C
  409. C ------ VERIFICATION DE LEUR PRESENCE et ...
  410. C ------ RECUPERATION DES CONTRAINTES INITIALES (s'il le faut) ...
  411. C
  412.          if(IPCHC1.gt.0) then
  413.  
  414.             NBTYPE=INFELL(16)
  415.             SEGINI NOTYPE
  416.             MOTYPE=NOTYPE
  417.             TYPE(1)='REAL*8'
  418.  
  419.             call komcha(ipchc1,ipmail,conm,mostrs,motype,1,infos,3,
  420.      &                  ivasi0)
  421.  
  422.             IF (IERR.NE.0) THEN
  423.                NSTRS=0
  424.                GOTO 9990
  425.             ENDIF
  426.  
  427.             SEGSUP NOTYPE
  428.          endif
  429.  
  430. C
  431. C ------ CREATION DU MCHAML DE LA SOUS ZONE
  432. C
  433.          N2=NSTRS
  434.          SEGINI MCHAML
  435.          ICHAML(ISOUS)=MCHAML
  436.  
  437.          NS=1
  438.          NCOSOU=NSTRS
  439.          SEGINI MPTVAL
  440.          IVASTR=MPTVAL
  441.  
  442. C
  443. C ------ RECHERCHE DE LA TAILLE DES MELVAL A ALLOUER
  444. C
  445.          N1PTEL=nbpgau
  446.          N1EL=NBELEM
  447.          N2PTEL=0
  448.          N2EL=0
  449. C
  450.          NOMID=MOSTRS
  451.          SEGACT NOMID
  452. C    initialisation   des melvals  des precontaintes
  453.          DO 100 ICOMP=1,NSTRS
  454.             NOMCHE(ICOMP)=LESOBL(ICOMP)
  455.             TYPCHE(ICOMP)='REAL*8'
  456.             SEGINI MELVAL
  457.             IELVAL(ICOMP)=MELVAL
  458.             IVAL(ICOMP)=MELVAL
  459.  100     CONTINUE
  460. C            ipacara  preparé dans splitag
  461. C  recuperation du module de young  dans ivamat ( position 1 )
  462. C  recuperation de la section       dans ivacar ( position 2 )
  463.          mchel6=ipcara
  464.          segact mchel6
  465.          mcham6=mchel6.ichaml(isous)
  466.          segact mcham6
  467.          melval= mcham6.ielval(1)
  468.          segact melval
  469.          ivamat=melval
  470.          melval= mcham6.ielval(2)
  471.          segact melval
  472.          ivacar=melval
  473.  
  474. C
  475. C ------ CALCUL lui même --------------------------
  476. C
  477. C ... cas du CERC ...
  478. C
  479.  
  480.       IF(MELE.EQ.95) THEN
  481.          DO 2004 IB=1,NBELEM
  482.             MPTVAL=IVASTR
  483.             DO 1701 I=1,NSTRS
  484.                MELVAL=IVAL(I)
  485.                VELCHE(1,IB)=PS1
  486.  1701       CONTINUE
  487.          iret=1
  488.  2004    CONTINUE
  489.  
  490.        ELSE
  491. C  ----------------valeurs par defaut ----------------------
  492.  
  493.       f1 = 0.18d0
  494.       f2 = 0.002d0
  495.       ganc = 0.d0
  496.       rmu0 = 0.43d0
  497.       fprg = 1.7d9
  498.       rh10 = 2.5d0
  499.  
  500.       if(iptab.ne.0) then
  501.       MTABLE = iptab
  502.       segact mtable
  503.       nbind= mlotab
  504.           IVALIN = 0
  505.          XVALIN = 0.D0
  506.          LOGIN  = .TRUE.
  507.          IOBIN  = 0
  508.          TAPIND = 'MOT'
  509. C
  510.       TYPOBJ = '        '
  511.       Call acctab(MTABLE,'MOT',IVALIN,XVALIN,'FF',LOGIN,0,
  512.      &             TYPOBJ,IVELRE,XVALRE,CHARRE,LOGRE,IOBRE)
  513.       IF(TYPOBJ.EQ.'FLOTTANT') F1 = XVALRE
  514.       IF(TYPOBJ.EQ.'ENTIER  ') F1 = IVELRE
  515.       TYPOBJ = '        '
  516.       Call acctab(MTABLE,'MOT',IVALIN,XVALIN,'PHIF',LOGIN,0,
  517.      &             TYPOBJ,IVELRE,XVALRE,CHARRE,LOGRE,IOBRE)
  518.       IF(TYPOBJ.EQ.'FLOTTANT') F2 = XVALRE
  519.       IF(TYPOBJ.EQ.'ENTIER  ') F2 = IVELRE
  520.       TYPOBJ = '        '
  521.       Call acctab(MTABLE,'MOT',IVALIN,XVALIN,'GANC',LOGIN,0,
  522.      &             TYPOBJ,IVELRE,XVALRE,CHARRE,LOGRE,IOBRE)
  523.       IF(TYPOBJ.EQ.'FLOTTANT') GANC = XVALRE
  524.       IF(TYPOBJ.EQ.'ENTIER  ') GANC = IVELRE
  525.       TYPOBJ = '        '
  526.       Call acctab(MTABLE,'MOT',IVALIN,XVALIN,'RMU0',LOGIN,0,
  527.      &             TYPOBJ,IVELRE,XVALRE,CHARRE,LOGRE,IOBRE)
  528.       IF(TYPOBJ.EQ.'FLOTTANT') RMU0 = XVALRE
  529.       IF(TYPOBJ.EQ.'ENTIER  ') RMU0 = IVELRE
  530.       TYPOBJ = '        '
  531.       Call acctab(MTABLE,'MOT',IVALIN,XVALIN,'FPRG',LOGIN,0,
  532.      &             TYPOBJ,IVELRE,XVALRE,CHARRE,LOGRE,IOBRE)
  533.       IF(TYPOBJ.EQ.'FLOTTANT') FPRG = XVALRE
  534.       IF(TYPOBJ.EQ.'ENTIER  ') RPRG = IVELRE
  535.       TYPOBJ = '        '
  536.       Call acctab(MTABLE,'MOT',IVALIN,XVALIN,'RH10',LOGIN,0,
  537.      &             TYPOBJ,IVELRE,XVALRE,CHARRE,LOGRE,IOBRE)
  538.       IF(TYPOBJ.EQ.'FLOTTANT') RH10 = XVALRE
  539.       IF(TYPOBJ.EQ.'ENTIER  ') RH10 = IVELRE
  540.       segdes MTABLE
  541.  
  542. 4321  format(6E12.5)
  543.        endif
  544. C----------------------------------------------------------
  545. C       write(6,4321) f1,f2,ganc,rmu0,fprg,rh10
  546. C
  547. C ... Boucle sur les éléments dans lesquels on fait le calcul ...
  548. C
  549.  
  550.          SLON=0.0D0
  551.          FAI=0.0D0
  552.          ids= idim+1
  553.             segact,minte
  554.             if(nbpgau.ne.2) then
  555.                call erreur(5)
  556.                return
  557.             endif
  558. *
  559. *    NEW   1-ERE BOUCLE SUR LES ELEMENTS POUR TROUVER
  560. *     LA LONGUEUR D'INFLUENCE DU RECUL D'ANCRAGE
  561. *
  562.             NAM =NBELEM
  563.             SEGINI ALTRAV
  564.             IAM=0
  565.             SUM =0.D0
  566. *****       PRINT *,'IVERS=',IVERS
  567. *
  568.             DO 5005 IC=1,NBELEM
  569.  
  570.             if(ivers) then
  571. C    ordre inverse
  572.                JC= NBELEM+1-IC
  573.                NC3=ipt2.NUM(1,JC)
  574.                NC2=ipt2.NUM(2,JC)
  575.               IF(JC.EQ.NBELEM) THEN
  576.                  NC1=NC2
  577.               ELSE
  578.                  NC1=IPT2.NUM(2,JC+1)
  579.               ENDIF
  580.             else
  581. C    ordre normal
  582.                JC=IC
  583.                NC3=ipt2.NUM(2,JC)
  584.                NC2=ipt2.NUM(1,JC)
  585.                IF(JC.EQ.1) THEN
  586.                   NC1=NC2
  587.                ELSE
  588.                   NC1=IPT2.NUM(1,JC-1)
  589.                ENDIF
  590.             endif
  591. *****         PRINT *,'NC1=',NC1,' NC2=',NC2,' NC3=',NC3
  592.               JS1=(NC1-1)*IDS
  593.               JS2=(NC2-1)*IDS
  594.               JS3=(NC3-1)*IDS
  595.  
  596.                DO     IW=1,IDIM
  597.                   X1(IW)=XCOOR(JS1+IW)
  598.                   X2(IW)=XCOOR(JS2+IW)
  599.                   X3(IW)=XCOOR(JS3+IW)
  600.                enddo
  601. C          --- Distance entre points
  602.                DS1=0.0D0
  603.                DS2=0.0D0
  604.                DO     IW=1,IDIM
  605.                   RL(IW)=X3(IW)-X2(IW)
  606.                   RS(IW)=X2(IW)-X1(IW)
  607.                   DS1=DS1+RL(IW)*RL(IW)
  608.                   DS2=DS2+RS(IW)*RS(IW)
  609.                enddo
  610.                DS1=SQRT(DS1)
  611.                DS2=SQRT(DS2)
  612.                CS1=0.0D0
  613.                IF(DS2.NE.0.D0) THEN
  614.                   DO     IW=1,IDIM
  615.                      RL(IW)=RL(IW)/DS1
  616.                      RS(IW)=RS(IW)/DS2
  617.                      CS1=CS1+RL(IW)*RS(IW)
  618.                   enddo
  619.                ELSE
  620.                      CS1=1.D0
  621.                ENDIF
  622.                IF(CS1.GT.1.0) CS1=1.0D0
  623.                ALFA=ACOS(CS1)
  624.                IAM=IAM+1
  625.                ANG(IAM)=FAI
  626.                DANG(IAM)=ALFA
  627.                ACUR(IAM)=SLON
  628.                DACUR(IAM)=DS1
  629.                SLON=SLON+DS1
  630.                FAI=FAI+ALFA
  631.  
  632. *          PRINT *,'IAM=',IAM, ' ACUR =',ACUR(IAM),
  633. *     &    ' ANG=',ANG(IAM),' DACUR=',DACUR(IAM),' DANG=',DANG(IAM)
  634. *
  635.                IF(SLON.NE.0.D0) THEN
  636.                 SUM = SUM + FAI/SLON
  637.                ENDIF
  638. *
  639.  5005       CONTINUE
  640.             ANG(NAM+1)=FAI
  641.             ACUR(NAM+1)=SLON
  642. *
  643.             SLONT=SLON
  644.             FAIT=FAI
  645.             XLMBDA=0.D0
  646. *
  647.             IF(GANC.EQ.0.D0) GO TO 7999
  648. *
  649. *------------------------------------------------------------
  650. *    CALCUL DE LA LONGUEUR D'INFLUENCE DU RECUL D'ANCRAGE
  651. *              LE CAS ECHEANT  ( GANC > 0 )
  652. *------------------------------------------------------------
  653. *
  654. *      ON PREND LES MODULE ET SECTION DU 1-ER ELEMENT
  655. *      POUR LA 1-ERE APPROXIMATION DE XLAMBDA
  656. *
  657.             MELVAL =IVAMAT
  658.             EA = VELCHE(1,1)
  659.             MELVAL =IVACAR
  660.             SECT1 = VELCHE(1,1)
  661. *        write(6,*) ' f1, f2, slont, ps1, sect1, ea,ganc,sum,nbelem,fai'
  662. *       write(6,*)  f1, f2, slont, ps1, sect1, ea,ganc,sum,nbelem,fai
  663.             IF(F1.EQ.0.D0.OR.FAI.LE.1.D-4) THEN
  664.               IF(F2.EQ.0.D0) THEN
  665.                 XLMBDA=SLONT
  666.               ELSE
  667.                 aaa= F2*GANC*EA*SECT1/PS1
  668.                  bbb = aaa ** 0.5
  669.                  ccc = 1.-bbb
  670.                  if( ccc .le.0.) then
  671.                      call erreur(768)
  672.                      return
  673.                  else
  674. *               write(6,*) ' aaa,bbb,ccc',aaa,bbb,ccc
  675.                     XLMBDA = -(LOG(1.D0-
  676.      &                     (F2*GANC*EA*SECT1/PS1)**0.5))/F2
  677.                  endif
  678.               ENDIF
  679.             ELSE
  680.               ALSURL=SUM/FLOAT(NBELEM)
  681.               PP= F1*ALSURL + F2
  682.               IF(PP.EQ.0.D0) THEN
  683.                 XLMBDA=SLONT
  684.               ELSE
  685.                 XLMBDA = SQRT( (GANC*EA*SECT1)/(PS1*PP) )
  686.               ENDIF
  687. *          PRINT *,'ALSURL=',ALSURL
  688. *          PRINT *,'F2=',F2, 'PP=',PP
  689.             ENDIF
  690. *
  691. *    TEST SI XLAMBBDA > SLON
  692. *
  693.             IF(XLMBDA.GT.SLONT) XLMBDA=SLONT
  694.  
  695. *       PRINT *,'INITIALISATION DE XLMBDA=',XLMBDA
  696. *
  697.             IF(XLMBDA.LT.0.D0) THEN
  698. *****           PRINT *,' ATTENTION XLAMBDA INITIAL EST < 0 '
  699.                 CALL ERREUR(460)
  700.                 SEGSUP ALTRAV
  701.                 GO TO 9990
  702.             ENDIF
  703. *
  704.             ITER=0
  705.             LAST=0
  706.             IDICH=0
  707. *
  708. *         ITERATIONS
  709. *
  710.  8000       CONTINUE
  711.             ITER=ITER+1
  712. *       PRINT *,'ITER = ',ITER
  713.             IF(ITER.GT.ITEMAX) THEN
  714.                 CALL ERREUR(460)
  715.                 SEGSUP ALTRAV
  716.                 GO TO 9990
  717.             ENDIF
  718.             PSL=0.D0
  719.             SLON=0.0D0
  720.             FAI=0.0D0
  721.             IAM=0
  722. *
  723.             IF(LAST.EQ.0) THEN
  724. *         write(6,*) ' nam + 1',nam+1
  725. *         write(6,*) ' xlmbda acur(1)  acur(2) ang(1) ang(2) '
  726. *          write(6,*)  xlmbda,acur(1), acur(2) ,ang(1),ang(2)
  727. *             CALCUL DE L'ANGLE XALFA ASSOCIE A XLMBDA
  728.               CALL PRECN0(XLMBDA,ACUR,ANG,NAM+1,ANGBDA,ISUC)
  729.               IF(ISUC.EQ.NAM+1) ISUC=NAM
  730. *
  731. *****         PRINT *,'ISUC=',ISUC
  732.               IF(ISUC.LT.0) THEN
  733.                  CALL ERREUR(-ISUC)
  734.                  SEGSUP ALTRAV
  735.                  GO TO 9990
  736.               ENDIF
  737.  
  738.               SD=F1*ANGBDA + F2*XLMBDA
  739.               DSDDD=F1*DANG(ISUC)/DACUR(ISUC)+F2
  740. *
  741.             ELSE
  742.               SD=F1*FAIT + F2 * SLONT
  743.             ENDIF
  744. *
  745.             GG=0.D0
  746.             DGDD=0.D0
  747.             DO 8005 IC=1,NBELEM
  748.                if(ivers) then
  749.                   JC= NBELEM+1-IC
  750.                else
  751.                  JC=IC
  752.                endif
  753. *
  754.                IAM=IAM+1
  755.                ALFA = DANG(IAM)
  756.                DS1  = DACUR(IAM)
  757.  
  758.                do 8006  iptg=1,nbpgau
  759.  
  760.                if(ivers) then
  761.                ig= nbpgau+1-iptg
  762.                else
  763.                ig = iptg
  764.                endif
  765.  
  766. *         ---- Récupération des caractéristiques matérielles et géométriques ----
  767.  
  768.                melval =ivamat
  769.                ea = velche(min(ig,velche(/1)),
  770.      &                       min(JC,velche(/2)))
  771.  
  772.                melval =ivacar
  773.                sect1 = velche(min(ig,velche(/1)),
  774.      &                        min(JC,velche(/2)))
  775.  
  776.                SABS=SLON + DS1*(1+QSIGAU(IPTG))/2
  777.                IF(SABS.GT.XLMBDA) GO TO 8007
  778.  
  779.                SP=F1*FAI + F2*SABS
  780.                IF(LAST.EQ.0) THEN
  781.                  GG = GG +
  782.      &           PS1 * (EXP(-SP) - EXP(SP-2.D0*SD))*DS1*0.5D0
  783.      &                                   /(EA*SECT1)
  784.  
  785.                  DGDD = DGDD +
  786.      &           PS1 * (EXP(SP-2.D0*SD))*DSDDD*DS1/(EA*SECT1)
  787.                ELSE
  788.                  GG = GG + PS1 * (EXP(-SP))*DS1*0.5D0
  789.      &                                   /(EA*SECT1)
  790.                  DGDD = DGDD + (EXP(SP-SD))*DS1*0.5D0
  791.      &                                   /(EA*SECT1)
  792.                ENDIF
  793.  
  794.  8006          continue
  795.                SLON=SLON+DS1
  796.                FAI=FAI+ALFA
  797.  8005       CONTINUE
  798. *
  799. *        TEST DE CONVERGENCE
  800. *
  801.  8007       CONTINUE
  802. *****       PRINT *,'GG=',GG, ' GANC=',GANC,' FG =',GG-GANC
  803. *****       PRINT *,'DGDD=', DGDD
  804.  
  805.             FGG(ITER)=GG-GANC
  806.             XLAM(ITER)=XLMBDA
  807.  
  808. *
  809. *        CAS OU L'EFFET DE L'ANCRAGE DEPASSE LA LONGUEUR
  810. *        DU CABLE
  811. *
  812. *****    PRINT *,'LAST=',LAST
  813.          IF(LAST.EQ.1) THEN
  814.  
  815.              IF(ABS(GG-GANC)/GANC.LT.1.D-3) THEN
  816.                 PSL=0.D0
  817.              ELSE
  818.                 PSL=(GG-GANC)/DGDD
  819.              ENDIF
  820.  
  821. *****        PRINT *,'PSL=',PSL
  822.              IF(PSL.GT.PS1*EXP(-SD)) THEN
  823. *****           PRINT *,' ATTENTION PSL EST TROP GRAND '
  824.              ENDIF
  825.              GO TO 7999
  826.          ENDIF
  827.  
  828.  
  829.          IF(ABS(GG-GANC)/GANC.LT.1.D-3) GO TO 7999
  830. *
  831.          IF(ITER.GE.6 .AND.IDICH.EQ.0.AND.
  832.      &                     FGG(ITER)*FGG(ITER-1).LT.0.D0) THEN
  833.             IDICH=1
  834.             XLAM1=XLAM(ITER-1)
  835.             XLAM2=XLAM(ITER)
  836.             FG1=FGG(ITER-1)
  837.             FG2=FGG(ITER)
  838.          ENDIF
  839. *****    PRINT *,'IDICH=',IDICH
  840. *
  841.          IF(IDICH.EQ.0) THEN
  842.             XLMBDA=XLMBDA - (GG-GANC)/DGDD
  843.          ELSE
  844. *****       PRINT *,'XLAM(ITER)=',XLAM(ITER)
  845. *****       PRINT *, '  FGG(ITER)=',FGG(ITER)
  846.  
  847.             IF(FG1*FGG(ITER).LT.0.D0) THEN
  848.                XLAM2=XLAM(ITER)
  849.                FG2=FGG(ITER)
  850.             ELSE IF (FG1*FGG(ITER).GT.0.D0) THEN
  851.                XLAM1=XLAM(ITER)
  852.                FG1=FGG(ITER)
  853.             ENDIF
  854. *
  855. *****       PRINT *,'XLAM1=',XLAM1,'  FG1=',FG1
  856. *****       PRINT *,'XLAM2=',XLAM2,'  FG2=',FG2
  857.             DLAM=XLAM2-XLAM1
  858.             PROF=FG1/(FG1-FG2)
  859.             XLMBDA=XLAM1+PROF*DLAM
  860.          ENDIF
  861. *
  862. *****      PRINT *,'NOUVELLE XLMBDA = ', XLMBDA
  863.  
  864.            IF(XLMBDA.GE.SLONT) THEN
  865.               IF(ITER.EQ.ITEMAX-1) THEN
  866.                  LAST=1
  867.               ELSE
  868.                  XLMBDA=SLONT
  869.               ENDIF
  870.            ENDIF
  871. *
  872. *          si on trouve XLMBDA < 0 , on reinitialise
  873. *
  874.            IF(XLMBDA.LT.0.D0) THEN
  875.                 IF(ITER.GE.2) THEN
  876.                    XLMBDA=(XLAM(ITER)+XLAM(1))/ITER
  877.                 ELSE
  878.                    XLMBDA=XLAM(1)/2.D0
  879.                 ENDIF
  880. *****           PRINT *,'XLMBDA  REINITIALISEE = ', XLMBDA
  881.            ENDIF
  882. *
  883. *          SORTIE SI IDICH=1 ET XLMBDA CONVERGE
  884. *
  885.            IF(IDICH.EQ.1.AND.ITER.EQ.ITEMAX) THEN
  886.               IF(XLMBDA.NE.0.D0) THEN
  887.                  IF(ABS(XLMBDA-XLAM(ITER))/XLMBDA.
  888.      &                 LT.1.D-3) GO TO 7999
  889.               ENDIF
  890.            ENDIF
  891. *
  892.            GO TO 8000
  893. *
  894. *    2-EME BOUCLE SUR LES ELEMENTS
  895. *
  896.  7999       CONTINUE
  897.             SLON=0.0D0
  898.             FAI=0.0D0
  899.             IAM=0
  900.             DO 4005 IC=1,NBELEM
  901.  
  902.             if(ivers) then
  903. C    ordre inverse
  904.                JC= NBELEM+1-IC
  905.             else
  906. C    ordre normal
  907.               JC=IC
  908.             endif
  909.  
  910. C     ---- ON RECUPERE LES GRANDEURS GEOMETRIQUES -----------
  911.  
  912.                IAM=IAM+1
  913.                ALFA = DANG(IAM)
  914.                DS1  = DACUR(IAM)
  915.  
  916.                do 12345 iptg=1,nbpgau
  917.  
  918.                if(ivers) then
  919.                ig= nbpgau+1-iptg
  920.                else
  921.                ig = iptg
  922.                endif
  923.  
  924. C         ---- Récupération des EFFX deja présents -----------
  925.  
  926.                IF (IPCHC1.EQ.0) THEN
  927.                   SINI = 0.D0
  928.                ELSE
  929.                   mptval = ivasi0
  930.                   MELVAL = IVAL(ipeffx)
  931.                   SINI =velche(min(ig,velche(/1)),
  932.      &                         min(JC,velche(/2)))
  933.                ENDIF
  934.  
  935. C         ---- Récupération des caractéristiques matérielles et géométriques ----
  936.  
  937.                melval =ivamat
  938.                ea = velche(min(ig,velche(/1)),
  939.      &                       min(JC,velche(/2)))
  940.  
  941.                melval =ivacar
  942.                sect1 = velche(min(ig,velche(/1)),
  943.      &                        min(JC,velche(/2)))
  944.  
  945. C          --- Calcul des pertes de précontrainte -------------------
  946.  
  947.                ICOMP = 0
  948.                call precn1(PS1,ea,f1,f2,ganc,RMU0,FPRG,RH10,
  949.      &                     SLON+DS1*(1+qsigau(iptg))/2,FAI,SINI,
  950.      &                     SECT1,XLMBDA,SD,LAST,PSL,PSOUT,ICOMP)
  951.                IF (ICOMP.EQ.1) THEN
  952.                   CALL ERREUR(768)
  953.                   goto 9990
  954.                   RETURN
  955.                ENDIF
  956.  
  957. C          --- On range les résultats
  958.  
  959.                iret=1
  960.                MPTVAL=IVASTR
  961.                DO 1700 I=1,NSTRS
  962.                   MELVAL=IVAL(I)
  963.                   VELCHE(ig,JC)=PSOUT
  964.  1700          CONTINUE
  965.  
  966. 12345          continue
  967.                SLON=SLON+DS1
  968.                FAI=FAI+ALFA
  969. C            write(6,*) 'slon et fai ',slon,fai
  970.  
  971.  4005       CONTINUE
  972.             SEGSUP ALTRAV
  973. C----------------------    endif  sur  CERC ou  ELEMENTS LINEAIRE
  974.       END IF
  975.  
  976. C
  977.  9990 CONTINUE
  978.  
  979.       SEGSUP WRK3
  980.       IF(IERR.NE.0)THEN
  981.          SEGSUP MCHAML
  982.       ENDIF
  983. c
  984. C      CALL DTMVAL(IVACAR,1)
  985.  
  986.       if(IPCHC1.gt.0) call dtmval(ivasi0,1)
  987. C
  988.       IF(IERR.NE.0)THEN
  989.           CALL DTMVAL(IVASTR,3)
  990.       ELSE
  991.           CALL DTMVAL(IVASTR,1)
  992.       ENDIF
  993. C
  994.       IF(MOMATR.NE.0)THEN
  995.          NOMID=MOMATR
  996.          SEGSUP NOMID
  997.       ENDIF
  998.       IF(MOCARA.NE.0)THEN
  999.          NOMID=MOCARA
  1000.          SEGSUP NOMID
  1001.       ENDIF
  1002.       IF(MOSTRS.NE.0)THEN
  1003.          NOMID=MOSTRS
  1004.          if(lsupco)SEGSUP NOMID
  1005.       ENDIF
  1006. C
  1007.       SEGSUP INFO
  1008.       if(iret.eq.0) return
  1009. C
  1010.   500 CONTINUE
  1011. c-------- Fin de la boucle sur les zones du MODELE  TEMPORAIRE
  1012. C
  1013. C =========== rangement dans le MCHELM  final   MCHEL5
  1014.       if(ifomod.ne.0.and.ifomod.ne.1) then
  1015. C           les cables standards
  1016.       segact siezo
  1017. C      segact mchel5
  1018.       MCHELM = IPSTRS
  1019.       segact MCHELM
  1020. C
  1021. C  boucle sur les zones du modele ou champ initial
  1022.  
  1023. C
  1024.         ibc= 0
  1025.       do 3010 isous=1,mchel5.imache(/1)
  1026.       sielc = iezon(isous)
  1027.       segact sielc
  1028. C
  1029.         mcham5= mchel5.ichaml(isous)
  1030.         segact mcham5
  1031.         melva5 =mcham5.ielval(1)
  1032.         segact melva5*mod
  1033. C  boucle sur les partitions provisoires de la sz
  1034. C
  1035.       inel = 0
  1036.         do 3020 ik=1,nbcz
  1037.         ibc=ibc+1
  1038.         mchaml=ichaml(ibc)
  1039.         segact mchaml
  1040.         melval=ielval(1)
  1041.         segact melval
  1042. C
  1043.       do iii=1,velche(/2)
  1044.       inel = inel+1
  1045.       iel=isens(1,inel)
  1046.         if(isens(2,inel).eq.1) then
  1047.       do ip=1,velche(/1)
  1048.       melva5.velche(ip,iel)=velche(ip,iii)
  1049.       enddo
  1050.         else
  1051.       do ip=1,velche(/1)
  1052.       melva5.velche(3-ip,iel)=velche(ip,iii)
  1053.       enddo
  1054.         endif
  1055.       enddo
  1056. 3020  continue
  1057.       segsup sielc
  1058. 3010  continue
  1059.       call DTCHAM(ipcar2)
  1060.       if(ipchc2.ne.0) call DTCHAM(ipchc2)
  1061.       segsup ipt5
  1062.       else
  1063. C          les cerc
  1064.       icha1 = MCHEL5
  1065.       MCHEL5 = MCHELM
  1066.       MCHELM = icha1
  1067.       endif
  1068. 333   format(i4,2E12.5)
  1069. C  maintenant les destructions d objets temporaires
  1070.       call DTCHAM(MCHELM)
  1071. C      call DTMODL(ipmodt)
  1072.       mmodel = ipmodt
  1073.       segsup mmodel
  1074. C
  1075.       IPSTRS =MCHEL5
  1076.       iret = 1
  1077.       END
  1078.  
  1079.  
  1080.  
  1081.  
  1082.  
  1083.  
  1084.  

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