C CAPA1     SOURCE    PV090527  26/04/30    21:15:11     12529          

C=======================================================================
C=                            C A P A 1                                =
C=                            ---------                                =
C=                                                                     =
C=  Fonction :                                                         =
C=  ----------                                                         =
C=  Calcul de la matrice de CAPACITE CALORIFIQUE (type RIGIDITE)       =
C=                                                                     =
C=  Parametres :  (E)=Entree   (S)=Sortie                              =
C=  ------------                                                       =
C=   IPMODE   (E)   Segment IMODEL pour un modele elementaire  (ACTIF) =
C=   IPCHEL   (E)   Segment MCHELM de CARACTERISTIQUES         (?)     =
C=   ISUPC    (E)   Support du champ de CARACTERISTIQUES               =
C=   ITABCP   (E)   TABLE pour le changement de PHASE                  =
C=   IPRIGI  (E/S)  Segment MRIGID : CAPACITE                  (ACTIF) =
C=                                                                     =
C=  Creation par Denis ROBERT le 15 fevrier 1988.                      =
C=  Modifications par DEGAY le 10 mai 1994 et ulterieurement.          =
C=======================================================================

      SUBROUTINE CAPA1 (IPMODE,IPCHEL,ISUPC,ITABCP, IPRIGI)

      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)


-INC PPARAM
-INC CCOPTIO
-INC CCHAMP
-INC CCREEL

-INC SMELEME
-INC SMMODEL
      POINTEUR NOMID1.NOMID,NOMID2.NOMID
-INC SMRIGID

      INTEGER OOOVAL

      SEGMENT NOTYPE
        CHARACTER*16 TYPE(NBTYPE)
      ENDSEGMENT

      CHARACTER*8  CMATE
      CHARACTER*16 MOFOR
      CHARACTER*(LCONMO) CONM

      PARAMETER ( NINF=3 )
      DIMENSION INFOS(NINF)

C=  LEFMAS   Liste des numeros d'elements finis MASSIFs implementes
C=  NEFMAS   Longueur de cette liste                                   =
C=  LEFCOQ   Liste des numeros d'elements finis COQUEs implementes
C=  NEFCOQ   Longueur de cette liste                                   =
      PARAMETER (NEFMAS=14, NEFCOQ=5,neftuy=2)
      DIMENSION LEFMAS(NEFMAS), LEFCOQ(NEFCOQ),leftuy(neftuy)

C Elements         TRI3 TRI6 QUA4 QUA8 CUB8 CU20 PRI6 PR15 TET4
C MASSIFs          TE10 PYR5 PY13 T1D2 T1D3
      DATA LEFMAS /   4,   6,   8,  10,  14,  15,  16,  17,  23,
     &               24,  25,  26, 191, 192 /
C COQUEs            COQ2  COQ3  COQ6  COQ4  COQ8
      DATA LEFCOQ /   44,   27,   56,   49,   41 /
C   Tuyau(pour advection)
C                    TUY2   TUY3
      DATA  leftuy/  269  , 270/

      IPINTE = 0
      IPINT2 = 0
      MOMATE = 0
      MOTYPE = 0
      IPCPHA = 0
      MOCPHA = 0

C- Matrice de capacite
      MRIGID = IPRIGI
c*      SEGACT,MRIGID
      NRIGE0 = IRIGEL(/2)

C- Recuperation du sous-modele et de la zone elementaire associee
      IMODEL = IPMODE
C*      SEGACT,IMODEL

      NEF = NEFMOD
      MFR = NUMMFR(NEF)
      IF (NEF.EQ.22) RETURN

      IMAS = 0
      ICOQ = 0
      ituy = 0
      CALL PLACE2(LEFMAS,NEFMAS,IMAS,NEF)
      CALL PLACE2(LEFCOQ,NEFCOQ,ICOQ,NEF)
      CALL PLACE2(LEFTUY,NEFTUY,ituy,NEF)
C- Recuperation d'informations sur le maillage elementaire
      IPT1 = IMAMOD
      SEGACT,IPT1
      NBNOE1 = IPT1.NUM(/1)
      NBELE1 = IPT1.NUM(/2)

C- Quelques informations sur le modele
      CONM  = CONMOD
      CMATE = CMATEE
      MATE  = IMATEE

      IRET = 1


      CALL IDENT(IPT1,CONM,IPCHEL,0, INFOS,IRET)
      IF (IRET.EQ.0) GOTO 9990

C- Recuperation d'informations sur l'element fini
      CALL TSHAPE(NEF,'GAUSS',IPINTE)
      IF (IERR.NE.0) GOTO 9990
      IF (ICOQ.NE.0) THEN
        IF (NEF.EQ.56 .OR. NEF.EQ.41 .OR. NEF.EQ.49) THEN
          CALL TSHAPE(NEF,'NOEUD',IPINT2)
        ENDIF
      ENDIF

C- Recuperation des caracteristiques materielles (obligatoires)
      nbrfac = 0
      nbrobl = 0
      MOFOR  = FORMOD(1)
      IF    (MOFOR .EQ. 'THERMIQUE') THEN
        INFOR = 1
        IF (NEF.EQ.46 .OR. ICOQ.NE.0 .or. ituy.ne.0) THEN
          nbrobl = 3
        ELSE
          nbrobl = 2
        ENDIF
        SEGINI,nomid

        IF(MFR .EQ. 75)THEN
C         Cas des JOI1 (MFR=75) ==> Ressorts THERMIQUES
C         ====================
          lesobl(1) = 'M   '
          lesobl(2) = 'C   '
        ELSE
          lesobl(1) = 'RHO '
          lesobl(2) = 'C   '
        ENDIF

      ELSEIF(MOFOR .EQ. 'DIFFUSION') THEN
        INFOR = 2
        IF (NEF.EQ.46 .OR. ICOQ.NE.0 .or. ituy.ne.0) THEN
          nbrobl = 2
        ELSE
          nbrobl = 1
        ENDIF
        SEGINI,nomid
        lesobl(1) = 'CDIF'
      ELSE
        CALL ERREUR(21)
        RETURN
      ENDIF

      IF (NEF.EQ.46.or.ituy.ne.0) THEN
        lesobl(nbrobl) = 'SECT'
      ELSE IF (ICOQ.NE.0) THEN
        lesobl(nbrobl) = 'EPAI'
      ENDIF

      NMATT   = nbrobl + nbrfac
      MOMATE  = nomid
c
      nbtype  = 1
      SEGINI,notype
      type(1) = 'REAL*8'
      MOTYPE  = notype

C- Recuperation de donnees dans le cas d'un CHANGEMENT DE PHASE
C  Dans le cas d'un changement de phase, on calcule une capacite
C  calorifique equivalente liee a la chaleur latente (MCHAML IPCPHA).
C  Le support de ce champ est celui des points de GAUSS (IPINTE).
      IF (ITABCP.NE.0) THEN
        CALL CAPA7(ITABCP,IPT1,ICOQ,IPINTE, IPCPHA)
        IF (IERR.NE.0) GOTO 9990
        nbrobl = 1
        nbrfac = 0
        SEGINI,nomid
        lesobl(1) = 'C   '
        MOCPHA = nomid
        NPHAT = nbrobl + nbrfac
      ENDIF

C- Definition du descripteur IDESCR pour le modele elementaire
      IDESCR = 0
      NOMPRI = LNOMID(1)
      NOMDUA = LNOMID(2)
      CALL TCOND2(ICOQ,NBNOE1,IDESCR,NOMPRI,NOMDUA)
      DESCR  = IDESCR
      SEGACT,DESCR
      NLIGRD = LISDUA(/2)
      NLIGRP = LISINC(/2)
      SEGDES,DESCR
      LRE    = NLIGRD

C- Partionnement si necessaire de la matrice de capacite
C- determinant ainsi le nombre d'objets elementaires de MRIGID
      LTRK = oooval(1,4)
      IF (LTRK.EQ.0) LTRK = oooval(1,1)
      LTRK=MAX(LTRK,2**24)
* Ajout a la taille en mots de la matrice des infos du segment
      LSEG = LRE*LRE*NBELE1 + 16
      NBLPRT = (LSEG-1)/LTRK + 1
      NBLMAX = (NBELE1-1)/NBLPRT + 1
      NBLPRT = (NBELE1-1)/NBLMAX + 1
*      write(ioimp,*) ' capa1 : nblprt nblmax = ',nblprt,nblmax,nbele1

C Ajout de la matrice de CAPACITE a la matrice globale
C ====================================================
      NRIGEL = NRIGE0 + NBLPRT
      SEGADJ,MRIGID

      descr  = IDESCR
      meleme = IPT1
      nbnn   = NBNOE1
      nbelem = NBELE1
      nbsous = 0
      nbref  = 0

C Boucle sur les PARTITIONS elementaires de la matrice
C ====================================================
      DO irige = 1, NBLPRT

        IVAMAT = 0
        IVAPHA = 0

        IF (NBLPRT.GT.1) THEN
C- Partitionnement du maillage support de la matrice elementaire
          SEGACT,IPT1
          ielem = (irige-1)*NBLMAX
          nbelem = MIN(NBLMAX,NBELE1-ielem)
*          write(ioimp,*) ' creation segment ',nbnn,nbelem
          SEGINI,meleme
          itypel = IPT1.itypel
          DO ielt = 1, nbelem
            jelt = ielt + ielem
            DO inoe = 1, nbnn
              num(inoe,ielt) = IPT1.NUM(inoe,jelt)
            ENDDO
            icolor(ielt) = IPT1.ICOLOR(jelt)
          ENDDO
C- Recopie du descripteur
          des1 = IDESCR
          SEGINI,descr=des1
          SEGDES,descr
        ENDIF
        ipmail = meleme
        ipdesc = descr

C- Initialisation de la matrice de rigidite elementaire (xmatri)
        NELRIG = nbelem
        rigrel=0
        SEGINI,xmatri
        ipmatr = xmatri

C- Recuperation des valeurs de caracteristiques sur la partition
        CALL KOMCHA(IPCHEL,ipmail,CONM,MOMATE,MOTYPE,1,INFOS,3,IVAMAT)
        IF (IERR.NE.0) GOTO 9991
        IF (ISUPC.EQ.1 .AND. NEF.NE.265) THEN
C           On ne change pas le support pour JOI1
          CALL VALCHE(IVAMAT,NMATT,IPINTE,0,MOMATE,NEF)
          IF (IERR.NE.0) THEN
            ISUPC = 0
            GOTO 9991
          ENDIF
        ENDIF

C- Idem pour capacite equivalente en cas de changement de phase
        IF (ITABCP.NE.0) THEN
          CALL KOMCHA(IPCPHA,ipmail,CONM,MOCPHA,MOTYPE,1,INFOS,3,IVAPHA)
          IF (IERR.NE.0) GOTO 9991
        ENDIF

C- Calcul de la matrice elementaire pour la paritition elementaire et
C  Remplissage de la matrice globale (IPRIGI)
C---> Elements MASSIFs a integration NUMERIQUE
        IF (IMAS.NE.0) THEN
          CALL CAPANU(NEF,ipmail,IPINTE,IVAMAT,NMATT,IVAPHA,NPHAT,
     &                ipmatr,LRE,INFOR)
C---> Elements de COQUEs
        ELSE IF (ICOQ.NE.0) THEN
C-----> Element COQ2 (axisymetrique)
          IF (NEF.EQ.44) THEN
            CALL CAPAC1(NEF,ipmail,IPINTE,IVAMAT,NMATT,IVAPHA,NPHAT,
     &                  ipmatr,LRE,INFOR)
C-----> Element COQ3
          ELSE IF (NEF.EQ.27) THEN
            CALL CAPAC3(NEF,ipmail,IPINTE,IVAMAT,NMATT,IVAPHA,NPHAT,
     &                  ipmatr,LRE,INFOR)
C-----> Elements COQ4, COQ6 et COQ8
C*          ELSE IF (NEF.EQ.56 .OR. NEF.EQ.41 .OR. NEF.EQ.49) THEN
          ELSE
            CALL CAPAC2(NEF,ipmail,IPINTE,IPINT2,IVAMAT,NMATT,
     &                  IVAPHA,NPHAT, ipmatr,LRE,INFOR)
          ENDIF
C---> Element BARRE (integration NUMERIQUE)
        ELSE IF (NEF.EQ.46.or.nef.eq.269.or.nef.eq.270) THEN
          CALL CAPABA(NEF,ipmail,IPINTE,IVAMAT,NMATT,IVAPHA,NPHAT,
     &                ipmatr,LRE,INFOR)
C---> Elements seg3, RAC2 et RAC3 : non implementes
        ELSE IF (NEF.EQ.3 .OR. NEF.EQ.12 .OR. NEF.EQ.13) THEN
          CALL ERREUR(251)
          RETURN
C---> Elements JOI1 : pas d'integration
        ELSE IF (NEF.EQ.265) THEN
          CALL CAPAJ1(IPMAIL,IVAMAT,NMATT,IPMATR,INFOR)
C---> Elements POI1 : pas d'integration
        ELSE IF (NEF.EQ.45) THEN
          CALL CAPAP1(IPMAIL,IVAMAT,IPMATR,INFOR)
C---> Option indisponible : ERREUR
        ELSE
          CALL ERREUR(19)
          RETURN
        ENDIF

C- Un peu de menage dans les segments
 9991   CONTINUE
        IF (ISUPC.EQ.1 .AND. NEF.NE.265 .OR. NBLPRT.NE.1) THEN
          CALL DTMVAL(IVAMAT,3)
        ELSE
          CALL DTMVAL(IVAMAT,1)
        ENDIF
        IF (ITABCP.NE.0) THEN
          IF (NBLPRT.NE.1) THEN
            CALL DTMVAL(IVAPHA,3)
          ELSE
            CALL DTMVAL(IVAPHA,1)
          ENDIF
        ENDIF
C- Sortie prematuree en cas d'erreur
        IF (IERR.NE.0) GOTO 9990

        xmatri = ipmatr
        SEGDES,xmatri
        IF (NBLPRT.NE.1) THEN
          meleme = ipmail
        ENDIF

C- Stockage de la matrice de CAPACITE du modele
        jrige = NRIGE0 + irige
        COERIG(jrige)   = 1.
        IRIGEL(1,jrige) = ipmail
        IRIGEL(2,jrige) = 0
        IRIGEL(3,jrige) = ipdesc
        IRIGEL(4,jrige) = ipmatr
        IRIGEL(5,jrige) = NIFOUR
        IRIGEL(6,jrige) = 0
        IRIGEL(7,jrige) = 0
        IRIGEL(8,jrige) = 0

      ENDDO

      IPRIGI = MRIGID

C MENAGE : DESACTIVATION/DESTRUCTION DE SEGMENTS
C ==============================================
 9990 CONTINUE
      IF (MOMATE.NE.0) THEN
        nomid = MOMATE
        SEGSUP,nomid
      ENDIF
      IF (MOTYPE.NE.0) THEN
        notype = MOTYPE
        SEGSUP,notype
      ENDIF
      IF (MOCPHA.NE.0) THEN
        nomid = MOCPHA
        SEGSUP,nomid
      ENDIF
      IF (IPCPHA.NE.0) THEN
        CALL DTCHAM(IPCPHA)
      ENDIF

      RETURN
      END