Télécharger divc.eso

Retour à la liste

Numérotation des lignes :

  1. C DIVC SOURCE PV 16/11/17 21:59:02 9180
  2. SUBROUTINE DIVC(IKAS,MTABX)
  3. IMPLICIT INTEGER(I-N)
  4. IMPLICIT REAL*8 (A-H,O-Z)
  5. C************************************************************************
  6. C Operateur KMAC
  7. C
  8. C OBJET : Cree un objet de type MATRIK
  9. C
  10. C
  11. C IKAS=1 KMAB calcul de B (DIV U)
  12. C IKAS=2 KMBT calcul de Bt uniquement (GRAD P)
  13. C IKAS=3 KBBT calcul de B assemblage pour B et Bt
  14. C
  15. C***********************************************************************
  16. -INC CCOPTIO
  17. -INC CCGEOME
  18. -INC CCREEL
  19. -INC SMCOORD
  20. -INC SIZFFB
  21. -INC SMCHPOI
  22. POINTEUR IZCH2.MCHPOI,IZCCH2.MPOVAL
  23. POINTEUR IZDV.MCHPOI,IZDDV.MPOVAL,IZTU1.MPOVAL,TETAN.MPOVAL
  24. POINTEUR IZTG1.MCHPOI,IZTGG1.MPOVAL,IZBETA.MPOVAL,VISCO.MPOVAL
  25. POINTEUR IPM.IZAFM
  26. SEGMENT IMATRS
  27. ENDSEGMENT
  28. POINTEUR IPMS.IZAFM,IPS1.IZAFM,IPS2.IZAFM,IPS3.IZAFM
  29.  
  30. PARAMETER (LRV=64)
  31. DIMENSION VISC(3)
  32.  
  33. -INC SMLENTI
  34. POINTEUR IZIPAD.MLENTI,MLENTI1.MLENTI,MLENTI2.MLENTI
  35. -INC SMLMOTS
  36. POINTEUR LINCO.MLMOTS
  37. -INC SMELEME
  38. POINTEUR MELEMZ.MELEME,MELEMB.MELEME,MELSTB.MELEME
  39. POINTEUR MELEM1.MELEME,MELES1.MELEME,MCTREI.MELEME
  40. POINTEUR IGEOM.MELEME,MELEMA.MELEME
  41. POINTEUR IZLEMC.MELEME,MELEMS.MELEME,MELEMC.MELEME
  42. POINTEUR MELEMI.MELEME,MELEMP.MELEME
  43. POINTEUR MACRO1.MELEME
  44.  
  45. CHARACTER*8 TYPE,TYPC,NOMZ,NOMP,NOMD,NOM0
  46. CHARACTER*8 NOMPP,NOMDD
  47. CHARACTER*4 NOM
  48. INTEGER IPAD,IPAD2,IK
  49. REAL*8 KAUX,TETA1
  50. DIMENSION IXV(5)
  51. C
  52. C*************************************************************************
  53. CKMIC
  54. C write(6,*)' Operateur KMIC MTABX=',MTABX
  55. C
  56. C- Récupération de la table EQEX (pointeur MTAB1)
  57. C
  58. CALL LEKTAB(MTABX,'EQEX',MTAB1)
  59. IF(MTAB1.EQ.0)THEN
  60. C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24
  61. MOTERR( 1: 8) = ' EQEX '
  62. MOTERR( 9:16) = ' EQEX '
  63. MOTERR(17:24) = ' KIZX '
  64. CALL ERREUR(786)
  65. RETURN
  66. ENDIF
  67.  
  68. CALL ACME(MTAB1,'DISCPRES',IDISCP)
  69.  
  70. CALL ACME(MTAB1,'NAVISTOK',NASTOK)
  71. IF(NASTOK.EQ.0)THEN
  72. CALL ZKMIC(IKAS,MTABX,MTAB1)
  73. RETURN
  74. ENDIF
  75. C
  76. C- Récupération de la table INCO (pointeur KINC)
  77. C
  78. CALL LEKTAB(MTAB1,'INCO',KINC)
  79. IF(KINC.EQ.0)THEN
  80. C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24
  81. MOTERR( 1: 8) = ' INCO '
  82. MOTERR( 9:16) = ' INCO '
  83. MOTERR(17:24) = ' EQEX '
  84. CALL ERREUR(786)
  85. RETURN
  86. ENDIF
  87.  
  88. C*************************************************************************
  89. C OPTIONS
  90. C KIMPL = 0 -> EXPL 1 -> IMPL 2 -> CN
  91. C KFORM = 0 -> EFSI 1 -> EF 2 -> VF 3 -> EFMC
  92. C KPRE=3 pression P0 KPRE=4 pression P1 KPRE=2 cas macro 1ère génération
  93. C KPRE=5 pression MSOMMET
  94.  
  95. IAXI=0
  96. IK=0
  97. IF(IFOMOD.EQ.0)IAXI=2
  98. C
  99. C- Récupération de la table des options KOPT (pointeur KOPTI)
  100. C
  101. CALL LEKTAB(MTABX,'KOPT',KOPTI)
  102. IF (KOPTI.EQ.0) THEN
  103. C Indice %m1:8 : Indice %m9:16 non trouvé dans la table %m17:24
  104. MOTERR( 1: 8) = ' KOPT '
  105. MOTERR( 9:16) = ' KOPT '
  106. MOTERR(17:24) = ' KIZX '
  107. CALL ERREUR(786)
  108. RETURN
  109. ENDIF
  110.  
  111. CALL ACME(KOPTI,'KIMPL',KIMPL)
  112. CALL ACME(KOPTI,'KPOIN',KPRE)
  113. C write(6,*)'KBBT KPRE=',KPRE
  114. CALL ACMF(KOPTI,'AIMPL',TETA1)
  115. CALL ACME(KOPTI,'KFORM',KFORM)
  116. C CALL ACME(KOPTI,'IKOMP',IKOMP)
  117.  
  118. IF (IERR.NE.0) RETURN
  119. C write(6,*)' Apres les options '
  120. C*************************************************************************
  121. C
  122. C- Récupération de la table DOMAINE associée au domaine local
  123. C
  124. CALL ACMM(MTABX,'NOMZONE',NOMZ)
  125. TYPE=' '
  126. CALL ACMO(MTABX,'DOMZ',TYPE,MMODEL)
  127. IF(TYPE.NE.'MMODEL')THEN
  128. C On attend un des objets : %m1:8 %m9:16 %m17:24 %m25:32 %m33:40
  129. MOTERR( 1: 8) = ' MMODEL '
  130. MOTERR( 8:16) = ' MMODEL '
  131. MOTERR(17:24) = ' MMODEL '
  132. MOTERR(25:32) = ' MMODEL '
  133. MOTERR(33:40) = ' MMODEL '
  134. CALL ERREUR(471)
  135. RETURN
  136. ENDIF
  137.  
  138. CALL LEKMOD(MMODEL,MTABZ,INEFMD)
  139.  
  140. CALL LEKTAB(MTABZ,'MAILLAGE',MELEME)
  141. CALL LEKTAB(MTABZ,'SOMMET',MELEMS)
  142. CALL LEKTAB(MTABZ,'CENTRE',MELEMC)
  143. C?? CALL LEKTAB(MTABZ,'MACRO',MACRO)
  144. MACRO1=0
  145. IF(INEFMD.EQ.2)CALL LEKTAB(MTABZ,'MACRO1',MACRO1)
  146. C?? IF(INEFMD.EQ.3)CALL LEKTAB(MTABZ,'QUADRATI',MQUAD)
  147. IF (IERR.NE.0) RETURN
  148. IF(INEFMD.EQ.4.AND.KPRE.NE.5)THEN
  149. C% Données incompatibles
  150. CALL ERREUR(21)
  151. RETURN
  152. ENDIF
  153.  
  154. MELEMI=MELEME
  155. IF(MACRO1.NE.0.AND.KPRE.NE.2)THEN
  156. MELEMI=MACRO1
  157. ENDIF
  158.  
  159. IF(KPRE.EQ.2.AND.INEFMD.EQ.3)KPRE=3
  160. IF(INEFMD.EQ.1.AND.KPRE.NE.5)KPRE=2
  161.  
  162. IF(KPRE.EQ.3)THEN
  163. CALL LEKTAB(MTABZ,'CENTREP0',MELEMC)
  164. MELEMP=MELEMC
  165. ELSEIF(KPRE.EQ.4)THEN
  166. CALL LEKTAB(MTABZ,'CENTREP1',MELEMC)
  167. CALL LEKTAB(MTABZ,'ELTP1NC ',MELEMP)
  168. ELSEIF(KPRE.EQ.2)THEN
  169. CALL LEKTAB(MTABZ,'CENTRE',MELEMC)
  170. MELEMP=MELEMC
  171. ELSEIF(KPRE.EQ.5)THEN
  172. CALL LEKTAB(MTABZ,'MSOMMET',MELEMC)
  173. CALL LEKTAB(MTABZ,'MMAIL ',MELEMP)
  174. ENDIF
  175. C write(6,*)' KMIC KPRE=',kpre
  176.  
  177. IF(IDISCP.EQ.0)THEN
  178. CALL ECME(MTAB1,'DISCPRES',KPRE)
  179. ELSEIF(IDISCP.NE.KPRE)THEN
  180. C% Données incompatibles
  181. CALL ERREUR(21)
  182. RETURN
  183. ENDIF
  184.  
  185. IKOMP=1
  186. IF(KPRE.EQ.5)IKOMP=0
  187.  
  188. C*************************************************************************
  189. C VERIFICATIONS SUR LES INCONNUES
  190.  
  191. C write(6,*)' Verification sur les inconnues '
  192. TYPE='LISTMOTS'
  193. CALL ACMO(MTABX,'LISTINCO',TYPE,LINCO)
  194. IF (IERR.NE.0) RETURN
  195. SEGACT LINCO
  196. NBINC=LINCO.MOTS(/2)
  197. IF(NBINC.NE.2)THEN
  198. WRITE(6,*)'Operateur KMAB '
  199. WRITE(6,*)'Nombre d''inconnues incorrecte : ',NBINC,' On attend 2'
  200. C Indice %m1:8 : contient plus de %i1 %m9:16
  201. MOTERR( 1:8) = 'LISTINCO'
  202. INTERR(1) = 2
  203. MOTERR(9:16) = ' MOTS '
  204. CALL ERREUR(799)
  205. RETURN
  206. ENDIF
  207.  
  208. C On recupere PHI n et TETA n pour Cranck-Nicholson
  209. NOMP=LINCO.MOTS(1)
  210. TYPE=' '
  211. CALL ACMO(KINC,NOMP,TYPE,MCHPOI)
  212. IF(TYPE.NE.'CHPOINT ')THEN
  213. WRITE(6,*)' Opérateur KMAB :'
  214. WRITE(6,*)' L objet CHPOINT ',NOMP,
  215. & ' n existe pas dans la table'
  216. C Indice %m1:8 : ne contient pas un objet de type %m9:16
  217. MOTERR( 1: 8) = 'INC '//NOMP
  218. MOTERR( 9:16) = 'CHPOINT '
  219. CALL ERREUR(800)
  220. RETURN
  221. ELSE
  222. CALL LICHT(MCHPOI,IZTU1,TYPC,IGEOM0)
  223. ENDIF
  224. C*************************************************************************
  225. C Le domaine de definition est donne par le SPG de la premiere inconnue
  226. C Les inconnues suivantes devront posseder ce meme pointeur
  227. C On verifie que les points de la zone sont tous inclus dans ce SPG
  228. C Inconnue Primale
  229.  
  230. C write(6,*)' Verification inconnue primale '
  231. CALL KRIPAD(IGEOM0,MLENTI)
  232. IF(IKAS.EQ.1.OR.IKAS.EQ.3)THEN
  233. MELEMK=MELEMS
  234. ELSE
  235. MELEMK=MELEMC
  236. ENDIF
  237.  
  238. C write(6,*)' MELEMK=',melemk
  239. CALL VERPAD(MLENTI,MELEMK,IRET)
  240. IF(IRET.NE.0)THEN
  241. WRITE(6,*)' Opérateur KMAB '
  242. WRITE(6,*)' La zone ',NOMZ,' n''est pas incluse dans le domaine'
  243. & , ' de définition de l''inconnue ',NOMP
  244. WRITE(6,*)' MELEMK=',melemk,' IGEOM0=',IGEOM0
  245. C Indice %m1:8 : L'objet %m9:16 n'a pas le bon support géométrique
  246. MOTERR(1: 8) = 'INC '//NOMP
  247. MOTERR(9:16) = 'CHPOINT '
  248. CALL ERREUR(788)
  249. IPAS=0
  250. RETURN
  251. ENDIF
  252.  
  253. SEGSUP MLENTI
  254.  
  255. C*************************************************************************
  256.  
  257. NOMD=LINCO.MOTS(2)
  258. TYPE=' '
  259. CALL ACMO(KINC,NOMD,TYPE,MCHPOI)
  260. IF(TYPE.NE.'CHPOINT ')THEN
  261. WRITE(6,*)' Opérateur KMAB :'
  262. WRITE(6,*)' L objet CHPOINT ',NOMD,
  263. & ' n existe pas dans la table'
  264. C Indice %m1:8 : ne contient pas un objet de type %m9:16
  265. MOTERR( 1: 8) = 'INC '//NOMD
  266. MOTERR( 9:16) = 'CHPOINT '
  267. CALL ERREUR(800)
  268. RETURN
  269. ELSE
  270. CALL LICHT(MCHPOI,TETAN,TYPC,IGEOM0)
  271. ENDIF
  272.  
  273. NC=TETAN.VPOCHA(/2)
  274. C*************************************************************************
  275. C Le domaine de definition est donne par le SPG de la premiere inconnue
  276. C Les inconnues suivantes devront posseder ce meme pointeur
  277. C On verifie que les points de la zone sont tous inclus dans ce SPG
  278. C Inconnue Duale
  279.  
  280. C write(6,*)' IGEOM0=',igeom0
  281. CALL KRIPAD(IGEOM0,MLENTI)
  282. IF(IKAS.EQ.1.OR.IKAS.EQ.3)THEN
  283. MELEMK=MELEMC
  284. ELSE
  285. MELEMK=MELEMS
  286. ENDIF
  287.  
  288. C write(6,*)' Verification inconnue duale ',MELEMK
  289. CALL VERPAD(MLENTI,MELEMK,IRET)
  290. IF(IRET.NE.0)THEN
  291. WRITE(6,*)' Opérateur KMAB '
  292. WRITE(6,*)' La zone ',NOMZ,' n''est pas incluse dans le domaine'
  293. & ,' de définition de l''inconnue ',NOMD
  294. WRITE(6,*)' MELEMK=',melemk,' IGEOM0=',IGEOM0
  295. C Indice %m1:8 : L'objet %m9:16 n'a pas le bon support géométrique
  296. MOTERR(1: 8) = 'INC '//NOMD
  297. MOTERR(9:16) = 'CHPOINT '
  298. CALL ERREUR(788)
  299. IPAS=0
  300. RETURN
  301. ENDIF
  302.  
  303. SEGSUP MLENTI
  304.  
  305. C*************************************************************************
  306. C Lecture du ou des coefficients
  307. C Type du coefficient :
  308. C IK1=0 CHPOINT IK1=1 scalaire IK1=2 vecteur
  309.  
  310. C write(6,*)' Verification sur les coefficients '
  311. CALL ACME(MTABX,'IARG',IARG)
  312. C write(6,*)' IARG=',iarg
  313.  
  314. IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.2)THEN
  315. TYPE=' '
  316. CALL ACMO(MTABZ,'MELSTB',TYPE,MELSTB)
  317. SEGACT MELSTB
  318. IF(IDIM.EQ.2)NBELEM=MELSTB.NUM(/2)/4
  319. IF(IDIM.EQ.3)NBELEM=MELSTB.NUM(/2)/8
  320. NBNN=MELSTB.NUM(/1)
  321. NBSOUS=0
  322. NBREF=0
  323. SEGINI MELEMA
  324. MELEMA.ITYPEL=MELSTB.ITYPEL
  325.  
  326. NKPE=4
  327. IF(IDIM.EQ.3)NKPE=8
  328. do 4878 k=1,nbelem
  329. mi=(k-1)*NKPE+1
  330. do 4879 i=1,nbnn
  331. MELEMA.num(i,k)=melstb.num(i,mi)
  332. 4879 continue
  333. C write(6,*)k,(MELEMA.num(i,k),i=1,nbnn)
  334. 4878 continue
  335.  
  336. TYPE=' '
  337. CALL ACMO(MTABZ,'MCHPOC',TYPE,MCHPOC)
  338. TYPE=' '
  339. CALL ACMO(MTABZ,'CENTRE',TYPE,MCTREI)
  340. ENDIF
  341.  
  342.  
  343. C 1er COEF
  344.  
  345. IXV(1)=MELEMC
  346. IXV(2)=1
  347. IXV(3)=0
  348. IXV(4)=MELEMS
  349. IXV(5)=-MELEMS
  350. IRET=5
  351. CALL LEKCOF('Opérateur KMAB :',
  352. & MTABX,KINC,1,IXV,IZTG1,IZTGG1,NPT1,NC1,IK1,IRET)
  353. IF(IRET.EQ.0)RETURN
  354. IF(IK1.GE.4)CALL KRIPAD(MELEMS,MLENTI)
  355.  
  356.  
  357. BETA0=1.D-1
  358. IF(IARG.EQ.2)THEN
  359. C 2ème COEF
  360. IXV(1)=0
  361. IXV(2)=1
  362. IXV(3)=0
  363. CALL LEKCOF('Opérateur KMAB :',
  364. & MTABX,KINC,2,IXV,IZTG2,IZBETA,NPT2,NC2,IK2,IRET)
  365. IF(IRET.EQ.0)RETURN
  366. BETA0=IZBETA.VPOCHA(1,1)
  367. C? IF(IZBETA.VPOCHA(1,1).LT.0.D0)THEN
  368. C% Données incompatibles
  369. C? CALL ERREUR(21)
  370. C? RETURN
  371. C? ENDIF
  372. ENDIF
  373.  
  374.  
  375. NOMP=LINCO.MOTS(1)
  376. NOMD=LINCO.MOTS(2)(1:4)
  377.  
  378. NRIGE=7
  379. NKID =9
  380. NKMT =7
  381. NMATRI=1
  382. IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.2)NMATRI=2
  383. IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.4.AND.IDIM.EQ.2)NMATRI=2
  384. IF(KPRE.EQ.5)NMATRI=NMATRI+1
  385. C write(6,*)' NMATRI=',nmatri,' KPRE=',kpre
  386. SEGINI MATRIK
  387.  
  388. C CAS Stabilisation via MACRO CENTRE
  389. IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.2)THEN
  390. IF(IKAS.EQ.3)BETA0=-BETA0
  391. C write(6,*)' CAS Stabilisation via MACRO CENTRE ',BETA0,IKAS
  392. NK=0
  393. I2=2
  394. NBME=1
  395. NBSOUS=1
  396. SEGINI IMATRI
  397. IRIGEL(4,i2)=IMATRI
  398. KSPGP=MCTREI
  399. KSPGD=MCTREI
  400. IRIGEL(1,i2)=MELEMA
  401. IRIGEL(2,i2)=MELEMA
  402. IRIGEL(7,i2)=0
  403. CALL LICHT(MCHPOC,MPOVAL,TYPC,IGEOM)
  404.  
  405. C write(6,*)' MELSTB=',melstb,' MCTREI=',MCTREI
  406.  
  407. SEGACT MELSTB
  408.  
  409. NBSOUS=MELSTB.LISOUS(/1)
  410. IF(NBSOUS.NE.0)THEN
  411. CALL ERREUR(5)
  412. ENDIF
  413.  
  414. NBEL=MELEMA.NUM(/2)
  415. NBCI=MELSTB.NUM(/2)
  416. NP =MELSTB.NUM(/1)
  417. MP =NP
  418.  
  419. C write(6,*)' nbel,nbci,np,mp=',nbel,nbci,np,mp
  420. SEGINI IZAFM
  421. LIZAFM(1,1)=IZAFM
  422. LISPRI(1)=NOMD
  423. LISDUA(1)=NOMD
  424.  
  425. CALL KRIPAD(MCTREI,MLENT1)
  426.  
  427. C write(6,*)' IK1=',ik1
  428. IF(IK1.LT.4)THEN
  429. DO 33 K=1,NBEL
  430.  
  431. NK=NK+1
  432. KC=1+(1-IK1)*(NK-1)
  433. C write(6,*)' K=',K,' NK=',nk,np,kc
  434.  
  435. DO 32 J=1,NP
  436. K1=MLENT1.LECT(MELEMA.NUM(J,K))
  437. II=J
  438. DO 321 I=1,NP
  439. U=VPOCHA(K1,I)*BETA0*IZTGG1.VPOCHA(KC,1)
  440. IF(I.EQ.1)U=ABS(VPOCHA(K1,I))*BETA0*IZTGG1.VPOCHA(KC,1)
  441. IF(II.LE.NP)THEN
  442. AM(K,II,J)=U
  443. ELSE
  444. AM(K,II-NP,J)=U
  445. ENDIF
  446. II=II+1
  447. 321 CONTINUE
  448. 32 CONTINUE
  449. 33 CONTINUE
  450. C write(6,*)' FIN BCL 33 '
  451. ELSE
  452.  
  453. DO 34 K=1,NBEL
  454.  
  455. NK=NK+1
  456. UA=1.D0
  457.  
  458. DO 35 J=1,NP
  459. K1=MLENT1.LECT(MELEMA.NUM(J,K))
  460. II=J
  461. DO 351 I=1,NP
  462. U=VPOCHA(K1,I)*BETA0*UA
  463. IF(I.EQ.1)U=ABS(VPOCHA(K1,I))*BETA0*UA
  464. IF(II.LE.NP)THEN
  465. AM(K,II,J)=U
  466. ELSE
  467. AM(K,II-NP,J)=U
  468. ENDIF
  469. II=II+1
  470. 351 CONTINUE
  471. 35 CONTINUE
  472. 34 CONTINUE
  473.  
  474.  
  475. ENDIF
  476. SEGSUP MLENT1
  477. ENDIF
  478. C write(6,*)' Apres 34 '
  479.  
  480. C CAS Stabilisation via MACRO CENTREP1
  481. DEUPI=1.D0
  482. IF(IAXI.NE.0)DEUPI=2.D0*XPI
  483. IF(IKAS.EQ.3)DEUPI=-DEUPI
  484.  
  485. IF(MACRO1.NE.0.AND.IKAS.NE.2.AND.KPRE.EQ.4.AND.IDIM.EQ.2)THEN
  486. I2=2
  487. NBME=1
  488. NBSOUS=1
  489. SEGINI IMATRI
  490. IRIGEL(4,i2)=IMATRI
  491. KSPGP=MELEMC
  492. KSPGD=MELEMC
  493. IRIGEL(1,i2)=MELEMP
  494. IRIGEL(2,i2)=MELEMP
  495. IRIGEL(7,i2)=0
  496.  
  497. SEGACT MELEMP
  498. NBSOUS=MELEMP.LISOUS(/1)
  499. IF(NBSOUS.NE.0)THEN
  500. CALL ERREUR(5)
  501. ENDIF
  502.  
  503. SEGACT MELEMI
  504. NBSOUM=MELEMI.LISOUS(/1)
  505. IF(NBSOUM.EQ.0)NBSOUM=1
  506.  
  507. NBEL=MELEMP.NUM(/2)
  508. NP =MELEMP.NUM(/1)
  509. MP =NP
  510. NOM0=NOMS(MELEMP.ITYPEL)//' '
  511.  
  512. SEGINI IZAFM
  513. LIZAFM(1,1)=IZAFM
  514. LISPRI(1)=NOMD
  515. LISDUA(1)=NOMD
  516.  
  517. CALL KALPBG(NOM0,'FONFORM ',IZFFM)
  518.  
  519. SEGACT IZFFM*MOD
  520. IZHR=KZHR(1)
  521. SEGACT IZHR*MOD
  522. NES=GR(/1)
  523. NPG=GR(/3)
  524.  
  525. NK=0
  526. DO 430 LM=1,NBSOUM
  527. IPT1=MELEMI
  528. IF(NBSOUM.NE.1)IPT1=MELEMI.LISOUS(LM)
  529. SEGACT IPT1
  530. ITYP=IPT1.ITYPEL
  531. BETA=0.D0
  532. IF(NOMS(ITYP).EQ.'TRI6')BETA=BETA0
  533. C write(6,*)' NOMS(ITYP)=',NOMS(ITYP),BETA
  534. NBEL=IPT1.NUM(/2)
  535.  
  536. C Cas coefficient FLOTTANT ou CENTRE
  537. IF(IK1.LT.4)THEN
  538. DO 43 KK=1,NBEL
  539. NK=NK+1
  540. DO 42 I=1,NP
  541. J=MELEMP.NUM(I,NK)
  542. DO 42 N=1,IDIM
  543. XYZ(N,I)=XCOOR((J-1)*(IDIM+1) +N)
  544. 42 CONTINUE
  545.  
  546. CALL CALJBR
  547. &(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI,AIRE,AJ,SGN)
  548. KC=1+(1-IK1)*(NK-1)
  549.  
  550. DO 44 I=1,NP
  551. DO 44 J=1,NP
  552. U=0.D0
  553. DO 433 L=1,NPG
  554. UH=0.D0
  555. DO 4333 KH=1,IDIM
  556. UH=UH+HR(KH,J,L)*HR(KH,I,L)
  557. 4333 CONTINUE
  558. U=U-UH*PGSQ(L)*DEUPI*BETA*IZTGG1.VPOCHA(KC,1)
  559. 433 CONTINUE
  560. AM(NK,I,J)=U
  561. 44 CONTINUE
  562.  
  563. 43 CONTINUE
  564.  
  565. ELSEIF(IK1.EQ.4.OR.IK1.EQ.5)THEN
  566. C Cas coefficient SOMMET
  567.  
  568. DO 45 KK=1,NBEL
  569. NK=NK+1
  570. DO 46 I=1,NP
  571. J=MELEMP.NUM(I,NK)
  572. DO 46 N=1,IDIM
  573. XYZ(N,I)=XCOOR((J-1)*(IDIM+1) +N)
  574. 46 CONTINUE
  575.  
  576. CALL CALJBR
  577. &(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI,AIRE,AJ,SGN)
  578.  
  579. DO 47 I=1,NP
  580. DO 47 J=1,NP
  581. U=0.D0
  582. DO 453 L=1,NPG
  583. UH=0.D0
  584. DO 4533 KH=1,IDIM
  585. UH=UH+HR(KH,J,L)*HR(KH,I,L)
  586. 4533 CONTINUE
  587.  
  588. UA=0.D0
  589. DO 454 JJ=1,NP
  590. J1=LECT(IPT1.NUM(JJ,NK))
  591. DO 454 KH=1,IDIM
  592. KG=(IK1-4)*(KH-1)+1
  593. UA=UA+FN(JJ,L)*IZTGG1.VPOCHA(J1,KG)
  594. 454 CONTINUE
  595. UA=UA/KG
  596.  
  597. U=U-UH*PGSQ(L)*DEUPI*BETA*UA
  598. 453 CONTINUE
  599. AM(NK,I,J)=U
  600. 47 CONTINUE
  601.  
  602. 45 CONTINUE
  603. ENDIF
  604.  
  605. SEGDES IPT1
  606. 430 CONTINUE
  607.  
  608. ENDIF
  609. C Fin Stabilisation de toutes sortes
  610. C write(6,*)'C Fin Stabilisation de toutes sortes'
  611.  
  612. NBME=IDIM
  613. CALL KRIPAD(MELEMI,MLENTI1)
  614. SEGACT MELEMI
  615. NBSOUS=MELEMI.LISOUS(/1)
  616. IF(NBSOUS.EQ.0)NBSOUS=1
  617. SEGINI IMATRI
  618.  
  619. IF(IKAS.EQ.2)THEN
  620. KSPGD=MELEMS
  621. KSPGP=MELEMC
  622. IRIGEL(2,1)=MELEMI
  623. IRIGEL(1,1)=MELEMP
  624. ELSE
  625. KSPGP=MELEMS
  626. KSPGD=MELEMC
  627. IRIGEL(1,1)=MELEMI
  628. IRIGEL(2,1)=MELEMP
  629. ENDIF
  630. SEGACT MELEMP
  631.  
  632. C write(6,*)' ds kmac melemp=',IRIGEL(1,1)
  633. C write(6,*)' ds kmac melemd=',IRIGEL(2,1)
  634.  
  635. IRIGEL(4,1)=IMATRI
  636. * GBM modifie le 22 nov 99 pour pas de multipl de lagrange
  637. IF(IKAS.EQ.1)IRIGEL(7,1)=3
  638. IF(IKAS.EQ.2)IRIGEL(7,1)=3
  639. IF(IKAS.EQ.3)IRIGEL(7,1)=4
  640.  
  641. SEGACT MELEMP
  642. NBSOUC=MELEMP.LISOUS(/1)
  643. IF(NBSOUC.EQ.0)NBSOUC=1
  644. C write(6,*)'NBSOUS,NBSOUC=',NBSOUS,NBSOUC
  645.  
  646. K0=0
  647. DO 11 L=1,NBSOUS
  648. IPT1=MELEMI
  649. IPT2=MELEMP
  650. IF(NBSOUS.NE.1)IPT1=MELEMI.LISOUS(L)
  651. IF(NBSOUC.NE.1)IPT2=MELEMP.LISOUS(L)
  652. SEGACT IPT1,IPT2
  653. NBEL=IPT1.NUM(/2)
  654.  
  655. IF(IKAS.EQ.2)THEN
  656. MP=IPT1.NUM(/1)
  657. C avt NP=MELEMP.NUM(/1)
  658. NP=IPT2.NUM(/1)
  659. ELSE
  660. NP=IPT1.NUM(/1)
  661. C avt MP=MELEMP.NUM(/1)
  662. MP=IPT2.NUM(/1)
  663. ENDIF
  664.  
  665. DO 12 I=1,NBME
  666. SEGINI IZAFM
  667. LIZAFM(L,I)=IZAFM
  668. IF(IKAS.EQ.2)THEN
  669. WRITE(NOM,FMT='(I1,A3)')I,NOMD(1:3)
  670. LISDUA(I)=NOM//' '
  671. LISPRI(I)=NOMP
  672. ELSE
  673. WRITE(NOM,FMT='(I1,A3)')I,NOMP(1:3)
  674. LISPRI(I)=NOM//' '
  675. LISDUA(I)=NOMD
  676. ENDIF
  677. 12 CONTINUE
  678. IPM1=LIZAFM(L,1)
  679. IPM2=LIZAFM(L,2)
  680. IPM3=LIZAFM(L,2)
  681. IF(IDIM.EQ.3)IPM3=LIZAFM(L,3)
  682.  
  683. IF(IK1.LT.4)THEN
  684. C write(6,*)' Appel KPRUSS'
  685. CALL KPRUSS(IPT1,IPT2,
  686. &IPM1,IPM2,IPM3,IAXI,IKAS,INEFMD,KPRE,IZTGG1,IK1,K0,IKOMP)
  687. C write(6,*)'Retour KPRUSS'
  688. ELSE
  689. C write(6,*)' Appel KPRJSS'
  690. CALL KPRJSS(IPT1,IPT2,IPM1,IPM2,IPM3,
  691. & IAXI,IKAS,INEFMD,KPRE,IZTGG1,MLENTI,IK1,IKOMP)
  692. ENDIF
  693.  
  694. K0=K0+NBEL
  695. SEGDES IPM1,IPM2,IPM3
  696. SEGDES IPT1
  697. 11 CONTINUE
  698. SEGDES MELEMI
  699. IF(IK1.EQ.4)SEGSUP MLENTI
  700.  
  701.  
  702. C%%%%%%%%%%%%%%%%% CAS DES PRESSIONS CONTINUES %%%%%%%%%%%%%%%%%%%%%%%
  703.  
  704.  
  705. IF(KPRE.EQ.5)THEN
  706. CALL LEKTAB(MTABZ,'MMAIL',MELEME)
  707. CALL LEKTAB(MTABZ,'MSOMMET',MELEMS)
  708.  
  709. SEGACT MELEME
  710. NBSOU1=LISOUS(/1)
  711. IF(NBSOU1.EQ.0)NBSOU1=1
  712.  
  713. IRIGEL(1,2)=MELEME
  714. IRIGEL(2,2)=MELEME
  715. IRIGEL(7,2)=0
  716.  
  717. c%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  718.  
  719. IHV=0
  720. AIMPL=1.D0
  721.  
  722. C --- CALCUL DU NOMBRE DE PAQUETS DE LRV ELEMENTS ---
  723. CALL KPLRVM(MELEME,LRV,NBSOUS)
  724. NBME=1
  725.  
  726. SEGINI IMATRI
  727. IRIGEL(4,2)=IMATRI
  728. KSPGP=MELEMS
  729. KSPGD=MELEMS
  730. LISPRI(1)=NOMD//' '
  731. LISDUA(1)=NOMD//' '
  732.  
  733. LL=0
  734. NUTOEL=0
  735. DO 101 L=1,NBSOU1
  736. IPT1=MELEME
  737. IF(NBSOU1.NE.1)IPT1=LISOUS(L)
  738. SEGACT IPT1
  739. NOM0=NOMS(IPT1.ITYPEL)//' '
  740. CALL KALPBG(NOM0,'FONFORM ',IZFFM)
  741. SEGACT IZFFM*MOD
  742. IZHR=KZHR(1)
  743. SEGACT IZHR*MOD
  744. NES=GR(/1)
  745. NPG=GR(/3)
  746.  
  747. NP = IPT1.NUM(/1)
  748. MP = NP
  749. NBELEM=IPT1.NUM(/2)
  750.  
  751. NNN=MOD(NBELEM,LRV)
  752. IF(NNN.EQ.0) NPACK=NBELEM/LRV
  753. IF(NNN.NE.0) NPACK=1+(NBELEM-NNN)/LRV
  754.  
  755. DO 701 KPACK=1,NPACK
  756.  
  757. C --- POUR CHAQUE PAQUET DE LRV ELEMENTS ou moins
  758. KDEB=1+(KPACK-1)*LRV
  759. KFIN=MIN(NBELEM,KDEB+LRV-1)
  760.  
  761. NBEL=KFIN-KDEB+1
  762. LL=LL+1
  763.  
  764. SEGINI IPM1,IPS1
  765. C write(6,*)' IPM1,LL=',IPM1,LL
  766. LIZAFM(LL,1)=IPM1
  767.  
  768. IPM2=IPM1
  769. IPM3=IPM1
  770. IPS2=IPS1
  771. IPS3=IPS1
  772.  
  773. KITT=2
  774. KJTT=IK1
  775. VISCO=IZTGG1
  776. CALL INITD(VISC,3,1.D0)
  777.  
  778. CALL XLAPL(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI,
  779. & VISC,VISC,VISC,KITT,KJTT,IK1,
  780. & IPT1.NUM(1,KDEB),NBEL,NUTOEL,XCOOR,AIMPL,IKOMP,
  781. & IPM1.AM,IPM2.AM,IPM3.AM,
  782. & IPS1.AM,IPS2.AM,IPS3.AM,
  783. & NINKO,IHV,IARG,VISC)
  784.  
  785. NUTOEL=NUTOEL+NBEL
  786. SEGDES IPM1
  787.  
  788. 701 CONTINUE
  789. SEGDES IPT1
  790. 101 CONTINUE
  791.  
  792.  
  793.  
  794.  
  795.  
  796.  
  797. c%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  798.  
  799. ENDIF
  800.  
  801. CALL ECROBJ('MATRIK',MATRIK)
  802.  
  803. NAT=2
  804. NSOUPO=0
  805. SEGINI MCHPOI
  806. JATTRI(1)=2
  807. CALL ECROBJ('CHPOINT',MCHPOI)
  808.  
  809. C write(6,*)' Fin operateur KMIC'
  810. SEGDES IMATRI,MATRIK
  811. RETURN
  812. 1001 FORMAT(20(1X,I5))
  813. 1002 FORMAT(10(1X,1PE11.4))
  814. END
  815.  
  816.  
  817.  
  818.  
  819.  
  820.  
  821.  
  822.  
  823.  
  824.  
  825.  
  826.  
  827.  
  828.  
  829.  
  830.  

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