Télécharger ksigmp.eso

Retour à la liste

Numérotation des lignes :

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

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