Sources Esope | Cast3M

crectr
C CRECTR    SOURCE    PV        20/03/24    21:16:35     10554          
      SUBROUTINE CRECTR
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)
C*************************************************************************
C
C OBJET   : Cree un point au centre de gravite des éléments d'un maillage
C  (On ignore si ce point existe déja dans le cas d'éléments quadratiques)
C SYNTAXE : OBJ2 = KCTR OBJ1 &lt;'INCL' TABDOM> ;
C
C  OBJ1 : objet 'MAILLAGE'
C  OBJ2 : objet 'MAILLAGE' constitué d'éléments POI1
C
C*************************************************************************
-INC SMELEME
      POINTEUR MP1.MELEME
      POINTEUR MELEMC.MELEME,MELEF1.MELEME
-INC SMTABLE
      POINTEUR MTABD.MTABLE
-INC SMCOORD
-INC PPARAM
-INC CCOPTIO
      CHARACTER*4 LISMO(1)
      PARAMETER (NTB=1)
      DIMENSION   KTAB(NTB)
      CHARACTER*8 LTAB(NTB),TYPE
      DATA LISMO /'INCL'/
      DATA LTAB /'DOMAINE '/
 
C***
 
      CALL LIROBJ('MAILLAGE',MELEME,1,IRET)
      IF(IRET.EQ.0)RETURN
 
      KINC=0
      CALL LIRMOT(LISMO,1,IP,0)
      IF(IP.NE.0)THEN
         KINC=1
         NTO=1
         CALL LITABS(LTAB,KTAB,NTB,NTO,IRET)
         IF(IRET.EQ.0)RETURN
         CALL LIRREE(XVAL,1,IRET)
         IF(IRET.EQ.0)RETURN
         MTABD=KTAB(1)
         TYPE=' '
         CALL ACMO(MTABD,'CENTRE',TYPE,MELEMC)
C      call ecrobj('MAILLAGE',MELEMC)
         IF(TYPE.NE.'MAILLAGE')RETURN
         TYPE=' '
         CALL ACMO(MTABD,'FACE',TYPE,MELEF1)
      ENDIF
 
      CALL KNBEL(MELEME,NELN)
 
      segact mcoord*mod
      NBPTI=nbpts
      NBPTS=NBPTI+NELN
      SEGADJ MCOORD
 
      NBSOUS=0
      NBREF=0
      NBNN=1
      NBELEM=NELN
      SEGINI MP1
      MP1.ITYPEL=1
      SEGACT MELEME
      NBSOUS=LISOUS(/1)
      IF(NBSOUS.EQ.0)NBSOUS=1
      K0=0
 
      DO 1 L=1,NBSOUS
         IF(NBSOUS.NE.1)THEN
            IPT1=LISOUS(L)
            SEGACT IPT1
         ELSE
            IPT1=MELEME
         ENDIF
         NP=IPT1.NUM(/1)
         NEL=IPT1.NUM(/2)
 
         IF(IDIM.EQ.2)THEN
            DO 2 K=1,NEL
               NK=K0+K
               XC=0.D0
               YC=0.D0
               DO 21 I=1,NP
                  IP=IPT1.NUM(I,K)
                  XC=XC+XCOOR((IP-1)*(IDIM+1)+1)
                  YC=YC+XCOOR((IP-1)*(IDIM+1)+2)
 21            CONTINUE
               XC=XC/DBLE(NP)
               YC=YC/DBLE(NP)
               XD=XCOOR(IP*(IDIM+1))
 
               IP=NBPTI+NK
               XCOOR((IP-1)*(IDIM+1)+1)=XC
               XCOOR((IP-1)*(IDIM+1)+2)=YC
               XCOOR(IP*(IDIM+1))=XD
               MP1.NUM(1,NK)=IP
 
 2          CONTINUE
 
 
         ELSEIF(IDIM.EQ.3)THEN
            DO 3 K=1,NEL
               NK=K0+K
               XC=0.D0
               YC=0.D0
               ZC=0.D0
               DO 31 I=1,NP
                  IP=IPT1.NUM(I,K)
                  XC=XC+XCOOR((IP-1)*(IDIM+1)+1)
                  YC=YC+XCOOR((IP-1)*(IDIM+1)+2)
                  ZC=ZC+XCOOR((IP-1)*(IDIM+1)+3)
 31            CONTINUE
               XC=XC/DBLE(NP)
               YC=YC/DBLE(NP)
               ZC=ZC/DBLE(NP)
               XD=XCOOR(IP*(IDIM+1))
 
               IP=NBPTI+NK
               XCOOR((IP-1)*(IDIM+1)+1)=XC
               XCOOR((IP-1)*(IDIM+1)+2)=YC
               XCOOR((IP-1)*(IDIM+1)+3)=ZC
               XCOOR(IP*(IDIM+1))=XD
               MP1.NUM(1,NK)=IP
 
 3          CONTINUE
         ENDIF
         K0=K0+NEL
         IF(NBSOUS.NE.1)THEN
            SEGDES IPT1
         ENDIF
 1    CONTINUE
      SEGDES MELEME
      IF(KINC.NE.0)THEN
         WRITE(6,1951)NELN
 1951    FORMAT(1X,'KCTR : Creation des points centre :',
     &        ' Nombre de points a eliminer :',I7)
         CALL ECMO(MTABD,'BIDON','MAILLAGE',MP1)
         CALL ECRREE(XVAL)
         CALL ECROBJ('MAILLAGE',MELEMC)
         CALL ECROBJ('MAILLAGE',MP1   )
         CALL PRELIM(0)
         CALL LIROBJ('MAILLAGE',IP,1,IRET)
C     write(6,*)' Retour prelim melemc,mp1,ip=',melemc,mp1,ip
         IF(MELEF1.NE.0)THEN
            WRITE(6,*)' Elimination avec les faces'
            CALL ECRREE(XVAL)
            CALL ECROBJ('MAILLAGE',MELEF1)
            CALL ECROBJ('MAILLAGE',MP1)
            CALL PRELIM(0)
            CALL LIROBJ('MAILLAGE',IP,1,IRET)
C     write(6,*)' Retour prelim melemc,mp1,ip=',melemc,mp1,ip
         ENDIF
      ELSE
C Car normalement, ELIM a sans doute desactive MP1
         SEGDES MP1
      ENDIF
      CALL ECROBJ('MAILLAGE',MP1)
      RETURN
 
      END