Télécharger la plaquette Cast3M


fr en

Télécharger ksigmp.eso

Retour à la liste

Numérotation des lignes : 

ksigmp

  1. C KSIGMP    SOURCE    PV090527  24/04/22    21:15:02     11913          
  2.  
  3.       SUBROUTINE KSIGMP(IPMODL,IPCHE1,IPCHE2,IFLAM, IPRIGG)
  4.  
  5. c_______________________________________________________________________
  6. c
  7. c
  8. c     construction de la matrice de raideur geometrique a partir d'un
  9. c     mchaml de contraintes
  10. c
  11. c     entr{es:
  12. c     ________
  13. c
  14. c       ipmodl  pointeur sur un mmodel
  15. c       ipche1  pointeur sur un mchaml de contraintes
  16. c       ipche2  pointeur sur un mchaml de caracteristiques
  17. c       iflam   flag de flambage
  18. c
  19. c     sorties:
  20. c     ________
  21. c
  22. c       iprigg  pointeur sur un objet rigidite
  23. c               = 0 en cas d'erreur
  24. c_______________________________________________________________________
  25. c
  26.       IMPLICIT INTEGER(I-N)
  27.       IMPLICIT REAL*8(A-H,O-Z)
  28.  
  29.  
  30. -INC PPARAM
  31. -INC CCOPTIO
  32. -INC CCHAMP
  33. -INC CCREEL
  34.  
  35. -INC SMCHAML
  36. -INC SMCOORD
  37. -INC SMELEME
  38. -INC SMINTE
  39. -INC SMMODEL
  40. -INC SMRIGID
  41. C
  42.       INTEGER oooval
  43.  
  44.       SEGMENT NOTYPE
  45.         CHARACTER*16 TYPE(NBTYPE)
  46.       ENDSEGMENT
  47. C
  48.       SEGMENT MPTVAL
  49.         INTEGER IPOS(NS)    ,NSOF(NS)
  50.         INTEGER      IVAL(NCOSOU)
  51.         CHARACTER*16 TYVAL(NCOSOU)
  52.       ENDSEGMENT
  53. c
  54.       SEGMENT MWRK1
  55.         REAL*8 REL(LRE,LRE) ,XE(3,NBBB) ,XSTRS(NSTRS)
  56.       ENDSEGMENT
  57. C
  58.       SEGMENT MWRK2
  59.         REAL*8 SHPWRK(6,NBNO) ,BGENE(NSTRS,LRE)
  60.       ENDSEGMENT
  61. C
  62.       SEGMENT MWRK3
  63.         REAL*8 WORK(LW)
  64.       ENDSEGMENT
  65. C
  66.       SEGMENT MWRK4
  67.         REAL*8 BPSS(3,3) ,XEL(3,NBBB)
  68.       ENDSEGMENT
  69. C
  70.       SEGMENT MWRK5
  71.         REAL*8 GEOM(20), tabw(6,9), tabrot(4,9), XX(3), YY(3)
  72.       ENDSEGMENT
  73. C
  74. C   segment pour shb8
  75.       SEGMENT MWRK7
  76.         REAL*8 PROPEL(1),out(1),d(1), work1(30)
  77.       ENDSEGMENT
  78. C
  79.       character*6 msorse
  80.       CHARACTER*8 CMATE
  81.       CHARACTER*(NCONCH) CONM
  82.       PARAMETER ( NINF=3 )
  83.       INTEGER INFOS(NINF)
  84.       LOGICAL lsupfo,lsupde,lsupco,BDPGE
  85.       INTEGER ISUP1,ISUP2
  86.  
  87.       ISUP1 = 0
  88.       ISUP2 = 0
  89.  
  90.       IPRIGG = 0
  91.  
  92.       IDIMP1 = IDIM+1
  93. C
  94. C      verification du lieu support du mchaml de contraintes
  95. C
  96.       CALL QUESUP(IPMODL,IPCHE1,3,0,ISUP1,IRET1C)
  97.       IF (ISUP1.GT.1) RETURN
  98. C
  99. C      verification du lieu support du mchaml de caracteristiques
  100. C
  101.       IF (IPCHE2.NE.0) THEN
  102.         CALL QUESUP(IPMODL,IPCHE2,3,0,ISUP2,iret2c)
  103.         IF (ISUP2.GT.1) RETURN
  104.       ENDIF
  105. c
  106. c_______________________________________________________________________
  107. c
  108. c     initialisation du chapeau de l objet rigidite
  109. c_______________________________________________________________________
  110. c
  111.       NRIGEL = 0
  112.       SEGINI,MRIGID
  113.       IFORIG = IFOUR
  114.       ICHOLE = 0
  115.       IMGEO1 = 0
  116.       IMGEO2 = 0
  117.       ISUPEQ = 0
  118.       IF (IFLAM.NE.0) THEN
  119.         MTYMAT = 'MASSE   '
  120.       ELSE
  121.         MTYMAT = 'RIGIDITE'
  122.       ENDIF
  123. c
  124. c_______________________________________________________________________
  125. c
  126. c     activation du modele
  127. c_______________________________________________________________________
  128. c
  129.       MMODEL = IPMODL
  130.       SEGACT,MMODEL
  131.       NSOUS = KMODEL(/1)
  132. c
  133. c      boucle sur les modeles elementaires
  134. c
  135.       DO 500 ISOUS = 1,NSOUS
  136. c
  137. c     traitement du modele
  138. c
  139.         IMODEL = KMODEL(ISOUS)
  140.         IF (FORMOD(1).EQ.'CHARGEMENT') GOTO 500
  141.         SEGACT,IMODEL
  142. C
  143. C     INITIALISATIONS
  144. C
  145.         IPMINT = 0
  146.         IPMIN1 = 0
  147.  
  148.         MOSTRS = 0
  149.         MOCARA = 0
  150.         MOTYPS = 0
  151.         MOTYPC = 0
  152.  
  153.         MODEPL = 0
  154.         MOFORC = 0
  155.         lsupde = .false.
  156.         lsupfo = .false.
  157.         lsupco = .false.
  158.  
  159.         IDESCR = 0
  160.  
  161. C- Recuperation d'informations sur le maillage elementaire
  162.         IPT1 = IMAMOD
  163.         SEGACT,IPT1
  164.         NBNOE1 = IPT1.NUM(/1)
  165.         NBELE1 = IPT1.NUM(/2)
  166.  
  167. C- Quelques informations sur le modele
  168.         IIPDPG = imodel.IPDPGE
  169.         IIPDPG = IPTPOI(IIPDPG)
  170.  
  171.         CONM  = CONMOD
  172.         CMATE = CMATEE
  173. C       MATE  = IMATEE
  174. c*      INAT  = INATUU
  175.  
  176.         IRTD = 1
  177.         CALL IDENT(IPT1,CONM,IPCHE1,IPCHE2, INFOS,IRTD)
  178.         IF (IRTD.EQ.0) GOTO 599
  179. C
  180. C- Recuperation d'informations sur l'element fini
  181.         MELE  = NEFMOD
  182. c  pour l'el. timo on utilise l'el. barr
  183. c         IF (MELE .EQ. 84) MELE = 46
  184. c  bp: comme il n y a plus elquoi, ce n'est pas ici que ca intervient...
  185.  
  186. c   coque integree ou pas ?
  187.         IF (INFMOD(/1).NE.0)THEN
  188.           NPINT = INFMOD(1)
  189.           IF (NPINT.NE.0) THEN
  190.             CALL ERREUR(615)
  191.             GOTO 599
  192.           ENDIF
  193.         ELSE
  194.           NPINT = 0
  195.         ENDIF
  196.  
  197. C       LHOOK = INFELE(10)
  198. c*      LHOO2 = LHOOK*LHOOK
  199.         NSTRS = INFELE(16)
  200.         MFR   = INFELE(13)
  201.         LW    = INFELE(7)
  202. C        NDDL  = INFELE(15)
  203.         LRE   = INFELE(9)
  204. C        IPORE = INFELE(8)
  205.         NHRM  = NIFOUR
  206.  
  207.         IPPORE = 0
  208.         IF (MFR.EQ.33) IPPORE = NBNOE1
  209.  
  210. c_______________________________________________________________________
  211. C   segments d'integration     *
  212. c_______________________________________________________________________
  213. C  minte : 1er  segment d'integration, il existe pour tous les e.f.
  214. C  minte1: 2eme segment d'integration, uniquement pour certains e.f.
  215. C          en particulier pour coq6 et coq8
  216. C  nbpg:nb de points de gauss = nbpgau du segment minte
  217. C  iele:no d'element geometrique associe a l'e.f. mele
  218. C  nbff:nb de fonctions de forme = nbno du segment minte
  219.         NBPGAU = INFELE(  6)
  220. C        IELE   = INFELE( 14)
  221. c*        ICARA  = INFELE(  5)
  222.         IPMINT = INFMOD(5)
  223. c*        IPMINT  = INFELE(11)
  224.         IPMIN1 = INFMOD(8)
  225.         MINTE = IPMINT
  226.         IF (IPMINT.NE.0) SEGACT,MINTE
  227.  
  228. c_______________________________________________________________________
  229. c
  230. C  initialisation du segment descr, segment descripteur des *
  231. C  des inconnues relatives a la matrice de rigidite         *
  232. c_______________________________________________________________________
  233.         if (lnomid(1).ne.0) then
  234.           MODEPL = lnomid(1)
  235.         else
  236.           lsupde = .true.
  237.           CALL IDPRIM(IMODEL,MFR,MODEPL,NDEPL,NDUM)
  238.         endif
  239.         nomid = MODEPL
  240.         segact,nomid
  241.         ndepl = lesobl(/2)
  242. c*        ndum  = lesfac(/2)
  243.  
  244.         if (lnomid(2).ne.0) then
  245.           moforc = lnomid(2)
  246.         else
  247.           lsupfo=.true.
  248.           CALL IDDUAL(IMODEL,MFR,MOFORC,NFORC,NDUM)
  249.         endif
  250.         nomid = MOFORC
  251.         segact,nomid
  252.         nforc = lesobl(/2)
  253. c*      ndum  = lesfac(/2)
  254.  
  255.         IF (NDEPL.EQ.0.OR.NFORC.EQ.0.OR.NDEPL.NE.NFORC) CALL ERREUR(5)
  256.  
  257. C CB215821 : On copie RIGI2 pour les deformations planes generalisees
  258.         IF (IIPDPG.GT.0) THEN
  259.           IF (IFOUR.EQ.-3) THEN
  260.             BDPGE=.TRUE.
  261.             IREF=(IIPDPG-1)*(IDIM+1)
  262. C            XDPGE=XCOOR(IREF+1)
  263. C            YDPGE=XCOOR(IREF+2)
  264.           ELSE IF (IFOUR.EQ. 7 .OR. IFOUR.EQ. 8 .OR. IFOUR.EQ. 9 .OR.
  265.      &             IFOUR.EQ.10 .OR. IFOUR.EQ.11 .OR. IFOUR.EQ.14) THEN
  266.             BDPGE=.TRUE.
  267. C            XDPGE=XZero
  268. C            YDPGE=XZero
  269.           ELSE
  270.             CALL ERREUR(21)
  271.             RETURN
  272.           ENDIF
  273.         ELSE
  274.           BDPGE=.FALSE.
  275. C          XDPGE=XZero
  276. C          YDPGE=XZero
  277.         ENDIF
  278.  
  279.         NCOMP=NDEPL
  280.         NBNNS=NBNOE1
  281.         IF (MFR.EQ.33)                NCOMP = NDEPL-1
  282.         IF (MFR.EQ.19 .OR. MFR.EQ.21) NBNNS = NBNOE1/2
  283. C
  284. C   REMPLISSAGE DU SEGMENT DESCRIPTEUR
  285.         IF(BDPGE)THEN
  286. C         CB215821 : En 'DPGE' on traite sans le Pt support
  287.           NCOMP = NCOMP - 3
  288.           LRE   = LRE   - 3
  289.         ENDIF
  290.         NLIGRP = LRE
  291.         NLIGRD = LRE
  292.         SEGINI,DESCR
  293.  
  294.         IDDL = 1
  295.         DO 1004 INOEUD=1,NBNNS
  296.           DO 1005 ICOMP=1,NCOMP
  297.             NOMID=MODEPL
  298.             LISINC(IDDL)=LESOBL(ICOMP)
  299.             NOMID=MOFORC
  300.             LISDUA(IDDL)=LESOBL(ICOMP)
  301.             NOELEP(IDDL)=INOEUD
  302.             NOELED(IDDL)=INOEUD
  303.             IDDL=IDDL+1
  304.  1005     CONTINUE
  305.  1004   CONTINUE
  306.  
  307. C           cas des milieux poreux
  308. C
  309. C         if (mfr.eq.33) then
  310. C           ipos = nspos(iele)
  311. C           do 1104 inoeud=1,nbsom(iele)
  312. C             nomid=modepl
  313. C             lisinc(iddl)=lesobl(ndepl)
  314. C             nomid=moforc
  315. C             lisdua(iddl)=lesobl(ndepl)
  316. C             i = ibsom(ipos+inoeud-1)
  317. C             noelep(iddl)=i
  318. C             noeled(iddl)=i
  319. C             iddl=iddl+1
  320. C 1104      continue
  321. C         endif
  322.  
  323. C           cas des element raccord
  324. C
  325.         IF (MFR.EQ.19.OR.MFR.EQ.21) THEN
  326.           CALL IDPRIM(IMODEL,MFR+1000,MODPL1,NDEPL,NDUM)
  327.           CALL IDDUAL(IMODEL,MFR+1000,MOFRC1,NFORC,NDUM)
  328.           DO 1106 INOEUD=NBNNS+1,NBNOE1
  329.             DO 1107 ICOMP=1,NDEPL
  330.               NOMID=MODPL1
  331.               LISINC(IDDL)=LESOBL(ICOMP)
  332.               NOMID=MOFRC1
  333.               LISDUA(IDDL)=LESOBL(ICOMP)
  334.               NOELEP(IDDL)=INOEUD
  335.               NOELED(IDDL)=INOEUD
  336.               IDDL=IDDL+1
  337. 1107        continue
  338. 1106      continue
  339.           NOMID=MODPL1
  340.           SEGSUP,NOMID
  341.           NOMID=MOFRC1
  342.           SEGSUP,NOMID
  343.         ENDIF
  344.  
  345.         SEGDES,DESCR
  346.         IDESCR = DESCR
  347. c_______________________________________________________________________
  348. c
  349. C  composantes de contraintes necessaires                      *
  350. c_______________________________________________________________________
  351.         if (lnomid(4).ne.0) then
  352.           MOSTRS = lnomid(4)
  353.         else
  354.           lsupco=.true.
  355.           CALL IDCONT(IMODEL,IFOUR,MOSTRS,NSTR,NFAC)
  356.         endif
  357.         nomid = MOSTRS
  358.         segact,nomid
  359.         nstr=lesobl(/2)
  360. c*        nfac=lesfac(/2)
  361. c*        write(6,*) 'mostrts',mostrs,nstr,nfac
  362.         nbtype = 1
  363.         SEGINI,notype
  364.         TYPE(1)='REAL*8'
  365.         MOTYPS = notype
  366.  
  367.         ifai = 1
  368.         if (mele.eq.260.and.IRET1C.eq.5) ifai = 0
  369.         ISUP1L = 0
  370.         IF (ISUP1.EQ.1.AND.ifai.eq.1) ISUP1L = 1
  371.  
  372. c____________________________________________________________________
  373. c
  374. C  traitement des champs de caracteristiques                   *
  375. c____________________________________________________________________
  376.         NBROBL = 0
  377.         NBRFAC = 0
  378.         IVECT  = 0
  379.         notype = 0
  380.         nomid  = 0
  381. C
  382. C       v1x v1y dans le cas de la coque dst orthotrope
  383. C
  384.         IF (MFR.EQ.9) THEN
  385.           IF (MELE.EQ.93.AND.CMATE.NE.'ISOTROPE')THEN
  386.             NBROBL=2
  387.             SEGINI NOMID
  388.             LESOBL(1)='V1X '
  389.             LESOBL(2)='V1Y '
  390. C
  391.             NBTYPE=1
  392.             SEGINI NOTYPE
  393.             TYPE(1)='REAL*8'
  394.           ENDIF
  395. C
  396. C  epaisseur dans le cas massif et coq2 en contraintes planes
  397. C
  398.         ELSE IF ( (MFR.EQ.1.OR.MFR.EQ.3.OR.MFR.EQ.31) .AND.
  399.      +            IFOUR.EQ.-2 .AND. IPCHE2.NE.0) THEN
  400.           NBRFAC=1
  401.           SEGINI NOMID
  402.           LESFAC(1)='DIM3'
  403. C
  404.           NBTYPE=1
  405.           SEGINI NOTYPE
  406.           TYPE(1)='REAL*8'
  407. C
  408. C  epaisseur et excentrement dans le cas des coques epaisses
  409. C
  410.         ELSE IF (MFR.EQ.5 .OR. (MFR.EQ.3.AND.IFOUR.NE.-2)) THEN
  411.           NBROBL=1
  412.           NBRFAC=1
  413.           SEGINI NOMID
  414.           LESOBL(1)='EPAI'
  415.           LESFAC(1)='EXCE'
  416.  
  417.           NBTYPE=1
  418.           SEGINI NOTYPE
  419.           TYPE(1)='REAL*8'
  420. C
  421. C  caracteristiques pour les poutres
  422. C
  423.         ELSE IF (MFR.EQ.7 ) THEN
  424.           IF (IFOUR.EQ.-2.OR.IFOUR.EQ.-1.OR.IFOUR.EQ.-3) THEN
  425.             NBROBL=2
  426.             NBRFAC=1
  427.             SEGINI NOMID
  428.             LESOBL(1)='SECT'
  429.             LESOBL(2)='INRZ'
  430.             LESFAC(1)='SECY'
  431. C
  432.             NBTYPE=1
  433.             SEGINI NOTYPE
  434.             TYPE(1)='REAL*8'
  435.           ELSE
  436.             NBROBL=4
  437.             NBRFAC=5
  438.             IVECT =1
  439.             SEGINI NOMID
  440.             LESOBL(1)='TORS'
  441.             LESOBL(2)='INRY'
  442.             LESOBL(3)='INRZ'
  443.             LESOBL(4)='SECT'
  444.             LESFAC(1)='SECY'
  445.             LESFAC(2)='SECZ'
  446.             LESFAC(3)='VX  '
  447.             LESFAC(4)='VY  '
  448.             LESFAC(5)='VZ  '
  449. C
  450.             NBTYPE=9 
  451.             SEGINI NOTYPE
  452.             TYPE(1)='REAL*8'
  453.             TYPE(2)='REAL*8'
  454.             TYPE(3)='REAL*8'
  455.             TYPE(4)='REAL*8'
  456.             TYPE(5)='REAL*8'
  457.             TYPE(6)='REAL*8'
  458.             TYPE(7)='REAL*8'
  459.             TYPE(8)='REAL*8'
  460.             TYPE(9)='REAL*8'
  461.           ENDIF
  462. C
  463. C  caracteristiques pour les tuyaux
  464. C
  465.         ELSE IF (MFR.EQ.13) THEN
  466.           NBROBL = 2
  467.           NBRFAC = 5
  468.           IVECT  = 1
  469.           SEGINI NOMID
  470.           LESOBL(1)='EPAI'
  471.           LESOBL(2)='RAYO'
  472.           LESFAC(1)='RACO'
  473.           LESFAC(2)='CISA'
  474.           LESFAC(3)='VX  '
  475.           LESFAC(4)='VY  '
  476.           LESFAC(5)='VZ  '
  477. C
  478.           NBTYPE = 7
  479.           SEGINI NOTYPE
  480.           TYPE(1)='REAL*8'
  481.           TYPE(2)='REAL*8'
  482.           TYPE(3)='REAL*8'
  483.           TYPE(4)='REAL*8'
  484.           TYPE(5)='REAL*8'
  485.           TYPE(6)='REAL*8'
  486.           TYPE(7)='REAL*8'
  487.         ENDIF
  488. C
  489.         MOCARA = NOMID
  490.         MOTYPC = NOTYPE
  491.         NCARA = NBROBL
  492.         NCARF = NBRFAC
  493.         NCARR = NCARA+NCARF
  494. **      write(6,*) 'mfr ncarr en 498',mfr,ncarr,ncara,ncarf
  495.  
  496.         IF (MOCARA.NE.0 .AND. IPCHE2.EQ.0) THEN
  497.           MOTERR(1:8)  = 'CARACTER'
  498.           MOTERR(9:12) = NOMTP(MELE)
  499.           MOTERR(13:20)= 'KSIGMA'
  500.           CALL ERREUR(145)
  501.           GOTO 598
  502.         ENDIF
  503.  
  504.         ifai = 1
  505.         IF (mele.EQ.260) ifai = 0
  506.         ISUP2L = 0
  507.         IF (ISUP2.EQ.1.AND.ifai.eq.1) ISUP2L = 1
  508.  
  509. C- Partionnement si necessaire de la matrice de capacite
  510. C- determinant ainsi le nombre d'objets elementaires de MRIGID
  511.         LTRK = oooval(1,4)
  512.         IF (LTRK.EQ.0) LTRK = oooval(1,1)
  513.       LTRK=MAX(LTRK,2**24)
  514. C  Ajout a la taille en mots de la matrice des infos du segment
  515.         LSEG = LRE*LRE*NBELE1 + 16
  516.         NBLPRT = (LSEG-1)/LTRK + 1
  517.         NBLMAX = (NBELE1-1)/NBLPRT + 1
  518.         NBLPRT = (NBELE1-1)/NBLMAX + 1
  519. C         write(ioimp,*) ' capa1 : nblprt nblmax = ',nblprt,nblmax,nbele1
  520.  
  521. C Ajout de la matrice a la matrice globale
  522. C ========================================
  523. C        NRIGE0 =  IRIGEL(/2)
  524. C         NRIGEL = NRIGE0 + NBLPRT
  525. C         SEGADJ,MRIGID
  526.  
  527.         descr  = IDESCR
  528.         meleme = IPT1
  529.         NBNN   = NBNOE1
  530.         nbelem = NBELE1
  531.         nbsous = 0
  532.         nbref  = 0
  533. C
  534. C ***********************************************************************
  535. C                           P H A S E   2
  536. C
  537. C  Boucle sur les PARTITIONS elementaires de la matrice
  538. C
  539. C ***********************************************************************
  540.         DO irige = 1, NBLPRT
  541.  
  542.           IF (NBLPRT.GT.1) THEN
  543. C- Partitionnement du maillage support de la matrice elementaire
  544. C- (IPT1 peut etre desactive suite a l'appel a KOMCHA !)
  545.             SEGACT,IPT1
  546.             ielem  = (irige-1)*NBLMAX
  547.             nbelem = MIN(NBLMAX,NBELE1-ielem)
  548. C           write(ioimp,*) ' creation segment ',nbnn,nbelem
  549.             SEGINI,meleme
  550.             itypel = IPT1.itypel
  551.             DO ielt = 1, nbelem
  552.               jelt = ielt + ielem
  553.               DO inoe = 1, NBNN
  554.                 num(inoe,ielt) = IPT1.NUM(inoe,jelt)
  555.               ENDDO
  556.               icolor(ielt) = IPT1.ICOLOR(jelt)
  557.             ENDDO
  558. C- Recopie du descripteur
  559.             des1 = IDESCR
  560.             SEGINI,descr=des1
  561.             SEGDES,descr
  562.           ENDIF
  563.  
  564.           ipmail = meleme
  565.           ipdesc = descr
  566.  
  567. C- Initialisation de la matrice de rigidite elementaire (xmatri)
  568.           NELRIG = nbelem
  569.           SEGINI,xmatri
  570.           ipmatr = xmatri
  571.  
  572. C- Recuperation des valeurs des contraintes et proprietes geometriques
  573.           IVASTR = 0
  574.           IVACAR = 0
  575.           IVECTL = IVECT
  576.           NCARR1 = NCARR
  577. C
  578.           CALL KOMCHA(IPCHE1,ipmail,CONM,MOSTRS,MOTYPS,1,INFOS,3,IVASTR)
  579.           IF (IERR.NE.0) GOTO 597
  580.           IF (ISUP1L.EQ.1) THEN
  581.             CALL VALCHE(IVASTR,NSTR,IPMINT,IPPORE,MOSTRS,MELE)
  582.             IF (IERR.NE.0) THEN
  583.               ISUP1L = 0
  584.               GOTO 597
  585.             ENDIF
  586.           ENDIF
  587.  
  588.           IF (MOCARA.NE.0)  THEN
  589.             CALL KOMCHA(IPCHE2,ipmail,CONM,MOCARA,MOTYPC,1,
  590.      &                  INFOS,3,IVACAR)
  591.             IF (IERR.NE.0) GOTO 597
  592.             IF (ISUP2L.EQ.1) THEN
  593.               CALL VALCHE(IVACAR,NCARR,IPMINT,IPPORE,MOCARA,MELE)
  594.               IF (IERR.NE.0) THEN
  595.                 ISUP2L = 0
  596.                 GOTO 597
  597.               ENDIF
  598.             ENDIF
  599.             IF (IVECT.EQ.1) THEN
  600. ***           MPTVAL = IVACAR
  601. *****         NCARR1 = NCARR - 3
  602.  
  603. ***           IF (IVAL(NCARR1).EQ.0) IVECTL = 2
  604.             ENDIF
  605.           ENDIF
  606. **        write(6,*) ' dans ksigmp ivect ivectl ',ivect,ivectl
  607.  
  608. c_______________________________________________________________________
  609. c
  610. c     numero des etiquettes      :
  611. c     etiquettes de 1 a 98 pour traitement specifique a l element
  612. c     dans la zone specifique a chaque element commencant par :
  613. c     5  continue
  614. c     element 5   etiquettes 1005 2005 3005 4005   ...
  615. c     44 continue
  616. c     element 44  etiquettes 1044 2044 3044 4044   ...
  617. c_______________________________________________________________________
  618.           IF (MELE.LE.100) THEN
  619.       GOTO (99,99,99, 4,99, 4,99, 4,99, 4,99,99,99, 4, 4, 4, 4,99,99,99,
  620.      1      99,99, 4, 4, 4, 4,27,28,29,99,99,99,99,99,99,99,99,99,99,99,
  621.      2      41,29,43,44,99,46,99,99,49,99,51,99,99,99,99,41,99,99,99,99,
  622.      3      99,99,99,99,99,99,99,99, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4,99,99,
  623.      4      99,99,99,29,99,99,99,99,99,99,99,99,28,99,46,99,99,99,99,99
  624.      5     ),MELE
  625.           ELSE IF (MELE.LE.200) THEN
  626.       GOTO (99,99,99,99,99,99,99,99,99,99, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4,
  627.      1       4, 4,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,
  628.      2      99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,
  629.      3      99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,
  630.      4      99,99, 4, 4,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99
  631.      5     ),MELE-100
  632.           ELSE IF (MELE.LE.300) THEN
  633.       GOTO (99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,
  634.      1      99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,
  635.      2      99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,
  636.      &      260,
  637.      3      99,99,99,99,99,99,99,99,99,99,99,99,99, 4, 4,99,99,99,99,99,
  638.      4      99,99, 4, 4,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99,99
  639.      5     ),MELE-200
  640.           ENDIF
  641. c
  642.  99       CONTINUE
  643. C           MOTERR(1:4)  = NOMTP(MELE)
  644. C           MOTERR(5:12) = 'KSIGMP  '
  645. C          CALL ERREUR(86)
  646.           GOTO 510
  647. c
  648. c_______________________________________________________________________
  649. c
  650. c     secteur de calcul pour les elements massifs
  651. c_______________________________________________________________________
  652.  4        CONTINUE
  653.           NBNO = NBNN
  654.           NBBB = NBNN
  655.           SEGINI,MWRK1,MWRK2
  656. c   recuperation de l'epaisseur
  657.           DIM3 = 1.D0
  658.           MEPDI3 = 0
  659. c*         IF (IFOUR.EQ.-2.AND.IPCHE2.NE.0) THEN
  660.           IF (IVACAR.NE.0) THEN
  661.             MPTVAL = IVACAR
  662.             MEPDI3 = IVAL(1)
  663.           ENDIF
  664.  
  665.           DO 3004 IB=1,NBELEM
  666. c
  667. c     on cherche  les coordonnees des noeuds de l element ib
  668. c
  669.             CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  670.             CALL ZERO(REL,LRE,LRE)
  671. c
  672. c     boucle sur les points de gauss
  673. c
  674.             ISDJC = 0
  675.             DO 4004 IGAU=1,NBPGAU
  676. c
  677. c   recuperation de l'epaisseur
  678.               IF (MEPDI3.NE.0) THEN
  679.                 MELVAL = MEPDI3
  680.                 IGMN=MIN(IGAU,VELCHE(/1))
  681.                 IBMN=MIN(  IB,VELCHE(/2))
  682.                 DIM3=VELCHE(IGMN,IBMN)
  683.               ENDIF
  684. c
  685.               DO 100 IA=1,NBNN
  686.                 DO 101 IO=1,IDIMP1
  687.                   SHPWRK(IO,IA)=SHPTOT(IO,IA,IGAU)
  688.  101            CONTINUE
  689.  100          CONTINUE
  690.               CALL DEVOLU(XE,SHPWRK,MFR,NBNN,IFOUR,NIFOUR,IDIM,DIM3,
  691.      &                    RR,DJAC)
  692. c
  693. c     verification du signe du jacobien
  694. c
  695.               IF (DJAC.LT.0.) ISDJC=ISDJC+1
  696.               DJAC = ABS(DJAC)
  697.               IF (DJAC.LT.XPETIT) THEN
  698.                 INTERR(1) = IB
  699.                 CALL ERREUR(259)
  700.                 GOTO 9004
  701.               ENDIF
  702.               DJAC = DJAC * POIGAU(IGAU)
  703. c
  704. c     on recupere les contraintes
  705. c
  706.               MPTVAL=IVASTR
  707.               DO 5004 ICOMP=1,NSTR
  708.                 MELVAL=IVAL(ICOMP)
  709.                 IGMN = MIN(IGAU,VELCHE(/1))
  710.                 IBMN = MIN(IB  ,VELCHE(/2))
  711.                 XSTRS(ICOMP)=VELCHE(IGMN,IBMN)
  712.  5004         CONTINUE
  713. c
  714.               IF (IFOUR.EQ.1) THEN
  715.                 IF (NIFOUR.EQ.0) THEN
  716.                   CALL THSIG1(SHPWRK,DJAC,XSTRS,NBNN,LRE,REL,RR)
  717.                 ELSE
  718.                   CALL THSIG2(SHPWRK,DJAC,XSTRS,NBNN,LRE,REL,NIFOUR,RR)
  719.                 ENDIF
  720.               ELSE IF (IFOUR.EQ.0) THEN
  721.                 CALL THSIG3(SHPWRK,DJAC,XSTRS,NBNN,LRE,REL,RR)
  722.               ELSE
  723.                 CALL THSIGH(SHPWRK,DJAC,XSTRS,NBNN,IDIM,LRE,REL)
  724.               ENDIF
  725. c
  726.  4004       CONTINUE
  727.  
  728.             IF (ISDJC.NE.0.AND.ISDJC.NE.NBPGAU) THEN
  729.               INTERR(1) = IB
  730.               CALL ERREUR(195)
  731.               GOTO 9004
  732.             ENDIF
  733. c
  734. c     remplissage de xmatri
  735. c
  736.             CALL REMPMT(REL,LRE,RE(1,1,ib))
  737. c
  738.  3004     CONTINUE
  739.  
  740.  9004     CONTINUE
  741.           SEGSUP MWRK1,MWRK2
  742.           GOTO 510
  743. c
  744. c_______________________________________________________________________
  745. c
  746. ccccccccccccccccccc element coq3
  747. c_______________________________________________________________________
  748.  27       CONTINUE
  749.           NBBB = NBNN
  750.           SEGINI,MWRK1,MWRK3,MWRK4
  751.  
  752.           DO 3027 IB = 1, NBELEM
  753. c
  754. c     on cherche les coordonnees des noeuds de l element ib
  755. c
  756.             CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  757. c
  758.             CALL ZERO(REL,LRE,LRE)
  759. c
  760. c     on cherche les contraintes
  761. c
  762.             MPTVAL=IVASTR
  763.             DO 5027 ICOMP=1,NSTR
  764.               MELVAL=IVAL(ICOMP)
  765.               IBMN=MIN(IB  ,VELCHE(/2))
  766.               XSTRS(ICOMP)=VELCHE(1,IBMN)
  767.  5027       CONTINUE
  768. c
  769. ccccccc   on calcule k(sigma)
  770. c
  771.             CALL COQ3KS(REL,XSTRS,XE,1.D0,WORK)
  772. c
  773. c     remplissage de xmatri
  774. c
  775.             CALL REMPMT(REL,LRE,RE(1,1,ib))
  776.  
  777.  3027     CONTINUE
  778.  
  779. C 9027     CONTINUE
  780.           SEGSUP,MWRK1,MWRK3,MWRK4
  781.           GOTO 510
  782. c
  783. c_______________________________________________________________________
  784. c
  785. c     element dkt , dst
  786. c_______________________________________________________________________
  787.  28       CONTINUE
  788.           DIM3 = 1.D0
  789.           NBNO = NBNN
  790.           IDI2=IDIM-1
  791.           NBBB=NBNN
  792.           SEGINI MWRK1,MWRK2,MWRK4,MWRK5
  793.           XX(1)=.5D0
  794.           XX(2)=.0D0
  795.           XX(3)=.5D0
  796.           YY(1)=.0D0
  797.           YY(2)=.5D0
  798.           YY(3)=.5D0
  799. C Pour la recuperation de l'epaisseur des elements DKT
  800.           IEPDKT = 0
  801.           IF (MFR.EQ.3 .AND. IFOUR.NE.-2) IEPDKT = IVACAR
  802. c*of 2011/06/22 : Quid de l'epaisseur pour les DST ????? EPAI = 0 ici !!
  803. C Pour la recuperation des axes d'orthotropie des elements DST
  804.           IAODST = 0
  805.           IF (MELE.EQ.93.AND.CMATE.NE.'ISOTROPE') IAODST = IVACAR
  806.  
  807.           DO 3028  IB=1,NBELEM
  808.  
  809. c     on cherche les coordonnees des noeuds de l element ib
  810. c
  811.             CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  812.  
  813.             CALL ZERO(REL,LRE,LRE)
  814.             CALL VPAST(XE,BPSS)
  815. c     bpss    stocke la matrice de passage
  816.             CALL VCORLC (XE,XEL,BPSS)
  817. c
  818. c     boucle sur les points de gauss
  819. c
  820.             DO 4028 IGAU=1,NBPGAU
  821. c
  822. c   recuperation de l'epaisseur (element DKT)
  823.               IF (IEPDKT.NE.0) THEN
  824.                 MPTVAL=IEPDKT
  825.                 MELVAL=IVAL(1)
  826.                 IGMN=MIN(IGAU,VELCHE(/1))
  827.                 IBMN=MIN(IB,VELCHE(/2))
  828.                 EPAI=VELCHE(IGMN,IBMN)
  829.               ELSE
  830.                 EPAI = XZERO
  831.               ENDIF
  832. c
  833.               call DKTSHP(IGAU,XEL,tabw,DJAC)
  834.               call GEOCST(XEL,GEOM)
  835.               call BBGFDK(XX(IGAU),YY(IGAU),GEOM,tabrot)
  836.  
  837.               DO 6028 IC=1,NBNN
  838.                 DO 60281 ID=1,6
  839.                   SHPWRK(ID,IC)=SHPTOT(ID,IC,IGAU)
  840. 60281           CONTINUE
  841.  6028         CONTINUE
  842.  
  843.               CALL DEVOLU(XEL,SHPWRK,MFR,NBNN,IFOUR,NIFOUR,IDI2,DIM3,
  844.      &                    RR,DJAC)
  845.               DJAC=DJAC*POIGAU(IGAU)
  846. c
  847. c     on cherche les contraintes
  848. c
  849.               MPTVAL=IVASTR
  850.               DO 5028 ICOMP=1,NSTRS
  851.                 MELVAL=IVAL(ICOMP)
  852.                 IGMN=MIN(IGAU,VELCHE(/1))
  853.                 IBMN=MIN(IB  ,VELCHE(/2))
  854.                 XSTRS(ICOMP)=VELCHE(IGMN,IBMN)
  855. C                 write(6,*)' xstrs(icomp)',icomp,XSTRS(ICOMP)
  856.  5028         CONTINUE
  857.  
  858. C Recuperation des axes d'orthotropie (element DST)
  859.               IF (IAODST.NE.0) THEN
  860.                 MPTVAL=IAODST
  861.                 MELVAL=IVAL(1)
  862.                 IBMN=MIN(IB  ,VELCHE(/2))
  863.                 IGMN=MIN(IGAU,VELCHE(/1))
  864.                 COSA=VELCHE(IGMN,IBMN)
  865.                 MELVAL=IVAL(2)
  866.                 IBMN=MIN(IB  ,VELCHE(/2))
  867.                 IGMN=MIN(IGAU,VELCHE(/1))
  868.                 SINA=VELCHE(IGMN,IBMN)
  869.                 CC=COSA*COSA
  870.                 SS=SINA*SINA
  871.                 CS=SINA*COSA
  872. C
  873. C        chgt d'axes
  874. C
  875.                 SIG1=CC*XSTRS(1)+SS*XSTRS(2)-2.D0*CS*XSTRS(3)
  876.                 SIG2=CC*XSTRS(2)+SS*XSTRS(1)+2.D0*CS*XSTRS(3)
  877.                 SIG3=CS*(XSTRS(1)-XSTRS(2))+(CC-SS)*XSTRS(3)
  878.                 XSTRS(1)=SIG1
  879.                 XSTRS(2)=SIG2
  880.                 XSTRS(3)=SIG3
  881.               ENDIF
  882. c
  883.               CALL DKTHSH(SHPWRK,tabw,tabrot,DJAC,XSTRS,REL,EPAI)
  884.  4028       CONTINUE
  885.  
  886.             CALL TRANSK(REL,BPSS,LRE,3,1)
  887. c
  888. c     remplissage de xmatri
  889. c
  890.             CALL REMPMT(REL,LRE,RE(1,1,ib))
  891. c
  892.  3028     CONTINUE
  893.  
  894. C 9028     CONTINUE
  895.           SEGSUP,MWRK1,MWRK2,MWRK4,MWRK5
  896.           GOTO 510
  897.  
  898. c_______________________________________________________________________
  899. c
  900. c  element poutre
  901. c_______________________________________________________________________
  902.  29       CONTINUE
  903.           NBBB = NBNN
  904.           SEGINI,MWRK1,MWRK3
  905.  
  906.           DO 3029 IB=1,NBELEM
  907. c
  908. c  on cherche les coordonnees des noeuds de l elementib
  909. c
  910.             CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  911. c
  912. c  il faudrait aussi modifier le vecteur local de la poutre
  913. c
  914. c  mise a zero de la raideur geometrique
  915. c
  916.             CALL ZERO(REL,LRE,LRE)
  917. c
  918. c     rangement des caracteristiques dans work
  919. c
  920.             MPTVAL=IVACAR
  921.             DO 6029 IC=1,NCARR
  922.               WORK(IC)=XZERO
  923.               IF (IVAL(IC).NE.0) THEN
  924.                 MELVAL=IVAL(IC)
  925.                 IBMN=MIN(IB,VELCHE(/2))
  926.                 DO 4029 IGAU=1,NBNN
  927.                   IGMN=MIN(IGAU,VELCHE(/1))
  928.                   IF (IGMN.GT.0.AND.IBMN.GT.0) THEN
  929.                     WORK(IC)=WORK(IC)+VELCHE(IGMN,IBMN)
  930.                   ENDIF
  931.  4029           CONTINUE
  932.                 WORK(IC)=WORK(IC)/NBNN
  933.               ENDIF
  934.  6029       CONTINUE
  935. c
  936. c  cas des tuyaux  -  on calcule les caracteristiques de la poutre
  937. c                     equivalente
  938. c
  939.             IF (MELE.EQ.42) THEN
  940. **            do ic=5,ncarr-1
  941. **            work(ic)=work(ic+1)
  942. **            enddo
  943. **            write (6,*) 'work isous ib',isous,ib
  944. **            write(6,*) (work(ic),ic=1,ncarr)
  945.  
  946.  
  947.  
  948.               CISA=WORK(4)
  949.               VX=WORK(5)
  950.               VY=WORK(6)
  951.               VZ=WORK(7)
  952. **            write(6,*) 'ksigmp vx vy vz',vx,vy,vz,isous,nbelem
  953.               CALL TUYCAR(WORK,CISA,VX,VY,VZ,KERRE,0)
  954.               IF (KERRE.EQ.77) THEN
  955.                 CALL ERREUR(77)
  956.                 GOTO 9029
  957.               ENDIF
  958.             ENDIF
  959. c
  960. c  on cherche les contraintes - on les met dans work
  961. c
  962.             IE = 9
  963.             MPTVAL=IVASTR
  964.             DO 7029 ID=1,2
  965.               ID2=ID
  966.               IF (NBPGAU.EQ.1.AND.ID.EQ.2) ID2=1
  967.               DO 70291 ICOMP=1,NSTR
  968.                 IE = IE+1
  969.                 MELVAL=IVAL(ICOMP)
  970.                 IGMN=MIN(ID2 ,VELCHE(/1))
  971.                 IBMN=MIN(IB  ,VELCHE(/2))
  972.                 WORK(IE)=VELCHE(IGMN,IBMN)
  973. 70291         CONTINUE
  974.  7029       CONTINUE
  975. c
  976. c  on calcule la rigidite geometrique
  977. c
  978.             IF (IFOUR.EQ.-2.OR.IFOUR.EQ.-1.OR.IFOUR.EQ.-3) THEN
  979.               CALL POUKS2(REL,LRE,WORK(10),WORK,XE,WORK(22),KERRE)
  980.             ELSE
  981.               CALL POUKSG(REL,LRE,WORK(10),WORK,XE,WORK(22),KERRE)
  982.             ENDIF
  983.  
  984.             IF (KERRE.NE.0) THEN
  985.               INTERR(1)=ISOUS
  986.               INTERR(2)=IB
  987.               CALL ERREUR(128)
  988.               GOTO 9029
  989.             ENDIF
  990. c
  991. c  remplissage de xmatri
  992. c
  993.             CALL REMPMT(REL,LRE,RE(1,1,ib))
  994. C
  995.  3029     CONTINUE
  996. c
  997.  9029     CONTINUE
  998.           SEGSUP,MWRK1,MWRK3
  999. c
  1000.           GOTO 510
  1001. c_______________________________________________________________________
  1002. c
  1003. c  elements coq8 et coq6
  1004. c_______________________________________________________________________
  1005.  41       CONTINUE
  1006.           NBBB=NBNN
  1007.           LRI =NBNN*5
  1008.           SEGINI,MWRK1,MWRK3
  1009. c
  1010.           MINTE1 = IPMIN1
  1011.           SEGACT,MINTE1
  1012. c
  1013.           DO 3041 IB=1,NBELEM
  1014. c
  1015. c  on cherche les coordonnees des noeuds de l elementib
  1016. c
  1017.             CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1018.  
  1019.             CALL ZERO(REL,LRE,LRE)
  1020. c
  1021. c  on cherche les caracteristiques de l element ib
  1022. c
  1023.             MPTVAL = IVACAR
  1024.             MELVAL = IVAL(1)
  1025.             IF (MELVAL.NE.0) THEN
  1026.               IBMN = MIN(IB  ,VELCHE(/2))
  1027.               DO IGAU = 1, NBNN
  1028.                 IGMN = MIN(IGAU,VELCHE(/1))
  1029.                 WORK(IGAU) = VELCHE(IGMN,IBMN)
  1030.               ENDDO
  1031.             ELSE
  1032.               DO IGAU = 1, NBNN
  1033.                 WORK(IGAU)=XZERO
  1034.               ENDDO
  1035.             ENDIF
  1036. c
  1037. c  on cherche les contraintes - on les met dans work
  1038. c
  1039.             IE = 9
  1040.             MPTVAL=IVASTR
  1041.             DO 7041 IGAU=1,NBPGAU
  1042.               DO 7042 ICOMP=1,NSTRS
  1043.                 MELVAL=IVAL(ICOMP)
  1044.                 IGMN=MIN(IGAU,VELCHE(/1))
  1045.                 IBMN=MIN(IB  ,VELCHE(/2))
  1046.                 WORK(IE)=VELCHE(IGMN,IBMN)
  1047.                 IE=IE+1
  1048.  7042         CONTINUE
  1049.  7041       CONTINUE
  1050. c
  1051. c  on calcule la rigidite geometrique
  1052. c
  1053.             CALL COQ8KS(REL,XE,SHPTOT,MINTE1.SHPTOT,
  1054.      &                  NBPGAU,POIGAU,DZEGAU,
  1055.      &                  WORK(1),WORK(9),NBNN,LRE,LRI,WORK(51))
  1056. c
  1057. c  remplissage de xmatri
  1058. c
  1059.             CALL REMPMT(REL,LRE,RE(1,1,ib))
  1060. C
  1061.  3041     CONTINUE
  1062.  
  1063. C 9041     CONTINUE
  1064.           SEGSUP,MWRK1,MWRK3
  1065.           GO TO 510
  1066. c_______________________________________________________________________
  1067. c
  1068. c     tuyau fissure
  1069. c_______________________________________________________________________
  1070.  43       CONTINUE
  1071. c     ksigma n a pas de sens evident pour cet element
  1072. c     on cree une matrice nulle
  1073. c          DO 3043 IB=1,NBELEM
  1074. c            do 4043 ic=1,lval
  1075. c              re(ic,ic,ib)=XZERO
  1076. c 4043       continue
  1077. c 3043     CONTINUE
  1078.           GOTO 510
  1079. c
  1080. c_______________________________________________________________________
  1081. c
  1082. c     element coq2
  1083. c_______________________________________________________________________
  1084. c
  1085.  44       CONTINUE
  1086.           DIM3=1.D0
  1087.           NBBB=NBNN
  1088.           SEGINI MWRK1,MWRK3,MWRK4
  1089. c
  1090.           DO 3044 IB=1,NBELEM
  1091. c
  1092. c     on cherche les coordonnees des noeuds de  l element ib
  1093. c
  1094.             CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1095.  
  1096.             CALL ZERO(REL,LRE,LRE)
  1097. c
  1098. c   recuperation de l'epaisseur
  1099. c
  1100.             IF (IFOUR.EQ.-2.AND.IPCHE2.NE.0) THEN
  1101.               MPTVAL=IVACAR
  1102.               MELVAL=IVAL(1)
  1103.               IF (MELVAL.NE.0) THEN
  1104.                 IBMN=MIN(IB,VELCHE(/2))
  1105.                 DIM3=VELCHE(1,IBMN)
  1106.               ELSE
  1107.                 DIM3 = 1.D0
  1108.               ENDIF
  1109.             ENDIF
  1110. c
  1111. c     on cherche les contraintes on les met dans work...
  1112. c
  1113.             JC = 0
  1114.             MPTVAL=IVASTR
  1115.             DO 5044 IGAU=1,NBPGAU
  1116.               DO 5045 ICOMP=1,NSTRS
  1117.                 MELVAL=IVAL(ICOMP)
  1118.                 IGMN=MIN(IGAU,VELCHE(/1))
  1119.                 IBMN=MIN(IB  ,VELCHE(/2))
  1120.                 JC=JC+1
  1121.                 WORK(JC)=VELCHE(IGMN,IBMN)
  1122.  5045         CONTINUE
  1123.  5044       CONTINUE
  1124. c
  1125. c     appel a coque2 ksigma...
  1126. c
  1127.             AN=NHRM
  1128.             CALL CQ2KSG(XE,1.D0,DIM3,IFOUR,AN,NBPGAU,WORK(1),WORK(19),
  1129.      1                  WORK(22),QSIGAU,POIGAU,WORK(25),WORK(30),
  1130.      2                  WORK(35),WORK(42),WORK(49),WORK(113),WORK(177),
  1131.      3                  WORK(241),WORK(305),LRE,REL)
  1132. c
  1133. c     remplissage de xmatri
  1134. c
  1135.             CALL REMPMT(REL,LRE,RE(1,1,ib))
  1136. c
  1137.  3044     CONTINUE
  1138. c
  1139. C 9044     CONTINUE
  1140.           SEGSUP,MWRK1,MWRK3,MWRK4
  1141.           GOTO 510
  1142. c
  1143. c_______________________________________________________________________
  1144. c
  1145. c  elements barre et cercle (et TIMO)
  1146. c_______________________________________________________________________
  1147.  46       CONTINUE
  1148. C Cas particulier :
  1149.           IF (MELE.EQ.95.AND.IFOUR.NE.0.AND.IFOUR.NE.1) GOTO 99
  1150. C
  1151.           NBBB = NBNN
  1152.           SEGINI MWRK1,MWRK3
  1153. c
  1154.           DO 3046 IB=1,NBELEM
  1155. c
  1156. c  on cherche les coordonnees des noeuds de l elementib
  1157. c
  1158.             CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1159. c
  1160. c  mise a zero de la raideur geometrique
  1161. c
  1162.             CALL ZERO(REL,LRE,LRE)
  1163. c
  1164. c  on cherche l'effort
  1165. c
  1166.             MPTVAL=IVASTR
  1167.             MELVAL=IVAL(1)
  1168.             NBPTEL=VELCHE(/1)
  1169.             IBMN=MIN(IB,VELCHE(/2))
  1170. c
  1171.             IF (NBPTEL.EQ.1) THEN
  1172.               EFFORT=VELCHE(1,IBMN)
  1173.             ELSE IF (NBPTEL.EQ.2) THEN
  1174.               EFF1=VELCHE(1,IBMN)
  1175.               EFF2=VELCHE(2,IBMN)
  1176.               EFFORT=0.5D0*(EFF1+EFF2)
  1177.             ENDIF
  1178. c
  1179. c  on calcule la rigidite geometrique
  1180. c
  1181.             IF (MELE.EQ.46.or.MELE.eq.84)
  1182.      &         CALL BARKSG(REL,LRE,EFFORT,XE,KERRE)
  1183.             IF (MELE.EQ.95) CALL CERKSG(REL,LRE,EFFORT,XE,KERRE)
  1184.             IF (KERRE.NE.0) THEN
  1185.               INTERR(1)=ISOUS
  1186.               INTERR(2)=IB
  1187.               CALL ERREUR(128)
  1188.               GO TO 9046
  1189.             ENDIF
  1190. c
  1191. c  remplissage de xmatri
  1192. c
  1193. c           cas particulier TIMO : on saute les ddls de rotation
  1194.             IF (MELE.EQ.84) THEN
  1195.               NCOMPU=NCOMP/2
  1196.               ii=0
  1197.               iii=0
  1198.               DO 841 INOEUD=1,NBNNS
  1199.               DO 842 ICOMP=1,NCOMP
  1200.                 ii=ii+1
  1201.                 if(ii.gt.NCOMPU) goto 842
  1202.                 iii=iii+1
  1203.                 jj=0
  1204.                 jjj=0
  1205.                 DO 843 JNOEUD=1,NBNNS
  1206.                 DO 844 JCOMP=1,NCOMP
  1207.                   jj=jj+1
  1208.                   if(jj.gt.ii) goto 842
  1209.                   if(jj.gt.NCOMPU) goto 844
  1210.                   jjj=jjj+1
  1211.                   RE(ii,jj,ib)=REL(iii,jjj)
  1212.                   RE(jj,ii,ib)=REL(iii,jjj)
  1213.  844            CONTINUE
  1214.  843            CONTINUE
  1215.  842          CONTINUE
  1216.  841          CONTINUE
  1217.             ELSE
  1218.               CALL REMPMT(REL,LRE,RE(1,1,ib))
  1219.             ENDIF
  1220. C
  1221.  3046     CONTINUE
  1222.  
  1223.  9046     CONTINUE
  1224.           SEGSUP,MWRK1,MWRK3
  1225.           GOTO 510
  1226.  
  1227. c_______________________________________________________________________
  1228. c
  1229. c     element coq4
  1230. c_______________________________________________________________________
  1231.  49       CONTINUE
  1232.           NBBB=NBNN
  1233.           NBNO=NBNN
  1234.           SEGINI,MWRK1,MWRK2,MWRK4
  1235.  
  1236.           DO 3049 IB=1,NBELEM
  1237. c
  1238. c     on cherche les coordonnees des noeuds de  l element ib
  1239. c
  1240.             CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1241.  
  1242.             CALL ZERO(REL,LRE,LRE)
  1243.             CALL CQ4LOC(XE,XEL,BPSS,IRRT,0)
  1244. C
  1245. C   attention  :  rien de prevu en cas d'excentrement
  1246. C
  1247.             CALL BCOQ4(5,XEL,SHPTOT,SHPWRK,BGENE,DJAC,XZERO,0,IRRT,0)
  1248. c
  1249.             MPTVAL=IVASTR
  1250.             DO 5049 ICOMP=1,NSTRS
  1251.               MELVAL=IVAL(ICOMP)
  1252.               IGMN=MIN(5,VELCHE(/1))
  1253.               IBMN=MIN(IB  ,VELCHE(/2))
  1254.               XSTRS(ICOMP)=VELCHE(IGMN,IBMN)
  1255.  5049       CONTINUE
  1256. c
  1257.             CALL CQ4KSG(DJAC,XSTRS,SHPWRK, REL)
  1258.             CALL TRANSK(REL,BPSS,LRE,4,0)
  1259. c
  1260. c     remplissage de xmatri
  1261. c
  1262.             CALL REMPMT(REL,LRE,RE(1,1,ib))
  1263. C
  1264.  3049     CONTINUE
  1265.  
  1266. C 9049     CONTINUE
  1267.           SEGSUP,MWRK1,MWRK2,MWRK4
  1268.           GOTO 510
  1269.  
  1270. c_______________________________________________________________________
  1271. c
  1272. c     element cof3
  1273. c_______________________________________________________________________
  1274. c
  1275.  51       CONTINUE
  1276. c
  1277.           NBBB=NBNN
  1278.           SEGINI,MWRK1,MWRK3,MWRK4
  1279. c
  1280.           CALL ERREUR(19)
  1281.           GOTO 9051
  1282.  
  1283. C          DO 3051 IB=1,NBELEM
  1284. Cc
  1285. Cc     on cherche les coordonnees des noeuds de  l element ib
  1286. Cc
  1287. C            CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1288. C
  1289. C            CALL ZERO(REL,LRE,LRE)
  1290. C
  1291. C            MPTVAL=IVACAR
  1292. C            MELVAL=IVAL(1)
  1293. C            IF (MELVAL.NE.0) THEN
  1294. C              IBMN=MIN(IB  ,VELCHE(/2))
  1295. C              EPAI=VELCHE(1,IBMN)
  1296. C            ELSE
  1297. C              EPAI=XZERO
  1298. C            ENDIF
  1299. Cc
  1300. Cc     on cherche les contraintes on les met dans work...
  1301. Cc
  1302. C            JC=0
  1303. C            MPTVAL=IVASTR
  1304. C            DO 5051 IGAU=1,NBPGAU
  1305. C              DO 5051 ICOMP=1,NSTRS
  1306. C                MELVAL=IVAL(ICOMP)
  1307. C                IGMN=MIN(IGAU,VELCHE(/1))
  1308. C                IBMN=MIN(IB  ,VELCHE(/2))
  1309. C                JC=JC+1
  1310. C                WORK(JC)=VELCHE(IGMN,IBMN)
  1311. C 5051       CONTINUE
  1312. Cc
  1313. Cc     appel a coque2 ksigma...
  1314. Cc
  1315. C            AN=NHRM
  1316. CC             call cq3ksg(xe,epai,an,nbpgau,work(1),work(19),work(22),
  1317. CC     1                   work(25),work(30),work(35),work(42),work(49),
  1318. CC     2                   work(113),work(177),work(241),work(305),rel)
  1319. Cc
  1320. Cc     remplissage de xmatri
  1321. Cc
  1322. C            CALL REMPMT(REL,LRE,RE(1,1,ib))
  1323. Cc
  1324. C 3051     CONTINUE
  1325.  
  1326.  9051     CONTINUE
  1327.           SEGSUP MWRK1,MWRK3,MWRK4
  1328.           GOTO 510
  1329.  
  1330. c_______________________________________________________________________
  1331. c
  1332. c     element shb8
  1333. c_______________________________________________________________________
  1334.  260      CONTINUE
  1335.           NBBB=NBNN
  1336.           SEGINI,MWRK1,MWRK7
  1337. C           write(6,*) ' nbnn nbpgau nstrs lre' , NBNN,nbpgau,nstrs,lre
  1338.  
  1339.           DO 3260 IB=1,NBELEM
  1340. c
  1341. c     on cherche les coordonnees des noeuds de  l element ib
  1342. c
  1343.             CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  1344.  
  1345.             MPTVAL=IVASTR
  1346.             IE=0
  1347.             do 3268 igau=1,nbpgau
  1348.               DO 3269 ICOMP=1,NSTRS
  1349.                 iE=IE+1
  1350.                 MELVAL=IVAL(ICOMP)
  1351.                 IGMN=MIN(IGAU,VELCHE(/1))
  1352.                 IBMN=MIN(IB  ,VELCHE(/2))
  1353.                 work1(ie)=VELCHE(IGMN,IBMN)
  1354. C                 write(6,*)' xstrs(icomp)',icomp,XSTRS(ICOMP)
  1355.  3269         CONTINUE
  1356.  3268       CONTINUE
  1357.             propel(1)=0.
  1358.             call shb8 (9,xe,D,propel,work1,rel,out)
  1359.  
  1360. C      remplissage de xmatri
  1361.             CALL REMPMT(REL,LRE,RE(1,1,ib))
  1362. C
  1363.  3260     CONTINUE
  1364. c
  1365. C 9260     CONTINUE
  1366.           SEGSUP,MWRK1,MWRK7
  1367.           GOTO 510
  1368.  
  1369. c_______________________________________________________________________
  1370. c
  1371. c     desactivation des segments propres a la zone geometrique isous
  1372. c_______________________________________________________________________
  1373. c
  1374.  510      CONTINUE
  1375.  597      CONTINUE
  1376.           IF (ISUP1L.EQ.1 .OR. nblprt.GT.1) THEN
  1377.             CALL DTMVAL(IVASTR,3)
  1378.           ELSE
  1379.             CALL DTMVAL(IVASTR,1)
  1380.           ENDIF
  1381.           IF (ISUP2L.EQ.1 .OR. nblprt.GT.1) THEN
  1382.             CALL DTMVAL(IVACAR,3)
  1383.           ELSE
  1384.             CALL DTMVAL(IVACAR,1)
  1385.           ENDIF
  1386.           xmatri = ipmatr
  1387.           SEGDES,xmatri
  1388.  
  1389. C- Sortie prematuree en cas d'erreur
  1390.           IF (IERR.NE.0) GOTO 598
  1391.  
  1392. C- Stockage de la matrice
  1393.           nrigel=irigel(/2) +1
  1394.           segadj,mrigid
  1395. C           jrige = NRIGE0 + irige
  1396.           jrige=nrigel
  1397.           COERIG(jrige)   = 1.
  1398.           IRIGEL(1,jrige) = ipmail
  1399.           IRIGEL(2,jrige) = 0
  1400.           IRIGEL(3,jrige) = ipdesc
  1401.           IRIGEL(4,jrige) = ipmatr
  1402.           IRIGEL(5,jrige) = NIFOUR
  1403.           IRIGEL(6,jrige) = 0
  1404.           IRIGEL(7,jrige) = 0
  1405.           IRIGEL(8,jrige) = 0
  1406.  
  1407.         ENDDO
  1408. C- Fin de la boucle sur les partitions
  1409. C
  1410.  598    CONTINUE
  1411.         IF (MOSTRS.NE.0) THEN
  1412.           nomid = MOSTRS
  1413.           IF (lsupco) SEGSUP,nomid
  1414.           notype = MOTYPS
  1415.           SEGSUP,notype
  1416.         ENDIF
  1417.         IF (MOCARA.NE.0) THEN
  1418.           NOMID = MOCARA
  1419.           SEGSUP,NOMID
  1420.           notype = MOTYPC
  1421.           SEGSUP,notype
  1422.         ENDIF
  1423. C
  1424.         NOMID=MODEPL
  1425.         IF (lsupde) SEGSUP,NOMID
  1426.         NOMID = MOFORC
  1427.         IF (lsupfo) SEGSUP,NOMID
  1428.  
  1429.  599    CONTINUE
  1430.  
  1431.         IF (IERR.NE.0) GOTO 999
  1432. C
  1433.  500  CONTINUE
  1434. C* Fin de la boucle sur les modeles elementaires
  1435.  
  1436.  999  CONTINUE
  1437.       IF (IERR.NE.0) THEN
  1438.         ktrace = -1
  1439. **  En situation d'erreur, on laisse le menage faire son travail
  1440. **      CALL DERIGI(MRIGID,ktrace,msorse)
  1441. **      SEGSUP,MRIGID
  1442.         IPRIGG = 0
  1443.       ELSE
  1444.         if(irigel(/2).eq.0) then
  1445.           call erreur (86)
  1446.           return
  1447.         endif
  1448.         SEGDES,MRIGID
  1449.         IPRIGG = MRIGID
  1450.       ENDIF
  1451.  
  1452.       END
  1453.  
  1454.  
  1455.  
  1456.  
  1457.  
  1458.  
  1459.  
  1460.  
  1461.  
  1462.  

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