Télécharger la plaquette Cast3M


fr en

Télécharger yfrt.eso

Retour à la liste

Numérotation des lignes : 

yfrt

  1. C YFRT      SOURCE    PV090527  26/04/30    21:16:49     12529          
  2.       SUBROUTINE YFRT(NOMPR,MTABX,IHV,MZIN,NOMI,MATRIK,MCHPO1)
  3.       IMPLICIT INTEGER(I-N)
  4.       IMPLICIT REAL*8 (A-H,O-Z)
  5. C***********************************************************************
  6. C
  7. C  CET OPERATEUR DISCRETISE Convection Diffusion + source
  8. C
  9. C  SYNTAXE :
  10. C  ---------
  11. C
  12. C
  13. C  COEFFICIENTS :
  14. C  --------------
  15. C
  16. C
  17. C  ALF  (SCAL DOMA)  DIFFUSIVITE THERMIQUE
  18. C       (SCAL ELEM)
  19. C       (SCAL NOEU)
  20. C
  21. C  INCONNUES :
  22. C  -----------
  23. C
  24. C   TN                 CHAMP DE TEMPERATURE
  25. C   UN                 CHAMP DE VITESSE TRANSPORTANT
  26. C
  27. C   NOMPR Nom de l'opérateur appelant (NS,TSCA,KONV ou LAPN)
  28. C         suivant l'opérateur les coefs sont différents
  29. C   Schéma en temps ou l'opérateur intervient :
  30. C    SEMI            CNG                         TVISQ
  31. C    teta - schéma   Crank Nicolson généralisé   Tenseur visqueux
  32. C   Ces schéma nécessitent que DFDT soit en CENTREE
  33. C   Le Schéma BDF2 ne fait intervenir que DFDT
  34. C
  35. C E/  MTABX Table de l'opérateur -> coefficients
  36. C E/  IHV=0 SCALAIRE
  37. C     IHV=1 VITESSE
  38. C E/  MZIN  pointeur du CHPOINT de l'inconnue
  39. C E/  NOMI  nom de l'inconnue
  40. C************************************************************************
  41.  
  42.  
  43. -INC PPARAM
  44. -INC CCOPTIO
  45. -INC CCREEL
  46. -INC CCGEOME
  47. -INC SMCOORD
  48. -INC SMLENTI
  49.       POINTEUR MLENT4.MLENTI
  50. -INC SMLMOTS
  51. -INC SMMODEL
  52. -INC SMELEME
  53.       POINTEUR IGEOM.MELEME,MELEMS.MELEME,MELEMC.MELEME,MELEMP.MELEME
  54.       POINTEUR MELEM2.MELEME
  55. -INC SMCHPOI
  56.       POINTEUR MZIN.MCHPOI
  57.       POINTEUR MTETA1.MPOVAL,MTETA2.MPOVAL,MTETA3.MPOVAL,MTETA4.MPOVAL
  58. -INC SMTABLE
  59.       POINTEUR MTABZ.MTABLE
  60. -INC SIZFFB
  61.       POINTEUR IZF1.IZFFM,IZH2.IZHR,IZW.IZFFM,IZWH.IZHR
  62.       SEGMENT SAJT
  63.       REAL*8 AJT(IDIM,IDIM,NPG),RF1(NP,MP,IDIM),SM1(NP,IDIM)
  64.       REAL*8 TN1(NP,IDIM),TN2(NP,IDIM)
  65.       REAL*8 WT(MP,NPG), WS(MP,NPG)
  66.       REAL*8 WTC(IDIM*NPG,MP),WSC(IDIM*NPG,MP),HRD(IDIM*NPG,NP)
  67.       REAL*8 FND(NPG,NP)
  68.       ENDSEGMENT
  69. -INC SMRIGID
  70. C-INC SMMATRIK
  71. C*******************************************************************
  72. C
  73.       SEGMENT MATRIK
  74.       REAL*8  COEMTK(NMATRI)
  75.       INTEGER JRIGEL(NRIGE,NMATRI)
  76.       INTEGER KSYM,KMINC,KMINCP,KMINCD,KIZM
  77.       INTEGER KISPGT,KISPGP,KISPGD
  78.       INTEGER KNTTT,KNTTP,KNTTD
  79.       INTEGER KIDMAT(NKID)
  80.       INTEGER KKMMT(NKMT)
  81.       ENDSEGMENT
  82.  
  83.       SEGMENT JMATRI
  84.       CHARACTER*8 LISPRI(NBMF),LJSDUA(NBMF)
  85.       INTEGER LIZAFM(NBSOUS,NBMF)
  86.       INTEGER KSPGP,KSPGD
  87.       ENDSEGMENT
  88.       POINTEUR JMATRS.JMATRI,JMATR1.JMATRI,JMATR2.JMATRI,JMATR3.JMATRI
  89.  
  90. C Stokage matrices elementaires non assemblees (valeurs)
  91.       SEGMENT IZAFM
  92.       REAL*8  AM(NBEL,NP,MP)
  93.       ENDSEGMENT
  94.       POINTEUR IPM1.IZAFM,IPM2.IZAFM,IPM3.IZAFM,IPM4.IZAFM
  95.  
  96. C*******************************************************************
  97.       POINTEUR IPM.IZAFM
  98. -INC SMCHAML
  99.  
  100.       LOGICAL TDFDT,TCONV,TSOUR,TSOURB,TMATRI
  101.       LOGICAL ICAL,XPG,XRIG
  102.       LOGICAL XDIAG,XTV,XTG,XBDF,XCONS
  103.       CHARACTER*(*) NOMPR,NOMI
  104.       CHARACTER*8 CHHH,TYPE,NOM,NOM0,CHAI,SCHT,NOMPER,MPRE
  105.       CHARACTER*4 INCOD
  106.       CHARACTER*(LOCOMP) NOMP(3),NOMD(3),NMACO
  107.       DIMENSION ICOEF(4),SU(3),XPOI(3)
  108. C*****************************************************************************
  109. CTCLSF
  110.       NOMPER=NOMPR
  111.       write(6,*)' DEBUT YFRT appele par ',NOMPER
  112.  
  113. C- Table de l'opérateur (pointeur MTABX)
  114.  
  115. C- Récupération de la table RV (pointeur MTAB1)
  116.       CALL LEKTAB(MTABX,'EQEX',MTAB1)
  117. C- Récupération de la table des options de l'opérateur (pointeur KOPTI)
  118.       CALL LEKTAB(MTABX,'KOPT',KOPTI)
  119. C- Récupération de la table DOMAINE de la zone (pointeur MTABZ)
  120.       CALL LEKTAB(MTABX,'DOMZ',MTABZ)
  121. C- Récupération de la table INCO (pointeur KINC)
  122.       CALL LEKTAB(MTAB1,'INCO',KINC)
  123.       IF(IERR.NE.0) RETURN
  124.  
  125. C*****************************************************************************
  126. C Traitement des options
  127. C KIMPL = 0 -> EXPL     1 -> IMPL     2 -> CN
  128. C                                   AIMPL=0 explicite
  129. C KKALE=1  Indice indiquant que l'on est en ALE
  130. C sinon KKALE=0
  131. C KFORM = 0 -> SI       1 -> EF       2 -> VF    3 -> EFMC
  132. C IDCEN  : entier indiquant le type de décentrement souhaité
  133. C IDCEN 1-> CENTREE  2-> SUPGDC  3-> SUPG  4-> TVISQUEU  5-> CNG
  134. C E/  CMD    : Coefficient multiplicateur du décentrement
  135. C Si IDCEN=4 ou =5 CMD=DT
  136.  
  137.       CALL ACME(KOPTI,'KFORM',KFORM)
  138.       CALL ACME(KOPTI,'KIMPL',KIMPL)
  139.       CALL ACMF(KOPTI,'AIMPL',AIMPL)
  140.       IKOMP=0
  141. *?    CALL ACME(KOPTI,'IKOMP',IKOMP)
  142.       CALL ACME(KOPTI,'KMACO',KMACO)
  143.       KKMACO=KMACO
  144.       CALL ACME(KOPTI,'IDCEN',IDCEN)
  145.       CALL ACMF(KOPTI,'CMD',CMD)
  146. *?    CALL ACME(KOPTI,'ALE',KKALE)
  147.  
  148. C Restrictions des options
  149.       IF(NOMPER.NE.'NS')KKALE=0
  150.       IF(NOMPER.EQ.'LAPN')IDCEN=1
  151.  
  152.       IF(KIMPL.EQ.0)AIMPL=0.D0
  153.       IF(KIMPL.EQ.1)AIMPL=1.D0
  154.  
  155. c       write(6,*)' YFRT KMACO=',KMACO,'IDCEN=',IDCEN
  156.        IF(KMACO.EQ.1.AND.IDCEN.NE.5)THEN
  157.          CALL ACMM(KOPTI,'NMACO',NMACO)
  158.          TYPE=' '
  159.          CALL ACMO(MTAB1,NMACO,TYPE,MATRIK)
  160.         IF(TYPE.EQ.'MATRIK')THEN
  161.           KKMACO=2
  162.         ELSE
  163.           KKMACO=1
  164.         ENDIF
  165.        ENDIF
  166.       IF(IDCEN.EQ.5)KKMACO=0
  167. c      write(6,*)' YFRT KMACO=',kmaco,' KKMACO=',KKMACO,'IKOMP=',ikomp
  168.  
  169. C*****************************************************************************
  170.  
  171.       KPOIND=0
  172.       NBMF=1
  173.       XRIG=.FALSE.
  174.       TMATRI=.TRUE.
  175.       XPG=.FALSE.
  176.       IF(IDCEN.GT.1)XPG=.TRUE.
  177.       IF(IDCEN.EQ.4)AIMPL=0.D0
  178.       IF(IDCEN.EQ.5)THEN
  179.        AIMPL=1.D0
  180.       ENDIF
  181.  
  182.        AIEX=AIMPL
  183.       IF(AIMPL.EQ.0.D0.AND.KKMACO.NE.0)THEN
  184.        TMATRI=.TRUE.
  185.        AIEX=1.D0
  186.       ENDIF
  187.  
  188.       IF(AIMPL.EQ.0.D0.AND.KKMACO.EQ.0)THEN
  189.        TMATRI=.FALSE.
  190.       ENDIF
  191.  
  192.       IAXI=0
  193.       IF(IFOMOD.EQ.0)IAXI=2
  194.       DEUPI=1.D0
  195.       IF(IAXI.NE.0)DEUPI=2.D0*XPI
  196.  
  197.         IF(XRIG)THEN
  198.          IF(IHV.EQ.0)THEN
  199.           NOMP(1)='T   '
  200.           NOMD(1)='Q   '
  201.          ELSEIF(IHV.EQ.1)THEN
  202.           NOMP(1)='UX  '
  203.           NOMP(2)='UY  '
  204.           NOMP(3)='UZ  '
  205.           NOMD(1)='FX  '
  206.           NOMD(2)='FY  '
  207.           NOMD(3)='FZ  '
  208.          ENDIF
  209.         ELSE
  210.          IF(IHV.EQ.0)THEN
  211.           NOMP(1)=NOMI
  212.           NOMD(1)=NOMI
  213.          ELSEIF(IHV.EQ.1)THEN
  214.           WRITE(NOMP(1),FMT='(I1,A3)')1,NOMI(1:LOCOMP-1)
  215.           WRITE(NOMD(1),FMT='(I1,A3)')1,NOMI(1:LOCOMP-1)
  216.           WRITE(NOMP(2),FMT='(I1,A3)')2,NOMI(1:LOCOMP-1)
  217.           WRITE(NOMD(2),FMT='(I1,A3)')2,NOMI(1:LOCOMP-1)
  218.           WRITE(NOMP(3),FMT='(I1,A3)')3,NOMI(1:LOCOMP-1)
  219.           WRITE(NOMD(3),FMT='(I1,A3)')3,NOMI(1:LOCOMP-1)
  220.          ENDIF
  221.         ENDIF
  222. c       write(6,*)' NOMP=',nomp
  223. c       write(6,*)' NOMD=',nomd
  224.  
  225.       CALL LEKTAB(MTABZ,'SOMMET',MELEMS)
  226.       CALL LEKTAB(MTABZ,'MAILLAGE',MELEME)
  227.       MELEM2=MELEME
  228.  
  229. C  On cree un second membre non vide
  230.             SEGACT MELEMS
  231.             N=MELEMS.NUM(/2)
  232.             IF(IHV.EQ.0)NC=1
  233.             IF(IHV.EQ.1)NC=IDIM
  234.             SEGDES MELEMS
  235.             NSOUPO=1
  236.             NAT=2
  237.             SEGINI MCHPO1,MSOUP1,MPOVA1
  238.             MCHPO1.JATTRI(1)=2
  239.             MCHPO1.IFOPOI=IFOUR
  240.             MCHPO1.MTYPOI='    '
  241.             MCHPO1.MOCHDE='                                            '
  242.             MCHPO1.IPCHP(1)=MSOUP1
  243.             MSOUP1.IGEOC=MELEMS
  244.             MSOUP1.IPOVAL=MPOVA1
  245.              IF(IHV.EQ.0)THEN
  246.               MSOUP1.NOCOMP(1)=NOMD(1)
  247.              ELSE
  248.               DO 91 N=1,NC
  249.               MSOUP1.NOCOMP(N)=NOMD(N)
  250.  91           CONTINUE
  251.              ENDIF
  252.             SEGDES MSOUP1,MCHPO1
  253.  
  254. C*************************************************************************
  255. C Lecture des Inco(s) aux temps precedents si Transitoire
  256. c       write(6,*)'Lecture des INCO MZIN=',MZIN
  257.  
  258.         CALL LRCHT(MZIN,MTETA1,TYPE,IGEOM)
  259.         CALL KRIPAD(IGEOM,MLENT1)
  260.         CALL VERPAD(MLENT1,MELEMS,IRET)
  261.          IF(IRET.NE.0)THEN
  262. C%    Indice %m1:8 : L'objet %m9:16 n'a pas le bon support géométrique
  263.           MOTERR(1: 8) = 'TETA1'
  264.           MOTERR(9:16) = 'CHPOINT '
  265.           WRITE(IOIMP,*)'Opérateur : ',NOMPER
  266.           CALL ERREUR(788)
  267.           RETURN
  268.          ENDIF
  269.  
  270.         NC=MTETA1.VPOCHA(/2)
  271.  
  272.         IF((IHV.EQ.0.AND.NC.NE.1).OR.
  273.      &     (IHV.EQ.1.AND.NC.NE.IDIM))THEN
  274. C% Indice %m1:8 : L'objet %m9:16 n'a pas le bon nombre de composantes
  275.          MOTERR( 1: 8) = 'Inconnue'
  276.          MOTERR( 9:16) = 'CHPOINT'
  277.          WRITE(IOIMP,*)'Operateur : ',NOMPER
  278.          CALL ERREUR(784)
  279.          RETURN
  280.         ENDIF
  281.  
  282.         SEGDES IGEOM,MELEMS
  283.  
  284. C*****************************************************************************
  285. C Lecture du coefficient
  286. c       write(6,*)'Lecture des COEF '
  287.  
  288.       CALL ACME(MTABX,'IARG',IARG)
  289. c      write(6,*)' YFRT IARG=',IARG,' KKALE=',KKALE
  290.  
  291. C--------- CAS FROT  Formulation non conservative ----------------------------
  292.        IF(NOMPER.EQ.'FROT')THEN
  293.  
  294. C Coef K={Kx Ky Kz}
  295.        NUCOEF=1
  296. c      write(6,*)' Coef 1 ',NUCOEF,' MTABZ=',MTABZ
  297.        CALL LEKMOF(MTABZ,NUCOEF,MTABX,KINC,1,MCHEL1,KPOIND,0,MCHELG)
  298.        IF (IERR.NE.0) RETURN
  299.  
  300. C Coef Beta={Betax Betay Betaz}
  301.        NUCOEF=NUCOEF+1
  302. c      write(6,*)' Coef 2 ',NUCOEF,' MTABZ=',MTABZ
  303.        CALL LEKMOF(MTABZ,NUCOEF,MTABX,KINC,1,MCHEL2,KPOIND,0,MCHELG)
  304.        IF (IERR.NE.0) RETURN
  305.  
  306. C Un
  307.        NUCOEF=NUCOEF+1
  308. c      write(6,*)' Coef 3 ',NUCOEF
  309.        CALL LEKMOF(MTABZ,NUCOEF,MTABX,KINC,1,MCHEL3,KPOIND,0,MCHELG)
  310.        IF (IERR.NE.0) RETURN
  311.  
  312.  
  313. c       write(6,*)'FIN CAS FROT  Formulation non conservative'
  314. C----- FIN CAS FROT  Formulation non conservative ----------------------------
  315.  
  316.  
  317.        ENDIF
  318. c      write(6,*)' Fin lect COEF '
  319. C*************************************************************************
  320. C******* CALCUL **********************************************************
  321. C*************************************************************************
  322.  
  323.        SEGACT MTETA1
  324.  
  325.       CALL KRIPAD(MELEMS,MLENTI)
  326.       SEGDES MELEMS
  327.  
  328.       SEGACT MELEME,MELEM2
  329.       NBSOUS=LISOUS(/1)
  330.       IF(NBSOUS.EQ.0)NBSOUS=1
  331.  
  332. C MATRIK . MATRIK . MATRIK . MATRIK . MATRIK . MATRIK  MATRIK . MATRIK .
  333.       IF(KKMACO.NE.2)THEN
  334. C ... Création / Ajustement des matrices élémentaires
  335.       IF(TMATRI)THEN
  336. C Cas RIGIDITE
  337.        IF(XRIG)THEN
  338.  
  339.          NRIGE=8
  340.          NRIGEL=0
  341.          SEGINI MRIGID
  342.  
  343.          NRIGE=8
  344.          NRI=IRIGEL(/2)
  345.          NRIGEL=NBSOUS+NRI
  346.          SEGADJ MRIGID
  347. c     write(6,*)' NRIGE,NRIGEL,MRIGID=',NRIGE,NRIGEL,MRIGID
  348.  
  349.        ELSE
  350. C Cas MATRIK
  351.  
  352.          NRIGE=7
  353.          NKID =9
  354.          NKMT =7
  355.          NMATRI=0
  356.          SEGINI MATRIK
  357. c        write(6,*)' YFRT creation MATRIK =',matrik
  358.  
  359.          NRIGE=7
  360.          NKID =9
  361.          NKMT =7
  362.          NMATR0=JRIGEL(/2)
  363.  
  364.          NBME=1
  365.          IF(IHV.EQ.1)NBME=IDIM
  366.          NMATRI=NBME
  367.          SEGADJ MATRIK
  368.  
  369.          DO 41 M=1,NBME
  370.          JRIGEL(1,M)=MELEME
  371.          JRIGEL(2,M)=MELEM2
  372.          JRIGEL(7,M)=0
  373. c?       IF(TCONV)JRIGEL(7,M)=2
  374.          SEGINI JMATRI
  375.          JRIGEL(4,M)=JMATRI
  376.          KSPGP=MELEMS
  377.          KSPGD=MELEMS
  378.          LISPRI(1)=NOMP(M)//'    '
  379.          LJSDUA(1)=NOMD(M)//'    '
  380.  41      CONTINUE
  381.        ENDIF
  382.       ENDIF
  383. C ... Fin Création / Ajustement des matrices élémentaires
  384.  
  385.  
  386.       ELSE
  387. C Activation RIGIDITE
  388.        IF(XRIG)THEN
  389.          SEGACT MRIGID
  390.        ELSE
  391. C Activation MATRIK
  392.          SEGACT MATRIK
  393.        ENDIF
  394.  
  395. C ... Activation des matrices élémentaires
  396.       IF(TMATRI)THEN
  397. C Cas RIGIDITE
  398.        IF(XRIG)THEN
  399.          NRIGE=8
  400.          NRI=IRIGEL(/2)
  401.          NRIGEL=NBSOUS+NRI
  402.          SEGAct MRIGID
  403. c     write(6,*)' NRIGE,NRIGEL,MRIGID=',NRIGE,NRIGEL,MRIGID
  404.  
  405.        ELSE
  406. C Cas MATRIK Pleine
  407.          NMATRI=JRIGEL(/2)
  408.          NBME=NMATRI
  409.  
  410.          DO 42  M=1,NBME
  411.          JMATRI=JRIGEL(4,M)
  412.          SEGACT JMATRI*MOD
  413.          LISPRI(1)=NOMP(M)//'    '
  414.          LJSDUA(1)=NOMD(M)//'    '
  415.  42      CONTINUE
  416.  
  417.        ENDIF
  418.       ENDIF
  419. C ... Fin Activation des matrices élémentaires
  420.  
  421.       ENDIF
  422. C FIN MATRIK . FIN MATRIK . FIN MATRIK . FIN MATRIK . FIN MATRIK . FIN M
  423. C ......................................................................
  424.  
  425.       SEGACT MCHEL1
  426.       SEGACT MCHEL2
  427.       SEGACT MCHEL3
  428.  
  429.        IF(MAX(1,MELEM2.LISOUS(/1)).NE.MAX(1,LISOUS(/1)))THEN
  430.         WRITE(6,*)' Geometries incompatibles dans ',nomper
  431. C% Données incompatibles
  432.         CALL ERREUR(22)
  433.         RETURN
  434.         ENDIF
  435.  
  436.       NKD=0
  437.       DO 101 L=1,MAX(1,LISOUS(/1))
  438.       IPT1=MELEME
  439.       IPT2=MELEM2
  440.       IF(LISOUS(/1).NE.0)IPT1=LISOUS(L)
  441.       IF(MELEM2.LISOUS(/1).NE.0)THEN
  442.       IPT2=MELEM2.LISOUS(L)
  443.       NKD=0
  444.       ENDIF
  445.       SEGACT IPT1,IPT2
  446.  
  447.       NOM0 = NOMS(IPT1.ITYPEL)//'    '
  448.       CALL KALPBG(NOM0,'FONFORM ',IZFFM)
  449.       IF(IZFFM.EQ.0)RETURN
  450.       SEGACT IZFFM*MOD
  451.       IZHR=KZHR(1)
  452.       SEGACT IZHR*MOD
  453.       IZF1 = KTP(1)
  454.       IZH2 = KZHR(2)
  455.       IZW  = IZFFM
  456.       IZWH = IZHR
  457.  
  458.       NES=GR(/1)
  459.       NPG=GR(/3)
  460.  
  461.       NP = IPT1.NUM(/1)
  462.       MP = IPT2.NUM(/1)
  463. C?    MP = IZW.FN(/1) ceci doit etre identique
  464.  
  465.       NBEL=IPT1.NUM(/2)
  466.  
  467.       SEGINI SAJT
  468.  
  469. c.......................................................................
  470. C MATRIK . MATRIK . MATRIK . MATRIK . MATRIK . MATRIK  MATRIK . MATRIK .
  471. c     write(6,*)' KKMACO=',KKMACO,'TMATRI=',TMATRI
  472.       IF(KKMACO.NE.2)THEN
  473.        IF(TMATRI)THEN
  474. C Cas RIGIDITE
  475.         IF(XRIG)THEN
  476.         IRIGEL(1,NRI+L)=MELEME
  477.         COERIG(L)=1.D0
  478.  
  479.         IRIGEL(7,NRI+L)=0
  480. c?      IF(TCONV)IRIGEL(7,NRI+L)=2
  481.  
  482.         NBME=1
  483.         IF(IHV.EQ.1)NBME=IDIM
  484.         NLIGRP=NP
  485.         NLIGRD=MP
  486.         SEGINI DESCR
  487.         IRIGEL(3,NRI+L)=DESCR
  488.          IF(NBME.EQ.1)THEN
  489.            DO 102 I=1,NLIGRP
  490.             LISINC(I)=NOMP(1)//'    '
  491.             NOELEP(I)=I
  492.  102       CONTINUE
  493.            DO 103 I=1,NLIGRD
  494.             LISDUA(I)=NOMD(1)//'    '
  495.             NOELED(I)=I
  496.  103       CONTINUE
  497.          ELSE
  498.          ENDIF
  499.         SEGDES DESCR
  500.  
  501.         NELRIG=NBEL
  502.         rigrel=0
  503.         SEGINI xMATRI
  504.         IRIGEL(4,NRI+L)=xMATRI
  505.         xmatri.symre=irigel(7,nri+l)
  506. c       write(6,*)'NELRIG,IMATRI=',NELRIG,IMATRI
  507.  
  508. *        DO 104 K=1,NELRIG
  509. *        SEGINI XMATRI
  510. c       write(6,*)'NLIGRD,NLIGRP,XMATRI=',NLIGRD,NLIGRP,XMATRI
  511. *        IMATTT(K)=XMATRI
  512. * 104    CONTINUE
  513.  
  514.         ELSE
  515. C Cas MATRIK
  516.         SEGINI IPM1
  517. c       mtail=(IPM1.AM(/1))*(ipm1.am(/2))*(ipm1.am(/3))
  518.         JMATR1=JRIGEL(4,NMATR0+1)
  519.         JMATR1.LIZAFM(L,1)=IPM1
  520.         IPM2=IPM1
  521.         IPM3=IPM1
  522.          IF(NBME.GE.2)THEN
  523.           JMATR2=JRIGEL(4,NMATR0+2)
  524.            SEGINI IPM2
  525. c       mtail=(IPM2.AM(/1))*(ipm2.am(/2))*(ipm2.am(/3))
  526.            JMATR2.LIZAFM(L,1)=IPM2
  527.            ICAL=.TRUE.
  528.          ENDIF
  529.          IF(NBME.GE.3)THEN
  530.           JMATR3=JRIGEL(4,NMATR0+3)
  531.            SEGINI IPM3
  532.            JMATR3.LIZAFM(L,1)=IPM3
  533.            ICAL=.TRUE.
  534.          ENDIF
  535.        ENDIF
  536.       ENDIF
  537. c       write(6,*)' YFRT appelé par ',NOMPER,' Taille matrice=',mtail
  538.  
  539.       ELSE
  540.  
  541. C Cas RIGIDITE
  542.        IF(XRIG)THEN
  543.  
  544.        ELSE
  545. C Cas MATRIK
  546.         JMATR1=JRIGEL(4,1)
  547.         SEGACT JMATR1
  548.         IPM1=JMATR1.LIZAFM(L,1)
  549.         IPM2=IPM1
  550.         IPM3=IPM1
  551.         SEGACT IPM1
  552.          IF(NBME.GE.2)THEN
  553.           JMATR2=JRIGEL(4,2)
  554.           SEGACT JMATR2
  555.           IPM2=JMATR2.LIZAFM(L,1)
  556.           SEGACT IPM2
  557.          ENDIF
  558.          IF(NBME.GE.3)THEN
  559.           JMATR3=JRIGEL(4,3)
  560.           IPM3=JMATR3.LIZAFM(L,1)
  561.           SEGACT IPM3
  562.          ENDIF
  563.        ENDIF
  564.       ENDIF
  565. C FIN MATRIK . FIN MATRIK . FIN MATRIK . FIN MATRIK . FIN MATRIK . FIN M
  566. C ......................................................................
  567.  
  568.  
  569. C----K
  570.       IK1=1
  571.       MCHAM1=MCHEL1.ICHAML(L)
  572.       SEGACT MCHAM1
  573.  
  574.       MELVA1=MCHAM1.IELVAL(1)
  575.       SEGACT MELVA1
  576.       N1PTEL=MELVA1.VELCHE(/1)
  577.       N1EL=MELVA1.VELCHE(/2)
  578.         IF(N1EL.EQ.1)THEN
  579.         IK1=1
  580.         ELSEIF(N1EL.EQ.NBEL)THEN
  581.         IK1=0
  582.         ENDIF
  583.  
  584. C----Beta
  585.       IK2=1
  586.       MCHAM2=MCHEL2.ICHAML(L)
  587.       SEGACT MCHAM2
  588.  
  589.       MELVA2=MCHAM2.IELVAL(1)
  590.       SEGACT MELVA2
  591.       N1PTEL=MELVA2.VELCHE(/1)
  592.       N1EL=MELVA2.VELCHE(/2)
  593.         IF(N1EL.EQ.1)THEN
  594.         IK2=1
  595.         ELSEIF(N1EL.EQ.NBEL)THEN
  596.         IK2=0
  597.         ENDIF
  598.  
  599. C----U
  600.       IK3=1
  601.       MCHAM3=MCHEL3.ICHAML(L)
  602.       SEGACT MCHAM3
  603.  
  604.       MELVA3=MCHAM3.IELVAL(1)
  605.       SEGACT MELVA3
  606.       N1PTEL=MELVA3.VELCHE(/1)
  607.       N1EL=MELVA3.VELCHE(/2)
  608.         IF(N1EL.EQ.1)THEN
  609.         IK3=1
  610.         ELSEIF(N1EL.EQ.NBEL)THEN
  611.         IK3=0
  612.         ENDIF
  613.  
  614.  
  615. c      write(6,*)' AVANT 108 NC=',NC,' NBEL=',NBEL,MP,NP,NC
  616. C===============================================
  617.         AI1=AIMPL-1.D0
  618.         IF(AIMPL.EQ.0.D0)THEN
  619.           AI2=-1.D0
  620.         ELSE
  621.           AI2=AI1/AIMPL
  622.         ENDIF
  623.  
  624.         DO 108 KE=1,NBEL
  625.  
  626.         NK1=KE + IK1*(1 - KE)
  627.         NK2=KE + IK2*(1 - KE)
  628.         NK3=KE + IK3*(1 - KE)
  629.  
  630.         DO     I=1,NP
  631.         J=IPT1.NUM(I,KE)
  632.         DO     N=1,IDIM
  633.         XYZ(N,I)=XCOOR((J-1)*(IDIM+1)+N)
  634.         ENDDO
  635.         ENDDO
  636.  
  637.         DO     I=1,NP
  638.          I1=MLENT1.LECT(IPT1.NUM(I,KE))
  639.         DO     N=1,NC
  640.          TN1(I,N)=MTETA1.VPOCHA(I1,N)
  641.         ENDDO
  642.         ENDDO
  643.  
  644.         CALL CALJBR(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,
  645.      *  IDIM,NP,NPG,IAXI,AIRE,AJ,SGN)
  646.  
  647.         IF (IDCEN.EQ.0.OR.IDCEN.EQ.1) THEN
  648.         CALL RSETD(WT,FN,(NP*NPG))
  649.         CALL RSETD(WS,FN,(NP*NPG))
  650.         ELSE
  651.         CALL CALDEC(WT,WS,XYZ,GR,HR,FN,NES,IDIM,NP,NPG,AJT,
  652.      &  IDCEN,CMD,MELVA1.VELCHE(1,NK1),MELVA2.VELCHE(1,NK2),
  653.      &  MELVA3.VELCHE(1,NK3),TN1,NC,IKOMP,XREF,Aire,KE)
  654.         ENDIF
  655.  
  656.         CALL INITD(SM1,(MP*IDIM),0.D0)
  657.  
  658. C:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
  659.         IF(KKMACO.EQ.2)THEN
  660.  
  661. C ...... Chargement Rigidite ou Matrik
  662. C Cas RIGIDITE
  663.           IF(XRIG)THEN
  664. *           XMATRI=IMATTT(KE)
  665.            DO     I=1,MP
  666.            DO     J=1,NP
  667.            RE(I,J,ke)=RF1(J,I,1)
  668.            ENDDO
  669.            ENDDO
  670. *           SEGDES XMATRI
  671.           ELSE
  672. C Cas MATRIK
  673.            DO     N=1,NBME
  674.            JMATR1=JRIGEL(4,N)
  675.            IPM4=JMATR1.LIZAFM(L,1)
  676.            DO     I=1,MP
  677.            DO     J=1,NP
  678.            RF1(J,I,N)=IPM4.AM(KE,J,I)
  679.            ENDDO
  680.            ENDDO
  681.            ENDDO
  682.           ENDIF
  683. C...... Multiplication Matrice Vecteur
  684.  
  685.         DO 710 I=1,MP
  686.         DO 711 J=1,NP
  687.  
  688.         DO 716 N=1,NC
  689.           SM1(I,N)   = SM1(I,N) + AI2*RF1(J,I,N)*TN1(J,N)
  690.  716    CONTINUE
  691.  
  692.  711    CONTINUE
  693.  710    CONTINUE
  694.  
  695.         ELSEIF(KKMACO.NE.2)THEN
  696. C..............................................................
  697. C...... Frot .... Formulation non conservative
  698. C préparation
  699.  
  700.         DO 416 N=1,NC
  701.         DO 410 I=1,MP
  702. c       DO 411 J=1,NP
  703.           U3=0.D0
  704.         DO 412 LG=1,NPG
  705.           PDR=PGSQ(LG)*DEUPI*RPG(LG)
  706.           C1=MELVA1.VELCHE((N-1)*NPG+LG,NK1)*PDR
  707.           C2=MELVA2.VELCHE((N-1)*NPG+LG,NK2)
  708.           C3=ABS ( MELVA3.VELCHE((N-1)*NPG+LG,NK3) ) + 1.e-30
  709.           CK=C1*(C3**(C2-1.D0))
  710.           U3=U3+CK * FN(I,LG)
  711.  412      CONTINUE
  712.  
  713.           SM1(I,N)   = SM1(I,N) + AI1*U3*TN1(I,N)
  714.           RF1(I,I,N) = AIEX*U3
  715. c?        SM1(I,N)   = SM1(I,N) + AI1*U3*TN1(J,N)
  716. c?        RF1(J,I,N) = AIEX*U3
  717.  
  718.  411    CONTINUE
  719.  410    CONTINUE
  720.  416    CONTINUE
  721.  
  722.  
  723. C..............................................................
  724.  
  725.  
  726.         ENDIF
  727.  
  728. C:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
  729. C=======================================================================
  730.  
  731. C>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
  732.       IF(TMATRI.AND.KKMACO.NE.2)THEN
  733. C ...... Chargement Rigidite ou Matrik
  734. C Cas RIGIDITE
  735.           IF(XRIG)THEN
  736. *           XMATRI=IMATTT(KE)
  737.            DO     I=1,NP
  738.            DO     J=1,NP
  739.            RE(I,J,ke)=RF1(J,I,1)
  740.            ENDDO
  741.            ENDDO
  742. *           SEGDES XMATRI
  743.           ELSE
  744. C Cas MATRIK
  745.            NBMM=1
  746.            IF(ICAL)NBMM=NBME
  747.            DO     N=1,NBMM
  748.            JMATR1=JRIGEL(4,NMATR0+N)
  749.            IPM4=JMATR1.LIZAFM(L,1)
  750.            DO     I=1,NP
  751.            DO     J=1,NP
  752.            IPM4.AM(KE,J,I)=RF1(J,I,N)
  753.            ENDDO
  754.            ENDDO
  755.            ENDDO
  756. c         write(6,*)' RF1 '
  757. c         write(6,1002)(((RF1(J,I,N),j=1,np),i=1,np),n=1,nbmm)
  758.           ENDIF
  759.       ENDIF
  760.  
  761. C ...... Chargement Second membre
  762.         DO     I=1,NP
  763.          I1=LECT(IPT1.NUM(I,KE))
  764.         DO     N=1,NC
  765.          MPOVA1.VPOCHA(I1,N)=MPOVA1.VPOCHA(I1,N)+SM1(I,N)
  766.         ENDDO
  767.         ENDDO
  768. C<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
  769.  
  770.  108    CONTINUE
  771.  
  772.          SEGSUP MCHAM1,MELVA1
  773.          SEGSUP MCHAM2,MELVA2
  774.          SEGSUP MCHAM3,MELVA3
  775.  
  776.  
  777.  
  778.       SEGDES IPT1,IPT2
  779.  
  780.       IF(TMATRI)THEN
  781. C Cas RIGIDITE
  782.        IF(XRIG)THEN
  783.         SEGDES xMATRI
  784.        ELSE
  785. C Cas MATRIK
  786.         SEGDES IPM1
  787.         IF(NBME.GE.2)SEGDES IPM2
  788.         IF(NBME.GE.3)SEGDES IPM3
  789.        ENDIF
  790.       ENDIF
  791.  
  792.       SEGSUP IZFFM,IZHR,IZF1,IZH2
  793.       SEGSUP SAJT
  794.  
  795.  101  CONTINUE
  796.  
  797.       IF(TMATRI)THEN
  798.        IF(.NOT.XRIG)THEN
  799.         NBMM=JRIGEL(/2)
  800.         IF(NBMM.NE.0)THEN
  801.         DO 141 M=1,NBMM
  802.         JMATRI=JRIGEL(4,M)
  803.         SEGDES JMATRI
  804.  141    CONTINUE
  805.         ENDIF
  806.        ENDIF
  807.       ENDIF
  808.  
  809. C"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
  810. C"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
  811. C"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
  812.  
  813.  
  814.          SEGSUP MCHEL1
  815.          SEGSUP MCHEL2
  816.          SEGSUP MCHEL3
  817.  
  818.       SEGDES MCHPO1,MPOVA1
  819.       SEGDES MELEME,MELEM2
  820.  
  821.       SEGSUP MLENTI
  822. c     IF(TLINCO)THEN
  823.       SEGSUP MLENT1
  824.        SEGDES MTETA1
  825. c     ENDIF
  826.  
  827.       IF(TMATRI)THEN
  828. C Cas RIGIDITE
  829.        IF(XRIG)THEN
  830.         SEGDES MRIGID
  831.        ELSE
  832. C Cas MATRIK
  833.         SEGDES MATRIK
  834.        ENDIF
  835.       ENDIF
  836.  
  837.  
  838. c     segact MCHPO1
  839. c     MSOUP1=MCHPO1.IPCHP(1)
  840. c     segact MSOUP1
  841. c     write(6,*)MCHPO1.IPCHP(/1),MSOUP1.NOCOMP(/2)
  842. c     write(6,*)(MSOUP1.NOCOMP(i),i=1,MSOUP1.NOCOMP(/2))
  843. c     if(ihv.eq.1)call ecmo(mtabz,'TOTO','CHPOINT',MCHPO1)
  844.  
  845. C*************************************************************************
  846.  
  847.       IF(KKMACO.EQ.1)THEN
  848.         TYPE='MATRIK'
  849.         CALL ECMO(MTAB1,NMACO,TYPE,MATRIK)
  850.       ENDIF
  851.  
  852.       IF(AIMPL.EQ.0.D0)THEN
  853.         NRIGE=7
  854.         NKID =9
  855.         NKMT =7
  856.         NMATRI=0
  857.         SEGINI MATRIK
  858.       ENDIF
  859.  
  860. c     write(6,*)' FIN YFRT'
  861.       RETURN
  862.  1001 FORMAT(20(1X,I5))
  863.  1002 FORMAT(10(1X,1PE11.4))
  864.       END
  865.  
  866.  
  867.  
  868.  
  869.  
  870.  
  871.  
  872.  
  873.  
  874.  
  875.  
  876.  
  877.  
  878.  
  879.  
  880.  
  881.  
  882.  
  883.  
  884.  
  885.  
  886.  
  887.  
  888.  
  889.  
  890.  
  891.  
  892.  
  893.  
  894.  
  895.  
  896.  
  897.  
  898.  
  899.  
  900.  

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