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

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