Télécharger la plaquette Cast3M


fr en

Télécharger capa1.eso

Retour à la liste

Numérotation des lignes : 

capa1

  1. C CAPA1     SOURCE    PV090527  26/04/30    21:15:11     12529          
  2.  
  3. C=======================================================================
  4. C=                            C A P A 1                                =
  5. C=                            ---------                                =
  6. C=                                                                     =
  7. C=  Fonction :                                                         =
  8. C=  ----------                                                         =
  9. C=  Calcul de la matrice de CAPACITE CALORIFIQUE (type RIGIDITE)       =
  10. C=                                                                     =
  11. C=  Parametres :  (E)=Entree   (S)=Sortie                              =
  12. C=  ------------                                                       =
  13. C=   IPMODE   (E)   Segment IMODEL pour un modele elementaire  (ACTIF) =
  14. C=   IPCHEL   (E)   Segment MCHELM de CARACTERISTIQUES         (?)     =
  15. C=   ISUPC    (E)   Support du champ de CARACTERISTIQUES               =
  16. C=   ITABCP   (E)   TABLE pour le changement de PHASE                  =
  17. C=   IPRIGI  (E/S)  Segment MRIGID : CAPACITE                  (ACTIF) =
  18. C=                                                                     =
  19. C=  Creation par Denis ROBERT le 15 fevrier 1988.                      =
  20. C=  Modifications par DEGAY le 10 mai 1994 et ulterieurement.          =
  21. C=======================================================================
  22.  
  23.       SUBROUTINE CAPA1 (IPMODE,IPCHEL,ISUPC,ITABCP, IPRIGI)
  24.  
  25.       IMPLICIT INTEGER(I-N)
  26.       IMPLICIT REAL*8 (A-H,O-Z)
  27.  
  28.  
  29. -INC PPARAM
  30. -INC CCOPTIO
  31. -INC CCHAMP
  32. -INC CCREEL
  33.  
  34. -INC SMELEME
  35. -INC SMMODEL
  36.       POINTEUR NOMID1.NOMID,NOMID2.NOMID
  37. -INC SMRIGID
  38.  
  39.       INTEGER OOOVAL
  40.  
  41.       SEGMENT NOTYPE
  42.         CHARACTER*16 TYPE(NBTYPE)
  43.       ENDSEGMENT
  44.  
  45.       CHARACTER*8  CMATE
  46.       CHARACTER*16 MOFOR
  47.       CHARACTER*(LCONMO) CONM
  48.  
  49.       PARAMETER ( NINF=3 )
  50.       DIMENSION INFOS(NINF)
  51.  
  52. C=  LEFMAS   Liste des numeros d'elements finis MASSIFs implementes
  53. C=  NEFMAS   Longueur de cette liste                                   =
  54. C=  LEFCOQ   Liste des numeros d'elements finis COQUEs implementes
  55. C=  NEFCOQ   Longueur de cette liste                                   =
  56.       PARAMETER (NEFMAS=14, NEFCOQ=5,neftuy=2)
  57.       DIMENSION LEFMAS(NEFMAS), LEFCOQ(NEFCOQ),leftuy(neftuy)
  58.  
  59. C Elements         TRI3 TRI6 QUA4 QUA8 CUB8 CU20 PRI6 PR15 TET4
  60. C MASSIFs          TE10 PYR5 PY13 T1D2 T1D3
  61.       DATA LEFMAS /   4,   6,   8,  10,  14,  15,  16,  17,  23,
  62.      &               24,  25,  26, 191, 192 /
  63. C COQUEs            COQ2  COQ3  COQ6  COQ4  COQ8
  64.       DATA LEFCOQ /   44,   27,   56,   49,   41 /
  65. C   Tuyau(pour advection)
  66. C                    TUY2   TUY3
  67.       DATA  leftuy/  269  , 270/
  68.  
  69.       IPINTE = 0
  70.       IPINT2 = 0
  71.       MOMATE = 0
  72.       MOTYPE = 0
  73.       IPCPHA = 0
  74.       MOCPHA = 0
  75.  
  76. C- Matrice de capacite
  77.       MRIGID = IPRIGI
  78. c*      SEGACT,MRIGID
  79.       NRIGE0 = IRIGEL(/2)
  80.  
  81. C- Recuperation du sous-modele et de la zone elementaire associee
  82.       IMODEL = IPMODE
  83. C*      SEGACT,IMODEL
  84.  
  85.       NEF = NEFMOD
  86.       MFR = NUMMFR(NEF)
  87.       IF (NEF.EQ.22) RETURN
  88.  
  89.       IMAS = 0
  90.       ICOQ = 0
  91.       ituy = 0
  92.       CALL PLACE2(LEFMAS,NEFMAS,IMAS,NEF)
  93.       CALL PLACE2(LEFCOQ,NEFCOQ,ICOQ,NEF)
  94.       CALL PLACE2(LEFTUY,NEFTUY,ituy,NEF)
  95. C- Recuperation d'informations sur le maillage elementaire
  96.       IPT1 = IMAMOD
  97.       SEGACT,IPT1
  98.       NBNOE1 = IPT1.NUM(/1)
  99.       NBELE1 = IPT1.NUM(/2)
  100.  
  101. C- Quelques informations sur le modele
  102.       CONM  = CONMOD
  103.       CMATE = CMATEE
  104.       MATE  = IMATEE
  105.  
  106.       IRET = 1
  107.  
  108.  
  109.       CALL IDENT(IPT1,CONM,IPCHEL,0, INFOS,IRET)
  110.       IF (IRET.EQ.0) GOTO 9990
  111.  
  112. C- Recuperation d'informations sur l'element fini
  113.       CALL TSHAPE(NEF,'GAUSS',IPINTE)
  114.       IF (IERR.NE.0) GOTO 9990
  115.       IF (ICOQ.NE.0) THEN
  116.         IF (NEF.EQ.56 .OR. NEF.EQ.41 .OR. NEF.EQ.49) THEN
  117.           CALL TSHAPE(NEF,'NOEUD',IPINT2)
  118.         ENDIF
  119.       ENDIF
  120.  
  121. C- Recuperation des caracteristiques materielles (obligatoires)
  122.       nbrfac = 0
  123.       nbrobl = 0
  124.       MOFOR  = FORMOD(1)
  125.       IF    (MOFOR .EQ. 'THERMIQUE') THEN
  126.         INFOR = 1
  127.         IF (NEF.EQ.46 .OR. ICOQ.NE.0 .or. ituy.ne.0) THEN
  128.           nbrobl = 3
  129.         ELSE
  130.           nbrobl = 2
  131.         ENDIF
  132.         SEGINI,nomid
  133.  
  134.         IF(MFR .EQ. 75)THEN
  135. C         Cas des JOI1 (MFR=75) ==> Ressorts THERMIQUES
  136. C         ====================
  137.           lesobl(1) = 'M   '
  138.           lesobl(2) = 'C   '
  139.         ELSE
  140.           lesobl(1) = 'RHO '
  141.           lesobl(2) = 'C   '
  142.         ENDIF
  143.  
  144.       ELSEIF(MOFOR .EQ. 'DIFFUSION') THEN
  145.         INFOR = 2
  146.         IF (NEF.EQ.46 .OR. ICOQ.NE.0 .or. ituy.ne.0) THEN
  147.           nbrobl = 2
  148.         ELSE
  149.           nbrobl = 1
  150.         ENDIF
  151.         SEGINI,nomid
  152.         lesobl(1) = 'CDIF'
  153.       ELSE
  154.         CALL ERREUR(21)
  155.         RETURN
  156.       ENDIF
  157.  
  158.       IF (NEF.EQ.46.or.ituy.ne.0) THEN
  159.         lesobl(nbrobl) = 'SECT'
  160.       ELSE IF (ICOQ.NE.0) THEN
  161.         lesobl(nbrobl) = 'EPAI'
  162.       ENDIF
  163.  
  164.       NMATT   = nbrobl + nbrfac
  165.       MOMATE  = nomid
  166. c
  167.       nbtype  = 1
  168.       SEGINI,notype
  169.       type(1) = 'REAL*8'
  170.       MOTYPE  = notype
  171.  
  172. C- Recuperation de donnees dans le cas d'un CHANGEMENT DE PHASE
  173. C  Dans le cas d'un changement de phase, on calcule une capacite
  174. C  calorifique equivalente liee a la chaleur latente (MCHAML IPCPHA).
  175. C  Le support de ce champ est celui des points de GAUSS (IPINTE).
  176.       IF (ITABCP.NE.0) THEN
  177.         CALL CAPA7(ITABCP,IPT1,ICOQ,IPINTE, IPCPHA)
  178.         IF (IERR.NE.0) GOTO 9990
  179.         nbrobl = 1
  180.         nbrfac = 0
  181.         SEGINI,nomid
  182.         lesobl(1) = 'C   '
  183.         MOCPHA = nomid
  184.         NPHAT = nbrobl + nbrfac
  185.       ENDIF
  186.  
  187. C- Definition du descripteur IDESCR pour le modele elementaire
  188.       IDESCR = 0
  189.       NOMPRI = LNOMID(1)
  190.       NOMDUA = LNOMID(2)
  191.       CALL TCOND2(ICOQ,NBNOE1,IDESCR,NOMPRI,NOMDUA)
  192.       DESCR  = IDESCR
  193.       SEGACT,DESCR
  194.       NLIGRD = LISDUA(/2)
  195.       NLIGRP = LISINC(/2)
  196.       SEGDES,DESCR
  197.       LRE    = NLIGRD
  198.  
  199. C- Partionnement si necessaire de la matrice de capacite
  200. C- determinant ainsi le nombre d'objets elementaires de MRIGID
  201.       LTRK = oooval(1,4)
  202.       IF (LTRK.EQ.0) LTRK = oooval(1,1)
  203.       LTRK=MAX(LTRK,2**24)
  204. * Ajout a la taille en mots de la matrice des infos du segment
  205.       LSEG = LRE*LRE*NBELE1 + 16
  206.       NBLPRT = (LSEG-1)/LTRK + 1
  207.       NBLMAX = (NBELE1-1)/NBLPRT + 1
  208.       NBLPRT = (NBELE1-1)/NBLMAX + 1
  209. *      write(ioimp,*) ' capa1 : nblprt nblmax = ',nblprt,nblmax,nbele1
  210.  
  211. C Ajout de la matrice de CAPACITE a la matrice globale
  212. C ====================================================
  213.       NRIGEL = NRIGE0 + NBLPRT
  214.       SEGADJ,MRIGID
  215.  
  216.       descr  = IDESCR
  217.       meleme = IPT1
  218.       nbnn   = NBNOE1
  219.       nbelem = NBELE1
  220.       nbsous = 0
  221.       nbref  = 0
  222.  
  223. C Boucle sur les PARTITIONS elementaires de la matrice
  224. C ====================================================
  225.       DO irige = 1, NBLPRT
  226.  
  227.         IVAMAT = 0
  228.         IVAPHA = 0
  229.  
  230.         IF (NBLPRT.GT.1) THEN
  231. C- Partitionnement du maillage support de la matrice elementaire
  232.           SEGACT,IPT1
  233.           ielem = (irige-1)*NBLMAX
  234.           nbelem = MIN(NBLMAX,NBELE1-ielem)
  235. *          write(ioimp,*) ' creation segment ',nbnn,nbelem
  236.           SEGINI,meleme
  237.           itypel = IPT1.itypel
  238.           DO ielt = 1, nbelem
  239.             jelt = ielt + ielem
  240.             DO inoe = 1, nbnn
  241.               num(inoe,ielt) = IPT1.NUM(inoe,jelt)
  242.             ENDDO
  243.             icolor(ielt) = IPT1.ICOLOR(jelt)
  244.           ENDDO
  245. C- Recopie du descripteur
  246.           des1 = IDESCR
  247.           SEGINI,descr=des1
  248.           SEGDES,descr
  249.         ENDIF
  250.         ipmail = meleme
  251.         ipdesc = descr
  252.  
  253. C- Initialisation de la matrice de rigidite elementaire (xmatri)
  254.         NELRIG = nbelem
  255.         rigrel=0
  256.         SEGINI,xmatri
  257.         ipmatr = xmatri
  258.  
  259. C- Recuperation des valeurs de caracteristiques sur la partition
  260.         CALL KOMCHA(IPCHEL,ipmail,CONM,MOMATE,MOTYPE,1,INFOS,3,IVAMAT)
  261.         IF (IERR.NE.0) GOTO 9991
  262.         IF (ISUPC.EQ.1 .AND. NEF.NE.265) THEN
  263. C           On ne change pas le support pour JOI1
  264.           CALL VALCHE(IVAMAT,NMATT,IPINTE,0,MOMATE,NEF)
  265.           IF (IERR.NE.0) THEN
  266.             ISUPC = 0
  267.             GOTO 9991
  268.           ENDIF
  269.         ENDIF
  270.  
  271. C- Idem pour capacite equivalente en cas de changement de phase
  272.         IF (ITABCP.NE.0) THEN
  273.           CALL KOMCHA(IPCPHA,ipmail,CONM,MOCPHA,MOTYPE,1,INFOS,3,IVAPHA)
  274.           IF (IERR.NE.0) GOTO 9991
  275.         ENDIF
  276.  
  277. C- Calcul de la matrice elementaire pour la paritition elementaire et
  278. C  Remplissage de la matrice globale (IPRIGI)
  279. C---> Elements MASSIFs a integration NUMERIQUE
  280.         IF (IMAS.NE.0) THEN
  281.           CALL CAPANU(NEF,ipmail,IPINTE,IVAMAT,NMATT,IVAPHA,NPHAT,
  282.      &                ipmatr,LRE,INFOR)
  283. C---> Elements de COQUEs
  284.         ELSE IF (ICOQ.NE.0) THEN
  285. C-----> Element COQ2 (axisymetrique)
  286.           IF (NEF.EQ.44) THEN
  287.             CALL CAPAC1(NEF,ipmail,IPINTE,IVAMAT,NMATT,IVAPHA,NPHAT,
  288.      &                  ipmatr,LRE,INFOR)
  289. C-----> Element COQ3
  290.           ELSE IF (NEF.EQ.27) THEN
  291.             CALL CAPAC3(NEF,ipmail,IPINTE,IVAMAT,NMATT,IVAPHA,NPHAT,
  292.      &                  ipmatr,LRE,INFOR)
  293. C-----> Elements COQ4, COQ6 et COQ8
  294. C*          ELSE IF (NEF.EQ.56 .OR. NEF.EQ.41 .OR. NEF.EQ.49) THEN
  295.           ELSE
  296.             CALL CAPAC2(NEF,ipmail,IPINTE,IPINT2,IVAMAT,NMATT,
  297.      &                  IVAPHA,NPHAT, ipmatr,LRE,INFOR)
  298.           ENDIF
  299. C---> Element BARRE (integration NUMERIQUE)
  300.         ELSE IF (NEF.EQ.46.or.nef.eq.269.or.nef.eq.270) THEN
  301.           CALL CAPABA(NEF,ipmail,IPINTE,IVAMAT,NMATT,IVAPHA,NPHAT,
  302.      &                ipmatr,LRE,INFOR)
  303. C---> Elements seg3, RAC2 et RAC3 : non implementes
  304.         ELSE IF (NEF.EQ.3 .OR. NEF.EQ.12 .OR. NEF.EQ.13) THEN
  305.           CALL ERREUR(251)
  306.           RETURN
  307. C---> Elements JOI1 : pas d'integration
  308.         ELSE IF (NEF.EQ.265) THEN
  309.           CALL CAPAJ1(IPMAIL,IVAMAT,NMATT,IPMATR,INFOR)
  310. C---> Elements POI1 : pas d'integration
  311.         ELSE IF (NEF.EQ.45) THEN
  312.           CALL CAPAP1(IPMAIL,IVAMAT,IPMATR,INFOR)
  313. C---> Option indisponible : ERREUR
  314.         ELSE
  315.           CALL ERREUR(19)
  316.           RETURN
  317.         ENDIF
  318.  
  319. C- Un peu de menage dans les segments
  320.  9991   CONTINUE
  321.         IF (ISUPC.EQ.1 .AND. NEF.NE.265 .OR. NBLPRT.NE.1) THEN
  322.           CALL DTMVAL(IVAMAT,3)
  323.         ELSE
  324.           CALL DTMVAL(IVAMAT,1)
  325.         ENDIF
  326.         IF (ITABCP.NE.0) THEN
  327.           IF (NBLPRT.NE.1) THEN
  328.             CALL DTMVAL(IVAPHA,3)
  329.           ELSE
  330.             CALL DTMVAL(IVAPHA,1)
  331.           ENDIF
  332.         ENDIF
  333. C- Sortie prematuree en cas d'erreur
  334.         IF (IERR.NE.0) GOTO 9990
  335.  
  336.         xmatri = ipmatr
  337.         SEGDES,xmatri
  338.         IF (NBLPRT.NE.1) THEN
  339.           meleme = ipmail
  340.         ENDIF
  341.  
  342. C- Stockage de la matrice de CAPACITE du modele
  343.         jrige = NRIGE0 + irige
  344.         COERIG(jrige)   = 1.
  345.         IRIGEL(1,jrige) = ipmail
  346.         IRIGEL(2,jrige) = 0
  347.         IRIGEL(3,jrige) = ipdesc
  348.         IRIGEL(4,jrige) = ipmatr
  349.         IRIGEL(5,jrige) = NIFOUR
  350.         IRIGEL(6,jrige) = 0
  351.         IRIGEL(7,jrige) = 0
  352.         IRIGEL(8,jrige) = 0
  353.  
  354.       ENDDO
  355.  
  356.       IPRIGI = MRIGID
  357.  
  358. C MENAGE : DESACTIVATION/DESTRUCTION DE SEGMENTS
  359. C ==============================================
  360.  9990 CONTINUE
  361.       IF (MOMATE.NE.0) THEN
  362.         nomid = MOMATE
  363.         SEGSUP,nomid
  364.       ENDIF
  365.       IF (MOTYPE.NE.0) THEN
  366.         notype = MOTYPE
  367.         SEGSUP,notype
  368.       ENDIF
  369.       IF (MOCPHA.NE.0) THEN
  370.         nomid = MOCPHA
  371.         SEGSUP,nomid
  372.       ENDIF
  373.       IF (IPCPHA.NE.0) THEN
  374.         CALL DTCHAM(IPCPHA)
  375.       ENDIF
  376.  
  377.       RETURN
  378.       END
  379.  
  380.  
  381.  
  382.  
  383.  

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