Sources Esope | Cast3M

ef0kon
C EF0KON    SOURCE    CB215821  20/11/25    13:27:07     10792          
C**************************************************************
C==============================================================
C     Cas de l'equation convection diffusion                  C
C     formulation EFM (elements finis mixte)                  C
C     Auteur : M. MARIN                                       C
C                                                             C
C     IERRKON : indicatif d'erreur = 1 si il y a une erreur   C
C                               = 0 sinon                     C
C     Cette routine cree les matrices elementaires pour un    C
C     terme u*teta ou :                                       C
C                         u : CHPOINT vectoriel au sommet     C
C                         teta : CHPOINT scalaire au centre   C
C                                                             C
C==============================================================
C**************************************************************
 
      SUBROUTINE EF0KON(KIZX,IZTU1,IZTGG1,IZTGG2,IZTGG3,TYPC,IERRKON,
     &     MELEME,MTABZ,NOMI,IKR,MLENTI,IAXI,NOMII,MZPHI,TYPCFI,
     &     IZTCO,NELZ,IKU,IKM,AIMPL,IDCEN,MLENT1,DT)
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)
 
 
-INC CCVQUA4
 
-INC PPARAM
-INC CCOPTIO
-INC CCGEOME
-INC SMCOORD
-INC SMLENTI
-INC SMELEME
      POINTEUR MELEMC.MELEME,MELEMS.MELEME
-INC SMCHAML
-INC SMCHPOI
      POINTEUR IZTU1.MPOVAL,MZPHI.MPOVAL,IZTCO.MPOVAL
      POINTEUR IZTGG1.MPOVAL,IZTGG2.MPOVAL,IZTGG3.MPOVAL
-INC SIZFFB
-INC SMMATRIK
-INC SMLMOTS
-INC SMTABLE
      POINTEUR KIZX.MTABLE,MTABZ.MTABLE
 
C===========================================
C     Variables locale ou argument         C
      CHARACTER*8 TYPC,NOMI,NOM,NOM0,NOMII,TYPCFI
      INTEGER IERRKON,IKR,IKU,IKM
C===========================================
 
c         WRITE(6,*) KIZX,IZTU1,IZTGG1,IZTGG2,IZTGG3,TYPC,IERRKON,
c     &     MELEM,MTABZ,NOMI,IKR,MLENTI,IAXI,NOMII,MZPHI,TYPCFI,
c     &     IZTCO,NELZ,IKU,IKM,AIMPL
C     Validation du passage d'argument
 
      IERRKON=0
 
C     Initialisation pour le nombre total final d'elements (TRI3+QUA4+...)
C     sur le domaine ou s'applique KONV
      NUTOEL=0
c      WRITE(6,*) 'MELEME=',MELEME
 
C      WRITE(6,*) 'Type du Coef 1:',IKR
 
C     On teste si la variable primale est de type CENTRE
      IF (TYPC.NE.'CENTRE  ') THEN
         CALL ERRKON(1,IERRKON)
         RETURN
      END IF
 
C     Test si la variable primale est de type scalaire
      IF (IZTU1.VPOCHA(/2).NE.1) THEN
         CALL ERRKON(5,IERRKON)
         RETURN
      END IF
 
C     On teste si la variable duale est de type SOMMET
      IF (TYPCFI.NE.'SOMMET  ') THEN
         CALL ERRKON(4,IERRKON)
         RETURN
      END IF
 
C     On teste si la variable duale est un vecteur
      IF ((MZPHI.VPOCHA(/2).EQ.1).AND.(IDIM.NE.1)) THEN
         CALL ERRKON(6,IERRKON)
         RETURN
      END IF
 
C     On teste pour savoir si RO est bien un scalaire
      IF (IKR.NE.1) THEN
         CALL ERRKON(7,IERRKON)
         RETURN
      END IF
C     On recupere le nombre de composante de la vitesse
      NBME=IDIM
c      WRITE(6,*) 'Nombre de composante de UN :',NBME
 
C     On active le segment contenant le maillage de la zone ou
C     KONV s'applique et on initialise le nombre de sous-objet.
C     ------------------------------------------------------
C     NBSOUS et NBME servent a initialiser le segment IMATRI
      SEGACT MELEME
      NBSOUS=MELEME.LISOUS(/1)
c      WRITE(6,*) 'NBSOUS=',NBSOUS
      IF (NBSOUS.EQ.0) NBSOUS=1
 
C     Initialisation de la matrice elementaire de l'operateur
      NRIGE=7
      NKID=9
      NKMT=7
      NMATRI=1
      SEGINI MATRIK
 
 
C     On recupere le maillage des centres
      CALL LEKTAB(MTABZ,'CENTRE',MELEMC)
 
C     On recupere le maillage des sommets
      CALL LEKTAB(MTABZ,'SOMMET',MELEMS)
 
C     Inconnue primale
      IRIGEL(1,1)=MELEMC
 
C     Inconnue duale
      IRIGEL(2,1)=MELEME
 
C     Matrice rectangulaire
      IRIGEL(7,1)=3
 
      SEGINI IMATRI
      IRIGEL(4,1)= IMATRI
 
C     Intialisation des supports inconnue primale et duale
C     pour le segment IMATRI
      KSPGP=MELEMC
      KSPGD=MELEMS
 
      DO I=1,NBME
         WRITE(NOM,FMT='(I1,A3)')I,NOMII(1:3)
         LISPRI(I)=NOMI(1:4)//'    '
         LISDUA(I)=NOM(1:4)//'    '
      END DO
 
C     On recupere le segment qui contient int(Ni)
         CALL LEKTAB(MTABZ,'XXPSOML',MCHELM)
 
         IF (MCHELM.EQ.0) THEN
            CALL ERRKON(2,IERRKON)
            RETURN
         END IF
 
         SEGACT MCHELM
 
C==========Boucle principale==============C
C     On effectue cette boucle pour tous les
C     sous-objets
 
 
      DO L=1,NBSOUS
         IPT1=MELEME
         SEGACT MELEME
 
C     Si NBSOUS > 1 on agit sur les sous-objets LISOUS
         IF (NBSOUS.NE.1) IPT1=LISOUS(L)
         SEGACT IPT1
 
c         WRITE(6,*) 'IPT1=',IPT1,'NBSOUS=',NBSOUS
C     On rempli NOM0 du veritable nom des elements traites (TRI3, QUA4,...)
         NOM0=NOMS(IPT1.ITYPEL)//'    '
 
C     Cette routine cree des objets contenant les fonctions
C     de forme lie aux elements traites
C     Le segment IZFFM contient les fonctions de forme et les gradients
C     a chaque point de Gauss
         CALL KALPBG(NOM0,'FONFORM ',IZFFM)
 
         SEGACT IZFFM
         IZHR=KZHR(1)
         SEGACT IZHR*MOD
 
C     Dimension de l'element de reference
         NES=GR(/1)
 
C     Nombre de points de Gauss
         NPG=GR(/3)
c         WRITE(6,*) 'Nbre de points de Gauss :' , NPG
 
C     NP est le nombre de points de l'element de reference
C     MP est la troisieme dimension du tableau contenant les matrices
C     elementaire de chaque elements (segment IZAFM). C'est le
C     nombre de point par element
         NP = 1
         MP = IPT1.NUM(/1)
C     NBEL est le nombre d'element de type NOM0
         NBEL=IPT1.NUM(/2)
 
C     On declare un segment IZAFM par composante
         SEGINI IPM1
         LIZAFM(L,1)=IPM1
 
         IF (IDIM.EQ.2) THEN
            SEGINI IPM2
            LIZAFM(L,2)=IPM2
         ELSE IF (IDIM.EQ.3) THEN
            SEGINI IPM2
            SEGINI IPM3
            LIZAFM(L,2)=IPM2
            LIZAFM(L,3)=IPM3
         END IF
 
C     Nombre total de points NTP
         NPT=IZTGG2.VPOCHA(/1)
 
C     On recupere le tableau contenant les integrales
C     de fonction test
         MCHAML=ICHAML(L)
         SEGACT MCHAML
         MELVAL=IELVAL(1)
         SEGACT MELVAL
         SEGACT MELEMC
C     ******
         CALL XCNEF0(NBEL,NP,MP,IPM1.AM,IZTGG1.VPOCHA,IZTGG2
     &        .VPOCHA,IZTGG3.VPOCHA,NPT,IDIM,IDCEN,XYZ,
     &        NUTOEL,XCOOR,IPT1.NUM,MLENTI.LECT,IPM2.AM,IPM3.AM,
     &        FN,GR,PG,HR,PGSQ,RPG,NES,NPG,IAXI,VELCHE,
     &        NBME,IZTGG3.VPOCHA,IZTCO.VPOCHA,NELZ,
     &        IKR,IKU,IKM,IZTU1.VPOCHA,AIMPL,MLENT1.LECT,DT,
     &        MELEMC)
C     ******
         SEGDES MELEMC
         SEGDES IPT1,IPM1,IPM2
         IF (IDIM.EQ.3) THEN
            SEGDES IPM3
         END IF
         SEGDES MCHAML,MELVAL
         SEGSUP IZHR,IZFFM
 
 
C     On met a jour le nombre d'element total (quelquesoit son type)
         NUTOEL=NUTOEL+NBEL
      END DO
 
C==========Fin Boucle principale==========C
 
C     On desactive les segments
 
      SEGDES MELEME
      SEGDES MCHELM
 
C     Si IKR = 0 alors IZTGG1 (coef 1) est un CHPOINT
C     que l'on desactive.
      IF(IKR.EQ.0)THEN
         SEGDES IZTGG1
      ENDIF
 
C     On ecrit dans la table KONV la matrice elementaire
C     de l'operateur
      CALL ECMO(KIZX,'MATELM','MATRIK',MATRIK)
 
      SEGDES IMATRI,MATRIK
      RETURN
 
      END