cgrhead.F

Go to the documentation of this file.

Copyright>        OpenRadioss
Copyright>        Copyright (C) 1986-2025 Altair Engineering Inc.
Copyright>
Copyright>        This program is free software: you can redistribute it and/or modify
Copyright>        it under the terms of the GNU Affero General Public License as published by
Copyright>        the Free Software Foundation, either version 3 of the License, or
Copyright>        (at your option) any later version.
Copyright>
Copyright>        This program is distributed in the hope that it will be useful,
Copyright>        but WITHOUT ANY WARRANTY; without even the implied warranty of
Copyright>        MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
Copyright>        GNU Affero General Public License for more details.
Copyright>
Copyright>        You should have received a copy of the GNU Affero General Public License
Copyright>        along with this program.  If not, see <https://www.gnu.org/licenses/>.
Copyright>
Copyright>
Copyright>        Commercial Alternative: Altair Radioss Software
Copyright>
Copyright>        As an alternative to this open-source version, Altair also offers Altair Radioss
Copyright>        software under a commercial license.  Contact Altair to discuss further if the
Copyright>        commercial version may interest you: https://www.altair.com/radioss/.
!||====================================================================
!||    cgrhead            ../starter/source/elements/shell/coque/cgrhead.F
!||--- called by ------------------------------------------------------
!||    lectur             ../starter/source/starter/lectur.F
!||--- calls      -----------------------------------------------------
!||    ancmsg             ../starter/source/output/message/message.F
!||--- uses       -----------------------------------------------------
!||    drape_mod          ../starter/share/modules1/drape_mod.F
!||    message_mod        ../starter/share/message_module/message_mod.F
!||    r2r_mod            ../starter/share/modules1/r2r_mod.F
!||    reorder_mod        ../starter/share/modules1/reorder_mod.F
!||    stack_mod          ../starter/share/modules1/stack_mod.F
!||====================================================================
      SUBROUTINE cgrhead(
     1       IXC     ,PM      ,GEO       ,INUM   ,ISEL   ,
     2       ITR1    ,EADD    ,INDEX     ,ITRI   ,XNUM   ,
     3       IPARTC  ,ND      ,THK       ,IGRSURF,IGRSH4N,
     4       CEP     ,XEP     ,IGEO      ,IPM    ,
     5       IPART   ,SH4TREE ,NOD2ELC   ,ISHEOFF,SH4TRIM,
     6       TAGPRT_SMS, LGAUGE,IWORKSH  ,MAT_PARAM,
     7       STACK   ,DRAPE  ,RNOISE ,SH4ANG,DRAPEG, PTSHEL,
     8       DAMP_RANGE_PART)
C-----------------------------------------------
C            A R G U M E N T S
C-----------------------------------------------
C     IXC(NIXC,NUMELC)  ARRAY MID(1)+CONECS(2-5)+PID(6)+           E
C                       N GLOBAL(7)                               E
C     PM(NPROPM,NUMMAT) ARRAY MATERIAL CHARACTERISTICS             E
C     GEO(NPROPG,NUMGEO)ARRAY PID CHARACTERISTICS                  E
C     INUM(9,NUMELC)    WORKING ARRAY                              E/S
C     ISEL(NSELC)       ARRAY SELECTED SHELLS FOR TH               E/S
C     ITR1(NSELC)       WORKING ARRAY                              E/S
C     EADD(NUMELC)      ARRAY ADDRESSES IN IDAM CHECKBOARD         S
C     INDEX(NUMELC)     WORKING ARRAY                              E/S
C     ITRI(7,NUMELC)    WORKING ARRAY                              E/S
C     IPARTC(NUMELC)    PART ARRAY                                 E/S
C     CEP(NUMELC)    PROC ARRAY                                    E/S
C     XEP(NUMELC)    PROC ARRAY                                    E/S
C     NOD2ELC(4*NUMELC)                                            E/S
C     ISHEOFF(NUMELC)   FLAG ELEM RBY ON/OFF                       E/S
C-----------------------------------------------
C   M o d u l e s
C-----------------------------------------------
      USE my_alloc_mod
      USE message_mod
      USE r2r_mod
      USE stack_mod
      USE message_mod
      USE reorder_mod
      USE groupdef_mod
      USE drape_mod
      USE matparam_def_mod
C-----------------------------------------------
C   I M P L I C I T   T Y P E S
C-----------------------------------------------
#include      "implicit_f.inc"
C-----------------------------------------------
C   C O M M O N   B L O C K S
C-----------------------------------------------
#include      "vect01_c.inc"
#include      "com04_c.inc"
#include      "com_xfem1.inc"
#include      "param_c.inc"
#include      "remesh_c.inc"
#include      "sms_c.inc"
#include      "scr17_c.inc"
#include      "r2r_c.inc"
#include      "drape_c.inc"
#include      "com01_c.inc"
C-----------------------------------------------
C   D U M M Y   A R G U M E N T S
C-----------------------------------------------
      INTEGER IXC(NIXC,*),ISEL(*),INUM(9,*),IPARTC(*), ISHEOFF(*),
     .        EADD(*),ITR1(*),INDEX(*),ITRI(8,*),
     .        ND, CEP(*), XEP(*),
     .        IGEO(NPROPGI,*),IPM(NPROPMI,*),IPART(LIPART1,*), 
     .        SH4TREE(KSH4TREE,*), NOD2ELC(*), SH4TRIM(*),
     .        TAGPRT_SMS(*) ,LGAUGE(3,*),
     .        IWORKSH(3,*)
      INTEGER , DIMENSION(NUMELC) , INTENT(INOUT):: PTSHEL
      INTEGER , INTENT(IN) :: DAMP_RANGE_PART(NPART)     !< flag to compute the damping range
C     REAL OR REAL*8
      my_real
     .    PM(NPROPM,*), GEO(NPROPG,*),XNUM(*),THK(*), RNOISE(NPERTURB,*),
     .    SH4ANG(*)
C-----------------------------------------------
      TYPE (STACK_PLY) :: STACK
      TYPE (DRAPE_)   , TARGET    :: DRAPE (NUMELC_DRAPE + NUMELTG_DRAPE)
      TYPE (DRAPEG_)                             :: DRAPEG
      TYPE (DRAPE_)  , DIMENSION(:), ALLOCATABLE :: XNUM_DRAPE
      TYPE (DRAPEG_)                             :: XNUM_DRAPEG
      TYPE (MATPARAM_STRUCT_) ,DIMENSION(NUMMAT),INTENT(IN) :: MAT_PARAM
C-----------------------------------------------
      TYPE (GROUP_)  , DIMENSION(NGRSHEL)  :: IGRSH4N
      TYPE (SURF_)   , DIMENSION(NSURF)   :: IGRSURF
C-----------------------------------------------
C   L O C A L   V A R I A B L E S
C-----------------------------------------------
      INTEGER, DIMENSION(:,:), ALLOCATABLE :: ISTOR,INUM_DRAPE
      INTEGER WORK(70000)
      INTEGER I, K,  MLN, NG, ISSN, NPN, IFIO, NN,L,IGTYP,
     .  mln0, issn0, ic, n, mid, mid0, pid, pid0, istr0,
     .  ihbe, ihbe0, ii, j, midn, pidn, nsg, nel, ne1,
     .  ithk, ithk0, ipla, ipla0,ii1,jj1,ii2,jj2,jj,ii3,jj3,ngrou,
     .  mskmln,msknpn,mskihb,mskisn,mskirb,mode,icsen,irb,
     .  mskist,mskipl,mskith,mskmid,mskpid,mskirp,msktyp,irep,
     .  ipt,imatly,ii0,jj0,ilev,prt,iadm,dir,ii4,jj4,n1,
     .  nfail,ifail,ixfem,inum_r2r(1+r2r_siu*numelc),
     .  ii5,jj5,ii6,jj6,
     .  isubstack,iigeo,iadi ,ippid,nb_law58,ipmat,
     .  ipert,stat,ip,nslice,kk,npt_drp,ie,ie0
      INTEGER, DIMENSION(:,:), ALLOCATABLE :: INUM_WORKC !(3,NUMELC)
      my_real, DIMENSION(:), ALLOCATABLE :: ANGLE !(NUMELC)
      EXTERNAL MY_SHIFTL,MY_SHIFTR,MY_AND
      INTEGER MY_SHIFTL,MY_SHIFTR,MY_AND,IPIDL
      INTEGER, DIMENSION(:), ALLOCATABLE :: INDEX2, INUM_PTSHEL
C     REAL OR REAL*8
      my_real, DIMENSION(:,:), ALLOCATABLE :: XNUM_RNOISE
C
      TYPE (DRAPE_PLY_)             , POINTER   :: DRAPE_PLY
C---+----1----+----2----+----3----+----4----+----5----+----6----+----7--
C1---------------------------------
      DATA mskmln /o'07770000000'/
      DATA msktyp /o'00007770000'/
      DATA mskihb /o'00000007000'/
      DATA mskisn /o'00000000700'/
      DATA mskist /o'00000000070'/
      DATA mskipl /o'00000000007'/
C2---------------------------------
      DATA mskith /o'10000000000'/
      DATA mskirp /o'07000000000'/
      DATA msknpn /o'00777000000'/
      DATA mskirb /o'00000000007'/
C3---------------------------------
      DATA mskmid /o'07777777777'/
C4---------------------------------
      DATA mskpid /o'07777777777'/
C======================================================================|
C   GLOBAL SORTING ON ALL CRITERIA FOR ALL ELEMENTS
C----------------------------------------------------------
      ALLOCATE(angle(numelc))
      ALLOCATE(inum_workc(3,numelc))
      IF(nadmesh /= 0)THEN
        ALLOCATE( istor(ksh4tree+1,numelc) )
      ELSE
        ALLOCATE( istor(0,0) )
      ENDIF
        IF (ndrape > 0 .AND. numelc_drape > 0) THEN
         ALLOCATE(xnum_drape(numelc))  
         ALLOCATE(xnum_drapeg%INDX(numelc))
         xnum_drapeg%INDX = 0                                             
         DO i =1, numelc 
           ie = drapeg%INDX(i)
           IF(ie == 0) cycle                                                                      
           npt_drp = drape(ie)%NPLY_DRAPE  
           npt     = drape(ie)%NPLY                                             
           ALLOCATE(xnum_drape(i)%INDX_PLY(npt))                                              
           ALLOCATE(xnum_drape(i)%DRAPE_PLY(npt_drp))        
           xnum_drape(i)%INDX_PLY = 0
           xnum_drape(i)%INDX_PLY = 0
           DO j = 1,npt_drp
                nslice = drape(ie)%DRAPE_PLY(j)%NSLICE    
                ALLOCATE(xnum_drape(i)%DRAPE_PLY(j)%RDRAPE(nslice,2))    
                ALLOCATE(xnum_drape(i)%DRAPE_PLY(j)%IDRAPE(nslice,2)) 
                xnum_drape(i)%DRAPE_PLY(j)%RDRAPE = 0             
                xnum_drape(i)%DRAPE_PLY(j)%IDRAPE = 0  
           ENDDO 
         ENDDO      
        ELSE
            ALLOCATE( xnum_drape(0) )
        ENDIF
      IF(abs(isigi) == 3 .OR. abs(isigi) == 4 .OR. abs(isigi) == 5) THEN
        ALLOCATE(inum_ptshel(numelc))
        inum_ptshel = 0
      ELSE
        ALLOCATE(inum_ptshel(0))
      ENDIF  
C
      IF (nperturb > 0) THEN
        ALLOCATE(xnum_rnoise(nperturb,numelc),stat=stat)
        IF (stat /= 0) CALL ancmsg(msgid=268,anmode=aninfo,
     .                           msgtype=msgerror,
     .                           c1='XNUM_RNOISE')
      ENDIF
C    
      CALL my_alloc(index2,numelc)
 
      IF(ndrape > 0 .AND. numelc_drape > 0) THEN
        DO i=1,numelc
            index2(i)=permutation%SHELL(i)
            eadd(i)=1
            itri(7,i)=i
            index(i)=i
            inum(1,i)=ipartc(i)
            inum(2,i)=isheoff(i)
            inum(3,i)=ixc(1,i)
            inum(4,i)=ixc(2,i)
            inum(5,i)=ixc(3,i)
            inum(6,i)=ixc(4,i)
            inum(7,i)=ixc(5,i)
            inum(8,i)=ixc(6,i)
            inum(9,i)=ixc(7,i)
            xnum(i)=thk(i)
            IF (nsubdom>0) inum_r2r(i) = tag_elcf(i)
            inum_workc(1,i) = iworksh(1,i)
            inum_workc(2,i) = iworksh(2,i)
            inum_workc(3,i) = iworksh(3,i)
            IF (nperturb > 0) THEN
              DO ipert = 1, nperturb
                xnum_rnoise(ipert,i) = rnoise(ipert,i) 
              ENDDO
            ENDIF
            angle(i)=sh4ang(i) 
            !drape structure
            ie = drapeg%INDX(i)
            xnum_drapeg%INDX(i) = drapeg%INDX(i)                                
            IF(ie == 0) cycle                                                   
            npt = drape(ie)%NPLY 
            xnum_drape(i)% NPLY  = npt                        
            xnum_drape(i)%INDX_PLY(1:npt) = drape(ie)%INDX_PLY(1:npt)
            npt = drape(ie)%NPLY_DRAPE                                      
            xnum_drape(i)%NPLY_DRAPE = npt   
            xnum_drape(i)%THICK  =  drape(ie)%THICK    
            DO jj = 1, npt                                                  
              drape_ply => drape(ie)%DRAPE_PLY(jj)                              
              nslice =  drape_ply%NSLICE                                        
              xnum_drape(i)%DRAPE_PLY(jj)%NSLICE =  nslice                      
              xnum_drape(i)%DRAPE_PLY(jj)%IPID   = drape_ply%IPID
              DO kk = 1,nslice                                                  
                xnum_drape(i)%DRAPE_PLY(jj)%IDRAPE(kk,1)=drape_ply%IDRAPE(kk,1) 
                xnum_drape(i)%DRAPE_PLY(jj)%IDRAPE(kk,2)=drape_ply%IDRAPE(kk,2) 
                xnum_drape(i)%DRAPE_PLY(jj)%RDRAPE(kk,1)=drape_ply%RDRAPE(kk,1) 
                xnum_drape(i)%DRAPE_PLY(jj)%RDRAPE(kk,2)=drape_ply%RDRAPE(kk,2)
              ENDDO    
              DEALLOCATE(drape_ply%IDRAPE, drape_ply%RDRAPE)
            ENDDO
             DEALLOCATE(drape(ie)%DRAPE_PLY)  
             DEALLOCATE(drape(ie)%INDX_PLY)                  
         ENDDO
       ELSE
        DO i=1,numelc
            index2(i)=permutation%SHELL(i)
            eadd(i)=1
            itri(7,i)=i
            index(i)=i
            inum(1,i)=ipartc(i)
            inum(2,i)=isheoff(i)
            inum(3,i)=ixc(1,i)
            inum(4,i)=ixc(2,i)
            inum(5,i)=ixc(3,i)
            inum(6,i)=ixc(4,i)
            inum(7,i)=ixc(5,i)
            inum(8,i)=ixc(6,i)
            inum(9,i)=ixc(7,i)
            xnum(i)=thk(i)
            IF (nsubdom>0) inum_r2r(i) = tag_elcf(i)
            inum_workc(1,i) = iworksh(1,i)
            inum_workc(2,i) = iworksh(2,i)
            inum_workc(3,i) = iworksh(3,i)          
            IF (nperturb > 0) THEN
              DO ipert = 1, nperturb
                xnum_rnoise(ipert,i) = rnoise(ipert,i) 
              ENDDO
            ENDIF
            angle(i)=sh4ang(i) 
         ENDDO
       ENDIF   
       IF(abs(isigi)==3.OR.abs(isigi)== 4.OR.abs(isigi)==5)THEN
         inum_ptshel(1:numelc) = ptshel(1:numelc)
       ENDIF
C
      IF(nadmesh /= 0)THEN
        DO  k=1,ksh4tree
          DO  i=1,numelc
            istor(k,i)=sh4tree(k,i)
          ENDDO
        ENDDO
        IF(lsh4trim/=0)THEN
          DO  i=1,numelc
            istor(ksh4tree+1,i)=sh4trim(i)
          ENDDO
        END IF
      END IF
C
      DO i=1,numelc
        xep(i)=cep(i)
      ENDDO
C      
      DO i = 1, numelc
        ii = i
C
        IF(nadmesh == 0)THEN
          itri(1,i)=0
        ELSE
C
C         ILEV must have strong weight on 1st key
          prt = ipartc(ii)
          iadm= ipart(10,prt)
          IF(iadm==0)THEN
C          not the same group as if adaptivity.
            itri(1,i)=0
          ELSE
            ilev= sh4tree(3,i)
            IF(ilev<0)ilev=-ilev-1
            itri(1,i)=ilev+1
          END IF
        END IF
C
        mid= ixc(1,ii)
        pid= ixc(6,ii)
        mln  = nint(pm(19,mid))
        igtyp= igeo(11,pid)
        jthe = nint(pm(71,mid))
        npn  = igeo(4,pid)
        ihbe = nint(geo(171,pid))
        ithk = nint(geo(35,pid))
        ipla = nint(geo(39,pid))
        irep = igeo(6,pid)
        ishxfem_ply = igeo(19,pid)
        nfail = 0
        ifail = 0
        ixfem = 0
        IF (igtyp == 11) THEN
          DO ipt = 1, npn
            imatly  = igeo(100+ipt,pid)
            nfail = max(nfail,mat_param(imatly)%NFAIL)
          ENDDO 
          IF(icrack3d > 0)THEN
C-          new multilayer -
            ixfem = mat_param(mid)%IXFEM
          ENDIF
        ELSEIF(igtyp == 17) THEN
          npn = iworksh(1,ii)
          isubstack =iworksh(3, ii)
!!          IIGEO   =  40 + 5*(ISUBSTACK - 1) 
!!          IADI    = IGEO(IIGEO + 3,PID)
!!          IPPID   = IADI    
          ippid = 2     
          DO ipt = 1, npn
!!            IPIDL = IGEO(IPPID+IPT,PID)
            ipidl = stack%IGEO(ippid  + ipt  ,isubstack)
            imatly  = igeo(101,ipidl)
            nfail = max(nfail,mat_param(imatly)%NFAIL)
          ENDDO   
        ELSEIF(igtyp == 51 ) THEN
C---
C new shell property (variable NPT through each layer)
C---
          nb_law58 = 0
          npn = iworksh(1,ii)
          isubstack = iworksh(3, ii)
          ippid = 2
          DO ipt = 1,npn  ! nb of plys
            ipidl  = stack%IGEO(ippid + ipt,isubstack)
            imatly = igeo(101,ipidl)
            nfail  = max(nfail,mat_param(imatly)%NFAIL)
C --- PID 51 combined with LAW58 ---
            IF (nint(pm(19,imatly)) == 58) nb_law58 = nb_law58 + 1
          ENDDO
C --- set IREP for sorting criteria:
          IF (nb_law58 == npn) THEN
            irep = 2
          ELSEIF (nb_law58 > 0) THEN
            irep = irep + 3
          ENDIF
        ELSEIF(igtyp == 52) THEN
C---
C new shell property (PCOMPP + STACK + PLY )
C---
          nb_law58 = 0
          npn = iworksh(1,ii)
          isubstack = iworksh(3, ii)
          ippid = 2
          ipmat = ippid + npn 
          DO ipt = 1,npn  ! nb of plys
            ipidl  = stack%IGEO(ippid + ipt,isubstack)
            imatly = stack%IGEO(ipmat + ipt,isubstack)
            nfail  = max(nfail,mat_param(imatly)%NFAIL)
C --- PID 51 combined with LAW58 ---
            IF (nint(pm(19,imatly)) == 58) nb_law58 = nb_law58 + 1
          ENDDO
C --- set IREP for sorting criteria:
          IF (nb_law58 == npn) THEN
            irep = 2
          ELSEIF (nb_law58 > 0) THEN
            irep = irep + 3
          ENDIF           
C
        ELSE ! IGTYP == 1
          nfail = mat_param(mid)%NFAIL
          IF(icrack3d > 0)THEN
C -         new monolayer -
            ixfem = mat_param(mid)%IXFEM
            IF (ixfem == 1) THEN
              ixfem = 2
              icrack3d = ixfem
            ENDIF
          END IF
        ENDIF
        IF (nfail > 0) ifail = 1
c        
C thermal material expansion
        iexpan  = ipm(218, mid)        
        icsen= igeo(3,pid)
        IF (icsen > 0) icsen=1
        IF(npn == 0.AND.(mln == 36.OR.mln == 86))THEN
          IF(ipla == 0) ipla=1
          IF(ipla == 2) ipla=0
C         IF(IPLA == 3) IPLA=2
        ELSEIF(npn == 0.AND.mln == 2)THEN
          IF(ipla == 2) ipla=0
        ELSE
         IF(ipla == 2) ipla=0
         IF(ipla == 3) ipla=2
        ENDIF
        IF(ithk == 2)THEN
          ithk = 0
        ELSEIF(mln == 32)THEN
          ithk = 1
        ENDIF
        ipla = iabs(ipla)
        ithk = iabs(ithk)
        istrain = nint(geo(11,pid))
        IF(mln == 19.OR.mln>=25.OR.mln == 15)istrain = 1
        issn = iabs(nint(geo(3,pid)))
C sorting on elem delete for rigidbody
C IRB = 0 : active elem
C IRB = 1 : inactive elem and optimized for SPMD
C IRB = 2 : inactive elem but optimized to be active in SPMD
        irb = isheoff(i)
C
C---     Key2 
        jsms = 0
        IF(isms/=0)THEN
          IF(idtgrs/=0)THEN
            IF(tagprt_sms(ipartc(ii))/=0)jsms=1
          ELSE
            jsms=1
          END IF
        END IF
C       JSMS=MY_SHIFTL(JSMS,0)
        itri(2,i) = jsms
C       NEXT=MY_SHIFTL(NEXT,1)
C
C---     Key3 
C       IPLA   = MY_SHIFTL(IPLA,0)
        istrain= my_shiftl(istrain,3)
        issn   = my_shiftl(issn,6)
        ihbe   = my_shiftl(ihbe,9)
        igtyp  = my_shiftl(igtyp,12)
        mln    = my_shiftl(mln,21)
        itri(3,i)=ipla+istrain+issn+ihbe+igtyp+mln
C
C---     Key4 
C
C       IRB    = MY_SHIFTL(IRB,0)
C
        ishxfem_ply  = my_shiftl(ishxfem_ply,10)
        ifail  = my_shiftl(ifail,11)
        iexpan = my_shiftl(iexpan,14)
        jthe   = my_shiftl(jthe,15)
        icsen= my_shiftl(icsen,16)
        npn  = my_shiftl(npn,17)
        irep = my_shiftl(irep,26)
        ithk = my_shiftl(ithk,30)
        IF(ixfem > 0) ixfem  = my_shiftl(ixfem,9)
C
        itri(4,i)=ithk+irep+npn+icsen+jthe+iexpan+irb+ifail+ishxfem_ply
     .           +ixfem
     
C---     Key5 
C       MID=MY_SHIFTL(MID,0)
        itri(5,i)=mid
C---     Key6 
C       PID=MY_SHIFTL(PID,0)
        itri(6,i)=pid
C --- key7 used for type17 iworkc=0 with/out type17 (or type51) PID  
        itri(7,i) =  iworksh(2,i) 
C --- key 8---------------------------------          
        itri(8,i )= damp_range_part(ipartc(ii))
      ENDDO
C
      mode=0
      CALL my_orders( mode, work, itri, index, numelc , 8)
C
       DO i=1,numelc
        ipartc(i) =inum(1,index(i))
        isheoff(i)=inum(2,index(i))
        IF (nsubdom>0) tag_elcf(i)=inum_r2r(index(i))
        thk(i)    =xnum(index(i))
      ENDDO
 
      DO i=1,numelc
        cep(i)=xep(index(i))
        permutation%SHELL(i)=index2(index(i))
      ENDDO
 
      DO  k=1,7
        DO  i=1,numelc
          ixc(k,i)=inum(k+2,index(i))
        ENDDO
      ENDDO
      IF(ndrape > 0 .AND. numelc_drape > 0 ) THEN
        ie = 0
        DO i=1,numelc
         iworksh(1,i)= inum_workc(1,index(i))
         iworksh(2,i)= inum_workc(2,index(i))
         iworksh(3,i)= inum_workc(3,index(i))
         IF (nperturb > 0) THEN
           DO ipert = 1, nperturb
             rnoise(ipert,i) = xnum_rnoise(ipert,index(i)) 
           ENDDO
         ENDIF
         sh4ang(i)=angle(index(i))
         !
         ie0 = xnum_drapeg%INDX(index(i))
         drapeg%INDX(i)= 0
         IF(ie0 == 0) cycle
         ie = ie + 1
         npt = xnum_drape(index(i))% NPLY    ! number of layer shell 
         drape(ie)%NPLY = npt
         drapeg%INDX(i)= ie
         ALLOCATE(drape(ie)%INDX_PLY(npt)) 
         drape(ie)%INDX_PLY(1:npt) =  xnum_drape(index(i))%INDX_PLY(1:npt)                                  
         npt = xnum_drape(index(i))%NPLY_DRAPE       ! NPT_DRP                                      
         ALLOCATE(drape(ie)%DRAPE_PLY(npt))  
         drape(ie)%NPLY_DRAPE= npt   
         drape(ie)%THICK = xnum_drape(index(i))%THICK     
         DO jj = 1, npt         
           drape_ply => drape(ie)%DRAPE_PLY(jj)                                  
           nslice = xnum_drape(index(i))%DRAPE_PLY(jj)%NSLICE                            
           drape_ply%NSLICE = nslice   
           drape_ply%IPID =  xnum_drape(index(i))%DRAPE_PLY(jj)%IPID
           ALLOCATE(drape_ply%IDRAPE(nslice,2), drape_ply%RDRAPE(nslice,2))
           DO kk = 1,nslice                                                               
            drape_ply%IDRAPE(kk,1) = xnum_drape(index(i))%DRAPE_PLY(jj)%IDRAPE(kk,1)      
            drape_ply%IDRAPE(kk,2) = xnum_drape(index(i))%DRAPE_PLY(jj)%IDRAPE(kk,2)      
            drape_ply%RDRAPE(kk,1) = xnum_drape(index(i))%DRAPE_PLY(jj)%RDRAPE(kk,1)      
            drape_ply%RDRAPE(kk,2) = xnum_drape(index(i))%DRAPE_PLY(jj)%RDRAPE(kk,2)      
           ENDDO                                                                          
         ENDDO 
        ENDDO
      ELSE
        DO i=1,numelc
         iworksh(1,i)= inum_workc(1,index(i))
         iworksh(2,i)= inum_workc(2,index(i))
         iworksh(3,i)= inum_workc(3,index(i))
         IF (nperturb > 0) THEN
           DO ipert = 1, nperturb
             rnoise(ipert,i) = xnum_rnoise(ipert,index(i)) 
           ENDDO
         ENDIF
          sh4ang(i)=angle(index(i))
        ENDDO
      ENDIF  
      IF(abs(isigi)==3.OR.abs(isigi)== 4.OR.abs(isigi)==5)THEN
         DO i=1,numelc
          ptshel(i) = inum_ptshel(index(i))
         ENDDO  
       ENDIF
C
      IF(nadmesh /= 0)THEN
        DO  k=1,ksh4tree
          DO  i=1,numelc
            sh4tree(k,i)=istor(k,index(i))
          ENDDO
        ENDDO
        IF(lsh4trim/=0)THEN
          DO  i=1,numelc
            sh4trim(i)=istor(ksh4tree+1,index(i))
          ENDDO
        END IF
      END IF
C
C INVERSION OF INDEX (IN ITR1)
C
      DO i=1,numelc
        itr1(index(i))=i
      ENDDO
C
C RENUMBERING OF THE TREE
      IF(nadmesh /= 0)THEN
        DO  i=1,numelc
          IF(sh4tree(1,i) /= 0)
     .       sh4tree(1,i)=itr1(sh4tree(1,i))
          IF(sh4tree(2,i) /= 0)
     .       sh4tree(2,i)=itr1(sh4tree(2,i))
        ENDDO
      END IF
C
C RENUMBERING FOR SURFACES
C
      DO i=1,nsurf
        nn=igrsurf(i)%NSEG
        DO j=1,nn
          IF(igrsurf(i)%ELTYP(j) == 3)
     .       igrsurf(i)%ELEM(j) = itr1(igrsurf(i)%ELEM(j))
        ENDDO
      ENDDO
C RENUMBERING FOR shell in Accel (gauge)
C
      DO i=1,nbgauge
         n1 = lgauge(1,i)
         IF(n1 <= 0) THEN
            n1=-lgauge(3,i)
            IF(n1 > 0) lgauge(3,i)=-itr1(n1)
         ENDIF
      ENDDO
C
C RENUMBERING FOR SHELL GROUPS
C
      DO i=1,ngrshel
        nn=igrsh4n(i)%NENTITY
        DO j=1,nn
          igrsh4n(i)%ENTITY(j) = itr1(igrsh4n(i)%ENTITY(j))
        ENDDO
      ENDDO
C
C renumbering INVERSE CONNECTIVITY
C
      DO i=1,4*numelc
        IF (nod2elc(i) /= 0) nod2elc(i)=itr1(nod2elc(i))
      END DO
C
C--------------------------------------------------------------
C         DETERMINATION OF SUPER_GROUPS
C--------------------------------------------------------------
      nd=1
      DO i=2,numelc
        ii0=itri(1,index(i))
        jj0=itri(1,index(i-1))
        ii =itri(2,index(i))
        jj =itri(2,index(i-1))
        ii1=itri(3,index(i))
        jj1=itri(3,index(i-1))
        ii2=itri(4,index(i))
        jj2=itri(4,index(i-1))
        ii3=itri(5,index(i))
        jj3=itri(5,index(i-1))
        ii4=itri(6,index(i))
        jj4=itri(6,index(i-1))
C for stack/ply pid
        ii5=itri(7,index(i))
        jj5=itri(7,index(i-1))
C damp freq range        
        ii6=itri(8,index(i))
        jj6=itri(8,index(i-1))       
        IF (ii0/=jj0 .or. 
     *   ii/=jj .or. 
     *   ii1/=jj1 .or.
     *   ii2/=jj2.OR.ii3 /= jj3.OR.ii4 /= jj4.OR.ii5 /= jj5 .or.
     *   ii6 /= jj6) THEN
          nd=nd+1
          eadd(nd)=i
        ENDIF
      ENDDO
      eadd(nd+1) = numelc+1
C-----------
c
      IF (nperturb > 0) THEN
        IF (ALLOCATED(xnum_rnoise)) DEALLOCATE(xnum_rnoise) 
      ENDIF
c
      DEALLOCATE(index2)
      DEALLOCATE( istor )
      IF(ndrape > 0 .AND. numelc_drape > 0) THEN
         DO i =1, numelc       
           ie = xnum_drapeg%INDX(i)
           IF(ie == 0) cycle                                
           npt_drp = xnum_drape(i)%NPLY_DRAPE  
           DO j = 1,npt_drp                                              
              DEALLOCATE(xnum_drape(i)%DRAPE_PLY(j)%RDRAPE)          
              DEALLOCATE(xnum_drape(i)%DRAPE_PLY(j)%IDRAPE)  
           ENDDO   
           DEALLOCATE(xnum_drape(i)%DRAPE_PLY,xnum_drape(i)%INDX_PLY) 
         ENDDO
         DEALLOCATE( xnum_drape ,xnum_drapeg%INDX)
      ELSE
         DEALLOCATE( xnum_drape )      
      ENDIF
      IF(ALLOCATED(inum_ptshel))DEALLOCATE(inum_ptshel)
 
      DEALLOCATE(angle,inum_workc)
      RETURN
      END