Télécharger la plaquette Cast3M


fr en

Télécharger ylap2c.eso

Retour à la liste

Numérotation des lignes : 

ylap2c

  1. C YLAP2C    SOURCE    CB215821  20/11/25    13:44:16     10792          
  2.       SUBROUTINE YLAP2C(ICOGRV,MPROC,MPVITC,
  3.      $     MPVOLU,MPNORM,MPSURF,MELEFL,
  4.      $     KRFACE,KRCENT,LCLIMV,KRVIMP,LCLITO,KRTOIM,
  5.      $     NOMINC,
  6.      $     MU,
  7.      $     IJACO,
  8.      $     IMPR,IRET)
  9.       IMPLICIT INTEGER(I-N)
  10.       IMPLICIT REAL*8 (A-H,O-Z)
  11. C***********************************************************************
  12. C NOM         : YLAP2C
  13. C DESCRIPTION : Calcul de la matrice jacobienne du résidu du laplacien
  14. C               VF 3D.
  15. C               Ici, on ne calcule que les contributions 'simples'
  16. C               à la matrice jacobienne faisant intervenir les
  17. C               coefficients pour le calcul des gradients de vitesse
  18. C               (ICOGRV).
  19. C               (contributions à (d Res_{\rho u} / d var) et (d
  20. C                Res_{\rho v} / d var)
  21. C                var prenant successivement les valeurs :
  22. C                \rho,   \rho u,    \rho v,     \rho e_t )
  23. C
  24. C
  25. C NOTE PERSO  : On pourrait programmer ça plus lisiblement en stockant
  26. C               les contributions dans un tableau de pointeurs (2
  27. C               indices, c'est possible ?) et en effectuant des produits
  28. C               matriciels (coeff. x matrices de dérivées).
  29. C               On n'y coupera pas si on veut regrouper 2D et 3D...
  30. C               On ne va pas le faire.
  31. C
  32. C
  33. C LANGAGE     : ESOPE
  34. C AUTEUR      : Stéphane GOUNAND (CEA/DEN/DM2S/SFME/LTMF)
  35. C               mél : gounand@semt2.smts.cea.fr
  36. C***********************************************************************
  37. C APPELES (UTIL)   : AJMTK : ajoute un objet de type MATSIM (non
  38. C                            standard) à un objet de type MATRIK.
  39. C APPELE PAR       : YLAP2A : Calcul de la matrice jacobienne du
  40. C                    résidu du laplacien VF 3D.
  41. C***********************************************************************
  42. C ENTREES            : ICOGRV (type MCHELM) : coefficients pour le
  43. C                        calcul du gradient de la vitesse aux
  44. C                        interfaces.
  45. C                      MPROC (type MPOVAL) : masse volumique par
  46. C                                            élément.
  47. C                      MPVITC (type MPOVAL) : vitesse par élément.
  48. C                      MPVOLU (type MPOVAL) : volume des éléments.
  49. C                      MPNORM (type MPOVAL) : normale aux faces.
  50. C                      MPSURF (type MPOVAL) : surface des faces.
  51. C                      MELEFL (type MELEME) : connectivités face-(centre
  52. C                                             gauche, centre droit).
  53. C                      KRFACE (type MLENTI) : tableau de repérage dans
  54. C                        le maillage des faces des éléments.
  55. C                      KRCENT (type MLENTI) : tableau de repérage dans
  56. C                        le maillage des centres des éléments.
  57. C                      LCLIMV (type logique) : .TRUE. => CL de Dirichlet
  58. C                        sur la vitesse.
  59. C                      KRVIMP (type MLENTI) : tableau de repérage dans
  60. C                        maillage des CL de Dirichlet sur la
  61. C                        vitesse.
  62. C                      LCLITO (type logique) : .TRUE. => CL de Dirichlet
  63. C                        sur le tenseur des contraintes.
  64. C                      KRTOIM (type MLENTI) : tableau de repérage dans
  65. C                        maillage des CL de Dirichlet sur le tenseur des
  66. C                        contraintes
  67. C                      NOMINC (type MLMOTS) : noms des inconnues.
  68. C                      MU (type réel) : viscosité dynamique (SI).
  69. C ENTREES/SORTIES    : IJACO (type MATRIK) : matrice jacobienne du
  70. C                        résidu du laplacien VF 3D.
  71. C SORTIES            : -
  72. C CODE RETOUR (IRET) : = 0 si tout s'est bien passé
  73. C***********************************************************************
  74. C VERSION    : v1, 28/08/2001, version initiale
  75. C HISTORIQUE : v1, 28/08/2001, création
  76. C HISTORIQUE :
  77. C HISTORIQUE :
  78. C***********************************************************************
  79. C Prière de PRENDRE LE TEMPS de compléter les commentaires
  80. C en cas de modification de ce sous-programme afin de faciliter
  81. C la maintenance !
  82. C***********************************************************************
  83.  
  84. -INC PPARAM
  85. -INC CCOPTIO
  86. -INC SMCHPOI
  87.       POINTEUR MPROC.MPOVAL ,MPVITC.MPOVAL
  88.       POINTEUR MPSURF.MPOVAL,MPNORM.MPOVAL,MPVOLU.MPOVAL
  89. -INC SMCHAML
  90.       POINTEUR ICOGRV.MCHELM,JCOGRV.MCHAML
  91.       POINTEUR KDUNDX.MELVAL,KDUNDY.MELVAL,KDUNDZ.MELVAL
  92. -INC SMELEME
  93.       POINTEUR MELEFL.MELEME
  94.       POINTEUR MCOGRV.MELEME
  95. -INC SMLENTI
  96.       POINTEUR KRVIMP.MLENTI,KRTOIM.MLENTI
  97.       POINTEUR KRCENT.MLENTI,KRFACE.MLENTI
  98. -INC SMLMOTS
  99.       POINTEUR NOMINC.MLMOTS
  101.       POINTEUR IJACO.MATRIK
  102. *
  103. * Objet matrice élémentaire simplifié
  104. *
  105.       SEGMENT GMATSI
  106.       INTEGER POIPR1(NPP1,NEL1)
  107.       INTEGER POIDU1(1,NEL1)
  108.       INTEGER POIPR2(NPP2,NEL2)
  109.       INTEGER POIDU2(2,NEL2)
  110.       POINTEUR LMATSI(0).MATSIM
  111.       ENDSEGMENT
  112. *  Contributions simples de la part du gradient de
  113. *  vitesse (CTSGRV)
  114.       POINTEUR CTSGRV.GMATSI
  115.       SEGMENT MATSIM
  116.       CHARACTER*8 NOMPRI,NOMDUA
  117.       REAL*8 VALMA1(1,NPP1,NEL1)
  118.       REAL*8 VALMA2(2,NPP2,NEL2)
  119.       ENDSEGMENT
  120.       POINTEUR ROURHO.MATSIM
  121.       POINTEUR ROVRHO.MATSIM
  122.       POINTEUR ROWRHO.MATSIM
  123.       POINTEUR ROUROU.MATSIM
  124.       POINTEUR ROVROU.MATSIM
  125.       POINTEUR ROWROU.MATSIM
  126.       POINTEUR ROUROV.MATSIM
  127.       POINTEUR ROVROV.MATSIM
  128.       POINTEUR ROWROV.MATSIM
  129.       POINTEUR ROUROW.MATSIM
  130.       POINTEUR ROVROW.MATSIM
  131.       POINTEUR ROWROW.MATSIM
  132. *
  133.       REAL*8 MU
  134. *
  135.       INTEGER IMPR,IRET
  136. *
  137.       LOGICAL LCLIMV,LCLITO
  138.       LOGICAL LMUR
  139.       LOGICAL LCTRB1
  140. *
  141.       INTEGER IELEM,IPD,IPP,ISOUCH,IEL1,IEL2
  142.       INTEGER NELEM,NPD,NPP,NSOUCH,NEL1,NEL2,NPP1,NPP2
  143.       INTEGER NGCDRO,NGCGAU,NGFACE,NPPRIM,NPDUAL
  144.       INTEGER NLCENP,NLCEND,NLFACE,NLCLTO,NLCLV
  145.       INTEGER NPTEL
  146. *
  147.       REAL*8 ALPHAX,ALPHAY,ALPHAZ,CNX,CNY,CNZ
  148.       REAL*8 SIGNOR,SURFFA,VOLUEL
  149.       REAL*8 RHOP,UP,VP,WP
  150.       REAL*8 FACTOR
  151.       REAL*8 DUDRHO,DUDROU
  152.       REAL*8 DVDRHO,DVDROV
  153.       REAL*8 DWDRHO,DWDROW
  154.       REAL*8 BROUDU,BROUDV,BROUDW
  155.       REAL*8 BROVDU,BROVDV,BROVDW
  156.       REAL*8 BROWDU,BROWDV,BROWDW
  157. C
  158. *
  159. * Executable statements
  160. *
  161.       IF (IMPR.GT.2) WRITE(IOIMP,*) 'Entrée dans ylap2c.eso'
  162. *   On calcule les contributions à (d Res_{\rho u} / d var), (d
  163. *   Res_{\rho v} / d var) et (d Res_{\rho w} / d var) ; var prenant
  164. *   successivement les valeurs :
  165. *      \rho,   \rho u,    \rho v,  \rho w,    \rho e_t
  166. *   Note :
  167. *     - (d Res_{\rho u} / d \rho e_t) = 0
  168. *     - (d Res_{\rho v} / d \rho e_t) = 0
  169. *     - (d Res_{\rho w} / d \rho e_t) = 0
  170. *   On dérive les termes : \tens{\tau} \prod \vect{n}
  171. *   Les noms de matrices élémentaires (type MATSIM) associées sont :
  172. *   ROURHO, ROUROU, ROUROV, ROUROW,
  173. *   ROVRHO, ROVROU, ROVROV, ROVROW,
  174. *   ROWRHO, ROWROU, ROWROV, ROWROW.
  175.       IF (LCLIMV) THEN
  176.          SEGACT KRVIMP
  177.       ENDIF
  178.       IF (LCLITO) THEN
  179.          SEGACT KRTOIM
  180.       ENDIF
  181.       SEGACT NOMINC
  182.       SEGACT KRCENT
  183.       SEGACT KRFACE
  184.       SEGACT MELEFL
  185.       SEGACT MPSURF
  186.       SEGACT MPNORM
  187.       SEGACT MPVOLU
  188.       SEGACT MPROC
  189.       SEGACT MPVITC
  190.       SEGACT ICOGRV
  191.       NSOUCH=ICOGRV.IMACHE(/1)
  192.       DO 1 ISOUCH=1,NSOUCH
  193.          MCOGRV=ICOGRV.IMACHE(ISOUCH)
  194.          JCOGRV=ICOGRV.ICHAML(ISOUCH)
  195.          SEGACT JCOGRV
  196.          KDUNDX=JCOGRV.IELVAL(1)
  197.          KDUNDY=JCOGRV.IELVAL(2)
  198.          KDUNDZ=JCOGRV.IELVAL(3)
  199.          SEGDES JCOGRV
  200.          SEGACT KDUNDX
  201.          SEGACT KDUNDY
  202.          SEGACT KDUNDZ
  203.          SEGACT MCOGRV
  204.          NELEM=MCOGRV.NUM(/2)
  205.          NPTEL=MCOGRV.NUM(/1)
  206.          NPP1=NPTEL-1
  207.          NPP2=NPTEL-1
  208.          NEL1=NELEM
  209.          NEL2=NELEM
  210.          IEL1=1
  211.          IEL2=1
  212.          SEGINI ROURHO
  213.          SEGINI ROVRHO
  214.          SEGINI ROWRHO
  215.          SEGINI ROUROU
  216.          SEGINI ROVROU
  217.          SEGINI ROWROU
  218.          SEGINI ROUROV
  219.          SEGINI ROVROV
  220.          SEGINI ROWROV
  221.          SEGINI ROUROW
  222.          SEGINI ROVROW
  223.          SEGINI ROWROW
  224.          SEGINI CTSGRV
  225.          ROURHO.NOMPRI(1:4)=NOMINC.MOTS(1)
  226.          ROURHO.NOMPRI(5:8)='    '
  227.          ROURHO.NOMDUA(1:4)=NOMINC.MOTS(2)
  228.          ROURHO.NOMDUA(5:8)='    '
  229.          ROVRHO.NOMPRI(1:4)=NOMINC.MOTS(1)
  230.          ROVRHO.NOMPRI(5:8)='    '
  231.          ROVRHO.NOMDUA(1:4)=NOMINC.MOTS(3)
  232.          ROVRHO.NOMDUA(5:8)='    '
  233.          ROWRHO.NOMPRI(1:4)=NOMINC.MOTS(1)
  234.          ROWRHO.NOMPRI(5:8)='    '
  235.          ROWRHO.NOMDUA(1:4)=NOMINC.MOTS(4)
  236.          ROWRHO.NOMDUA(5:8)='    '
  237.          ROUROU.NOMPRI(1:4)=NOMINC.MOTS(2)
  238.          ROUROU.NOMPRI(5:8)='    '
  239.          ROUROU.NOMDUA(1:4)=NOMINC.MOTS(2)
  240.          ROUROU.NOMDUA(5:8)='    '
  241.          ROVROU.NOMPRI(1:4)=NOMINC.MOTS(2)
  242.          ROVROU.NOMPRI(5:8)='    '
  243.          ROVROU.NOMDUA(1:4)=NOMINC.MOTS(3)
  244.          ROVROU.NOMDUA(5:8)='    '
  245.          ROWROU.NOMPRI(1:4)=NOMINC.MOTS(2)
  246.          ROWROU.NOMPRI(5:8)='    '
  247.          ROWROU.NOMDUA(1:4)=NOMINC.MOTS(4)
  248.          ROWROU.NOMDUA(5:8)='    '
  249.          ROUROV.NOMPRI(1:4)=NOMINC.MOTS(3)
  250.          ROUROV.NOMPRI(5:8)='    '
  251.          ROUROV.NOMDUA(1:4)=NOMINC.MOTS(2)
  252.          ROUROV.NOMDUA(5:8)='    '
  253.          ROVROV.NOMPRI(1:4)=NOMINC.MOTS(3)
  254.          ROVROV.NOMPRI(5:8)='    '
  255.          ROVROV.NOMDUA(1:4)=NOMINC.MOTS(3)
  256.          ROVROV.NOMDUA(5:8)='    '
  257.          ROWROV.NOMPRI(1:4)=NOMINC.MOTS(3)
  258.          ROWROV.NOMPRI(5:8)='    '
  259.          ROWROV.NOMDUA(1:4)=NOMINC.MOTS(4)
  260.          ROWROV.NOMDUA(5:8)='    '
  261.          ROUROW.NOMPRI(1:4)=NOMINC.MOTS(4)
  262.          ROUROW.NOMPRI(5:8)='    '
  263.          ROUROW.NOMDUA(1:4)=NOMINC.MOTS(2)
  264.          ROUROW.NOMDUA(5:8)='    '
  265.          ROVROW.NOMPRI(1:4)=NOMINC.MOTS(4)
  266.          ROVROW.NOMPRI(5:8)='    '
  267.          ROVROW.NOMDUA(1:4)=NOMINC.MOTS(3)
  268.          ROVROW.NOMDUA(5:8)='    '
  269.          ROWROW.NOMPRI(1:4)=NOMINC.MOTS(4)
  270.          ROWROW.NOMPRI(5:8)='    '
  271.          ROWROW.NOMDUA(1:4)=NOMINC.MOTS(4)
  272.          ROWROW.NOMDUA(5:8)='    '
  273.          DO 12 IELEM=1,NELEM
  274. *     Le premier point du support de ICOGRV est un point FACE
  275.             NGFACE=MCOGRV.NUM(1,IELEM)
  276.             NLFACE=KRFACE.LECT(NGFACE)
  277.             IF (NLFACE.EQ.0) THEN
  278.                WRITE(IOIMP,*) 'Erreur de programmation n°1'
  279.                GOTO 9999
  280.             ENDIF
  281. *     On calcule la contribution à la matrice jacobienne IJACO de la face
  282. *     NGFAC (points duaux     : centres à gauche et à droite de la face)
  283. *           (points primaux   : une partie (bicoz conditions aux limites)
  284. *            de ceux du stencil pour le calcul du gradient
  285. *            à la face, ils doivent être des points centres)
  286. *       Si le tenseur des contraintes sur la face est imposé par les
  287. *       conditions aux limites, la contribution de la face à IJACO est
  288. *       nulle.
  289.             LCTRB1=.TRUE.
  290.             IF (LCLITO) THEN
  291.                NLCLTO=KRTOIM.LECT(NGFACE)
  292.                IF (NLCLTO.NE.0) THEN
  293.                   LCTRB1=.FALSE.
  294.                ENDIF
  295.             ENDIF
  296.             IF (LCTRB1) THEN
  297.                NGCGAU=MELEFL.NUM(1,NLFACE)
  298.                NGCDRO=MELEFL.NUM(3,NLFACE)
  299.                LMUR=(NGCGAU.EQ.NGCDRO)
  300. *     On distingue le cas où la face est un bord du maillage (mur)
  301. *                  du cas où la face est interne au maillage
  302.                IF (.NOT.LMUR) THEN
  303.                   NPD=2
  304.                ELSE
  305.                   NPD=1
  306.                ENDIF
  307.                NPP=NPTEL-1
  308. *     IPD=1 : point à gauche du point NGFACE
  309. *     IPD=2 : point à droite du point NGFACE
  310.                DO 122 IPD=1,NPD
  311.                   NPDUAL=MELEFL.NUM((2*IPD)-1,NLFACE)
  312.                   IF (.NOT.LMUR) THEN
  313.                      CTSGRV.POIDU2(IPD,IEL2)=NPDUAL
  314.                   ELSE
  315.                      CTSGRV.POIDU1(IPD,IEL1)=NPDUAL
  316.                   ENDIF
  317.                   NLCEND=KRCENT.LECT(NPDUAL)
  318.                   IF (NLCEND.EQ.0) THEN
  319.                      WRITE(IOIMP,*) 'Erreur grave n°1'
  320.                      GOTO 9999
  321.                   ENDIF
  322.                   DO 124 IPP=1,NPP
  323.                      NPPRIM=MCOGRV.NUM(IPP+1,IELEM)
  324.                      NLCENP=KRCENT.LECT(NPPRIM)
  325. C
  326. C AB We do not check BC anymore
  327. C
  328.                      IF(NLCENP .EQ. 0)THEN
  329. *     Lorsque une contribution est nulle, on fixe artificiellement le
  330. *     point primal égal au point dual.
  331.                         IF (.NOT.LMUR) THEN
  332.                            CTSGRV.POIPR2(IPP,IEL2)=NPDUAL
  333.                            ROURHO.VALMA2(IPD,IPP,IEL2)=0.D0
  334.                            ROVRHO.VALMA2(IPD,IPP,IEL2)=0.D0
  335.                            ROWRHO.VALMA2(IPD,IPP,IEL2)=0.D0
  336.                            ROUROU.VALMA2(IPD,IPP,IEL2)=0.D0
  337.                            ROVROU.VALMA2(IPD,IPP,IEL2)=0.D0
  338.                            ROWROU.VALMA2(IPD,IPP,IEL2)=0.D0
  339.                            ROUROV.VALMA2(IPD,IPP,IEL2)=0.D0
  340.                            ROVROV.VALMA2(IPD,IPP,IEL2)=0.D0
  341.                            ROWROV.VALMA2(IPD,IPP,IEL2)=0.D0
  342.                            ROUROW.VALMA2(IPD,IPP,IEL2)=0.D0
  343.                            ROVROW.VALMA2(IPD,IPP,IEL2)=0.D0
  344.                            ROWROW.VALMA2(IPD,IPP,IEL2)=0.D0
  345.                         ELSE
  346.                            CTSGRV.POIPR1(IPP,IEL1)=NPDUAL
  347.                            ROURHO.VALMA1(IPD,IPP,IEL1)=0.D0
  348.                            ROVRHO.VALMA1(IPD,IPP,IEL1)=0.D0
  349.                            ROWRHO.VALMA1(IPD,IPP,IEL1)=0.D0
  350.                            ROUROU.VALMA1(IPD,IPP,IEL1)=0.D0
  351.                            ROVROU.VALMA1(IPD,IPP,IEL1)=0.D0
  352.                            ROWROU.VALMA1(IPD,IPP,IEL1)=0.D0
  353.                            ROUROV.VALMA1(IPD,IPP,IEL1)=0.D0
  354.                            ROVROV.VALMA1(IPD,IPP,IEL1)=0.D0
  355.                            ROWROV.VALMA1(IPD,IPP,IEL1)=0.D0
  356.                            ROUROW.VALMA1(IPD,IPP,IEL1)=0.D0
  357.                            ROVROW.VALMA1(IPD,IPP,IEL1)=0.D0
  358.                            ROWROW.VALMA1(IPD,IPP,IEL1)=0.D0
  359.                         ENDIF
  360.                      ELSE
  361. *     Les contributions valent :
  362. * 1.  (d Res_{\rho u})_d / (d var)_p =
  363. *     +/-1 (normale sortante, rentrante)  (1/V_d) * (S_f) * \mu
  364. *     * [ [ ( 4/3 (n_x * \alpha_x) + (n_y * \alpha_y)
  365. *                                  + (n_z * \alpha_z)) *
  366. *           ((du)_p / (d var)_p)]
  367. *       + [ (-2/3 (n_x * \alpha_y) + (n_y * \alpha_x)) *
  368. *           ((dv)_p / (d var)_p)]
  369. *       + [ (-2/3 (n_x * \alpha_z) + (n_z * \alpha_x)) *
  370. *           ((dw)_p / (d var)_p)]
  371. *       ]
  372. * 2.  (d Res_{\rho v})_d / (d var)_p =
  373. *     +/-1 (normale sortante, rentrante)  (1/V_d) * (S_f) * \mu
  374. *     * [ [ (-2/3 (n_y * \alpha_x) + (n_x * \alpha_y)) *
  375. *           ((du)_p / (d var)_p)]
  376. *       + [ ( 4/3 (n_y * \alpha_y) + (n_x * \alpha_x)
  377. *                                  + (n_z * \alpha_z)) *
  378. *           ((dv)_p / (d var)_p)]
  379. *       + [ (-2/3 (n_y * \alpha_z) + (n_z * \alpha_y)) *
  380. *           ((dw)_p / (d var)_p)]
  381. *       ]
  382. * 3.  (d Res_{\rho w})_d / (d var)_p =
  383. *     +/-1 (normale sortante, rentrante)  (1/V_d) * (S_f) * \mu
  384. *     * [ [ (-2/3 (n_z * \alpha_x) + (n_x * \alpha_z)) *
  385. *           ((du)_p / (d var)_p)]
  386. *       + [ (-2/3 (n_z * \alpha_y) + (n_y * \alpha_z)) *
  387. *           ((dv)_p / (d var)_p)]
  388. *       + [ ( 4/3 (n_z * \alpha_z) + (n_x * \alpha_x)
  389. *                                  + (n_y * \alpha_y)) *
  390. *           ((dw)_p / (d var)_p)]
  391. *       ]
  392. *
  393. *     avec :
  394. *          (du)_p / (d \rho)_p   = - u / \rho_p
  395. *          (du)_p / (d \rho u)_p =   1 / \rho_p
  396. *          (du)_p / (d \rho v)_p =   0
  397. *          (du)_p / (d \rho w)_p =   0
  398. *          (dv)_p / (d \rho)_p   = - v / \rho_p
  399. *          (dv)_p / (d \rho u)_p =   0
  400. *          (dv)_p / (d \rho v)_p =   1 / \rho_p
  401. *          (dv)_p / (d \rho w)_p =   0
  402. *          (dv)_p / (d \rho)_p   = - w / \rho_p
  403. *          (dv)_p / (d \rho u)_p =   0
  404. *          (dv)_p / (d \rho v)_p =   0
  405. *          (dv)_p / (d \rho w)_p =   1 / \rho_p
  406. *
  407. *       normale sortante pour IPD=1, rentrante pour IPD=2
  408.                         SIGNOR=(-1.D0)**(IPD+1)
  409.                         VOLUEL=MPVOLU.VPOCHA(NLCEND,1)
  410.                         SURFFA=MPSURF.VPOCHA(NLFACE,1)
  411.                         CNX   =MPNORM.VPOCHA(NLFACE,1)
  412.                         CNY   =MPNORM.VPOCHA(NLFACE,2)
  413.                         CNZ   =MPNORM.VPOCHA(NLFACE,3)
  414.                         ALPHAX=KDUNDX.VELCHE(IPP+1,IELEM)
  415.                         ALPHAY=KDUNDY.VELCHE(IPP+1,IELEM)
  416.                         ALPHAZ=KDUNDZ.VELCHE(IPP+1,IELEM)
  417.                         RHOP  =MPROC.VPOCHA(NLCENP,1)
  418.                         UP    =MPVITC.VPOCHA(NLCENP,1)
  419.                         VP    =MPVITC.VPOCHA(NLCENP,2)
  420.                         WP    =MPVITC.VPOCHA(NLCENP,3)
  421.                         FACTOR=SIGNOR*(1.D0/VOLUEL)*SURFFA*MU
  422.                         BROUDU=(( 4.D0/3.D0)*(CNX*ALPHAX))
  423.      $                       + (CNY*ALPHAY) + (CNZ*ALPHAZ)
  424.                         BROUDV=((-2.D0/3.D0)*(CNX*ALPHAY))
  425.      $                       + (CNY*ALPHAX)
  426.                         BROUDW=((-2.D0/3.D0)*(CNX*ALPHAZ))
  427.      $                       + (CNZ*ALPHAX)
  428.                         BROVDU=((-2.D0/3.D0)*(CNY*ALPHAX))
  429.      $                       + (CNX*ALPHAY)
  430.                         BROVDV=(( 4.D0/3.D0)*(CNY*ALPHAY))
  431.      $                       + (CNX*ALPHAX) + (CNZ*ALPHAZ)
  432.                         BROVDW=((-2.D0/3.D0)*(CNY*ALPHAZ))
  433.      $                       + (CNZ*ALPHAY)
  434.                         BROWDU=((-2.D0/3.D0)*(CNZ*ALPHAX))
  435.      $                       + (CNX*ALPHAZ)
  436.                         BROWDV=((-2.D0/3.D0)*(CNZ*ALPHAY))
  437.      $                       + (CNY*ALPHAZ)
  438.                         BROWDW=(( 4.D0/3.D0)*(CNZ*ALPHAZ))
  439.      $                       + (CNX*ALPHAX) + (CNY*ALPHAY)
  440.                         DUDRHO=-UP /RHOP
  441.                         DUDROU=1.D0/RHOP
  442.                         DVDRHO=-VP /RHOP
  443.                         DVDROV=1.D0/RHOP
  444.                         DWDRHO=-WP /RHOP
  445.                         DWDROW=1.D0/RHOP
  446.                         IF (.NOT.LMUR) THEN
  447.                            CTSGRV.POIPR2(IPP,IEL2)=NPPRIM
  448.                            ROURHO.VALMA2(IPD,IPP,IEL2)=
  449.      $                          FACTOR*((BROUDU*DUDRHO)
  450.      $                          +(BROUDV*DVDRHO)+(BROUDW*DWDRHO))
  451.                            ROVRHO.VALMA2(IPD,IPP,IEL2)=
  452.      $                          FACTOR*((BROVDU*DUDRHO)
  453.      $                          +(BROVDV*DVDRHO)+(BROVDW*DWDRHO))
  454.                            ROWRHO.VALMA2(IPD,IPP,IEL2)=
  455.      $                          FACTOR*((BROWDU*DUDRHO)
  456.      $                          +(BROWDV*DVDRHO)+(BROWDW*DWDRHO))
  457.                            ROUROU.VALMA2(IPD,IPP,IEL2)=
  458.      $                          FACTOR* (BROUDU*DUDROU)
  459.                            ROVROU.VALMA2(IPD,IPP,IEL2)=
  460.      $                          FACTOR* (BROVDU*DUDROU)
  461.                            ROWROU.VALMA2(IPD,IPP,IEL2)=
  462.      $                          FACTOR* (BROWDU*DUDROU)
  463.                            ROUROV.VALMA2(IPD,IPP,IEL2)=
  464.      $                          FACTOR* (BROUDV*DVDROV)
  465.                            ROVROV.VALMA2(IPD,IPP,IEL2)=
  466.      $                          FACTOR* (BROVDV*DVDROV)
  467.                            ROWROV.VALMA2(IPD,IPP,IEL2)=
  468.      $                          FACTOR* (BROWDV*DVDROV)
  469.                            ROUROW.VALMA2(IPD,IPP,IEL2)=
  470.      $                          FACTOR* (BROUDW*DWDROW)
  471.                            ROVROW.VALMA2(IPD,IPP,IEL2)=
  472.      $                          FACTOR* (BROVDW*DWDROW)
  473.                            ROWROW.VALMA2(IPD,IPP,IEL2)=
  474.      $                          FACTOR* (BROWDW*DWDROW)
  475.                         ELSE
  476.                            CTSGRV.POIPR1(IPP,IEL1)=NPPRIM
  477.                            ROURHO.VALMA1(IPD,IPP,IEL1)=
  478.      $                          FACTOR*((BROUDU*DUDRHO)
  479.      $                          +(BROUDV*DVDRHO)+(BROUDW*DWDRHO))
  480.                            ROVRHO.VALMA1(IPD,IPP,IEL1)=
  481.      $                          FACTOR*((BROVDU*DUDRHO)
  482.      $                          +(BROVDV*DVDRHO)+(BROVDW*DWDRHO))
  483.                            ROWRHO.VALMA1(IPD,IPP,IEL1)=
  484.      $                          FACTOR*((BROWDU*DUDRHO)
  485.      $                          +(BROWDV*DVDRHO)+(BROWDW*DWDRHO))
  486.                            ROUROU.VALMA1(IPD,IPP,IEL1)=
  487.      $                          FACTOR* (BROUDU*DUDROU)
  488.                            ROVROU.VALMA1(IPD,IPP,IEL1)=
  489.      $                          FACTOR* (BROVDU*DUDROU)
  490.                            ROWROU.VALMA1(IPD,IPP,IEL1)=
  491.      $                          FACTOR* (BROWDU*DUDROU)
  492.                            ROUROV.VALMA1(IPD,IPP,IEL1)=
  493.      $                          FACTOR* (BROUDV*DVDROV)
  494.                            ROVROV.VALMA1(IPD,IPP,IEL1)=
  495.      $                          FACTOR* (BROVDV*DVDROV)
  496.                            ROWROV.VALMA1(IPD,IPP,IEL1)=
  497.      $                          FACTOR* (BROWDV*DVDROV)
  498.                            ROUROW.VALMA1(IPD,IPP,IEL1)=
  499.      $                          FACTOR* (BROUDW*DWDROW)
  500.                            ROVROW.VALMA1(IPD,IPP,IEL1)=
  501.      $                          FACTOR* (BROVDW*DWDROW)
  502.                            ROWROW.VALMA1(IPD,IPP,IEL1)=
  503.      $                          FACTOR* (BROWDW*DWDROW)
  504.                         ENDIF
  505.                      ENDIF
  506.  124              CONTINUE
  507.  122           CONTINUE
  508.                IF (.NOT.LMUR) THEN
  509.                   IEL2=IEL2+1
  510.                ELSE
  511.                   IEL1=IEL1+1
  512.                ENDIF
  513.             ENDIF
  514.  12      CONTINUE
  515.          NPP1=NPTEL-1
  516.          NPP2=NPTEL-1
  517.          NEL1=IEL1-1
  518.          NEL2=IEL2-1
  519.          SEGADJ ROURHO
  520.          SEGADJ ROVRHO
  521.          SEGADJ ROWRHO
  522.          SEGADJ ROUROU
  523.          SEGADJ ROVROU
  524.          SEGADJ ROWROU
  525.          SEGADJ ROUROV
  526.          SEGADJ ROVROV
  527.          SEGADJ ROWROV
  528.          SEGADJ ROUROW
  529.          SEGADJ ROVROW
  530.          SEGADJ ROWROW
  531.          SEGADJ CTSGRV
  532.          CTSGRV.LMATSI(**)=ROURHO
  533.          CTSGRV.LMATSI(**)=ROVRHO
  534.          CTSGRV.LMATSI(**)=ROWRHO
  535.          CTSGRV.LMATSI(**)=ROUROU
  536.          CTSGRV.LMATSI(**)=ROVROU
  537.          CTSGRV.LMATSI(**)=ROWROU
  538.          CTSGRV.LMATSI(**)=ROUROV
  539.          CTSGRV.LMATSI(**)=ROVROV
  540.          CTSGRV.LMATSI(**)=ROWROV
  541.          CTSGRV.LMATSI(**)=ROUROW
  542.          CTSGRV.LMATSI(**)=ROVROW
  543.          CTSGRV.LMATSI(**)=ROWROW
  544. * On accumule les matrices résultantes dans IJACO
  545.          CALL AJMTK(CTSGRV,IJACO,IMPR,IRET)
  546.          IF (IRET.NE.0) GOTO 9999
  547.          SEGSUP ROURHO
  548.          SEGSUP ROVRHO
  549.          SEGSUP ROWRHO
  550.          SEGSUP ROUROU
  551.          SEGSUP ROVROU
  552.          SEGSUP ROWROU
  553.          SEGSUP ROUROV
  554.          SEGSUP ROVROV
  555.          SEGSUP ROWROV
  556.          SEGSUP ROUROW
  557.          SEGSUP ROVROW
  558.          SEGSUP ROWROW
  559.          SEGSUP CTSGRV
  560. *
  561.          SEGDES MCOGRV
  562.          SEGDES KDUNDZ
  563.          SEGDES KDUNDY
  564.          SEGDES KDUNDX
  565.  1    CONTINUE
  566.       SEGDES ICOGRV
  567.       SEGDES MPVITC
  568.       SEGDES MPROC
  569.       SEGDES MPVOLU
  570.       SEGDES MPNORM
  571.       SEGDES MPSURF
  572.       SEGDES MELEFL
  573.       SEGDES KRFACE
  574.       SEGDES KRCENT
  575.       SEGDES NOMINC
  576.       IF (LCLITO) THEN
  577.          SEGDES KRTOIM
  578.       ENDIF
  579.       IF (LCLIMV) THEN
  580.          SEGDES KRVIMP
  581.       ENDIF
  582. *
  583. * Normal termination
  584. *
  585.       IRET=0
  586.       RETURN
  587. *
  588. * Format handling
  589. *
  590. *
  591. * Error handling
  592. *
  593.  9999 CONTINUE
  594.       IRET=1
  595.       WRITE(IOIMP,*) 'An error was detected in subroutine ylap2c'
  596.       RETURN
  597. *
  598. * End of subroutine YLAP2C
  599. *
  600.       END
  601.  
  602.  
  603.  
  604.  
  605.  
  606.  
  607.  
  608.  
  609.  
  610.  
  611.  
  612.  
  613.  

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