i21tri.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/.
!||====================================================================
!||    i21tri                 ../starter/source/interfaces/inter3d1/i21tri.F
!||--- called by ------------------------------------------------------
!||    i21buc1                ../starter/source/interfaces/inter3d1/i21buc3.F
!||--- calls      -----------------------------------------------------
!||    ancmsg                 ../starter/source/output/message/message.f
!||    i21cor3t               ../starter/source/interfaces/inter3d1/i21cor3t.F
!||    i7cmp3                 ../starter/source/interfaces/inter3d1/i7cmp3.F
!||    i7dst3                 ../starter/source/interfaces/inter3d1/i7dst3.F
!||    i7dstk                 ../starter/source/interfaces/inter3d1/i7dstk.F
!||    i7pen3                 ../starter/source/interfaces/inter3d1/i7pen3.F
!||--- uses       -----------------------------------------------------
!||    message_mod            ../starter/share/message_module/message_mod.F
!||====================================================================
      SUBROUTINE i21tri(
     1      BPE    ,PE    ,BPN   ,PN    ,ADD   ,
     2      IRECT  ,X     ,NB_NC ,NB_EC ,XYZM  ,
     3      I_ADD  ,NSV   ,I_AMAX,XMAX  ,YMAX  ,
     4      ZMAX   ,MAXSIZ,I_STOK,I_MEM ,NB_N_B,
     5      CAND_N ,CAND_E,NSN   ,NOINT ,TZINF ,
     6      MAXBOX ,MINBOX,J_STOK,MSR   ,XM0   ,
     7      MULTIMP,ITAB  ,GAP   ,GAP_S ,IGAP  ,
     8      GAPMIN,GAPMAX,MARGE  ,DEPTH ,DRAD  ,
     9      ID,TITR,
     1      IX1 ,IX2    ,IX3    ,IX4    ,NSVG   ,
     2      X1  ,X2     ,X3     ,X4     ,Y1     ,
     3      Y2  ,Y3     ,Y4     ,Z1     ,Z2     ,
     4      Z3  ,Z4     ,XI     ,YI     ,ZI     ,
     5      X0  ,Y0     ,Z0     ,STIF   ,NX1    ,
     6      NY1 ,NZ1    ,NX2    ,NY2    ,NZ2    ,
     7      NX3 ,NY3    ,NZ3    ,NX4    ,NY4    ,
     8      NZ4 ,P1     ,P2     ,P3     ,P4     ,
     9      LB1 ,LB2    ,LB3    ,LB4    ,LC1    ,
     1      LC2 ,LC3    ,LC4    ,PENE   ,PROV_N ,
     2      PROV_E,N11  ,N21    ,N31    ,DGAPLOAD)
      USE message_mod
      USE names_and_titles_mod , ONLY : nchartitle
C-----------------------------------------------
C   I m p l i c i t   T y p e s
C-----------------------------------------------
#include      "implicit_f.inc"
C-----------------------------------------------
C   G l o b a l   P a r a m e t e r s
C-----------------------------------------------
#include      "mvsiz_p.inc"
C-----------------------------------------------
C   C o m m o n   B l o c k s
C-----------------------------------------------
#include      "com01_c.inc"
#include      "com04_c.inc"
#include      "param_c.inc"
#include      "vect07_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER NB_NC,NB_EC,I_ADD,MAXSIZ,I_STOK,J_STOK,I_MEM
      INTEGER I_BID, I_AMAX,NB_N_B, NOINT, NSN,MULTIMP, IGAP
      INTEGER ADD(2,0:*),IRECT(4,*),BPE(*),PE(*),BPN(*),PN(*)
      INTEGER NSV(*),CAND_N(*),CAND_E(*), ITAB(*), MSR(*)
      INTEGER, DIMENSION(MVSIZ), INTENT(INOUT) ::PROV_N,PROV_E
      INTEGER, DIMENSION(MVSIZ), INTENT(INOUT) :: IX1,IX2,IX3,IX4,NSVG
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: X1,X2,X3,X4
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: Y1,Y2,Y3,Y4
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: Z1,Z2,Z3,Z4
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: XI,YI,ZI
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: X0,Y0,Z0,STIF
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: N11,N21,N31,PENE
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: NX1,NY1,NZ1
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: NX2,NY2,NZ2
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: NX3,NY3,NZ3
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: NX4,NY4,NZ4
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: P1,P2,P3,P4
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: LB1,LB2,LB3,LB4
      my_real, DIMENSION(MVSIZ), INTENT(INOUT) :: LC1,LC2,LC3,LC4
C     REAL
      my_real
     .   x(3,*),xyzm(6,*),tzinf,
     .   maxbox,minbox, xmax, ymax, zmax,
     .   gap, gap_s(*),
     .   gapmin, gapmax, marge, gapsmx, bgapsmx, depth, xm0(3,*), drad
      my_real , INTENT(IN) :: dgapload
      INTEGER ID
      CHARACTER(LEN=NCHARTITLE) :: TITR
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER NB_NCN, NB_ECN, ADDNN, ADDNE, I, J
      INTEGER INF,SUP,DIR,N1,N2,N3,N4,NN,NE
C     REAL
      my_real
     .   dx,dy,dz,dsup,seuil,xmx,xmn,
     .   gapv(mvsiz)
C-----------------------------------------------
C   ROLE OF THE ROUTINE:
C   ===================
C   CLASSIFIES BPE ELEMENTS AND BPN NODES INTO TWO ZONES
C   > OR < TO A HERE DEFINED BOUNDARY AND OUTPUTS EVERYTHING
C   INTO bpe,hpe, and bpn,hpn
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C
C     NOM          DESCRIPTION                       E/S
C
C     BPE          ARRAY OF FACETTES TO SORT      E/S
C                  AND OF THE RESULT ON THE MAX SIDE
C     PE           ARRAY OF FACETTES              S
C                  RESULTAT COTE MIN
C     BPN          SORTED NODES ARRAY        E/S
C                  AND OF THE RESULT ON THE MAX SIDE
C     PN           NODES ARRAY S
C                  RESULTAT COTE MIN
C     ADD(2,*)     ARRAY OF ADRESSES              E/S 
C          1.......ADRESSES NODES C          2.......ADRESSES ELEMENTS
C     ZYZM(6,*)     ARRAY OF XYZMIN               E/S 
C          1.......XMIN BOITE
C          2.......YMIN BOITE
C          3.......ZMIN BOITE
C          4.......XMAX BOITE
C          5.......YMAX BOITE
C          6.......ZMAX BOITE
C     IRECT(4,*)   ARRAY OF CONEC FACETTES        E
C     X(3,*)       COORDONNEES NODALES               E
C     NB_NC        NUMBER OF CANDIDATE NODES        E/S
C     NB_EC        NUMBER OF CANDIDATE ELEMENTS           E/S
C     I_ADD        POSITION IN THE E/S ADDRESS TABLE
C     NSV          NOS SYSTEMES DES NODES E
C     Xmax larger abcisse existing e
C     XMAX largest order.existing E
C     Xmax larger existing side E
C     MAXSIZ       TAILLE MEMOIRE MAX POSSIBLE       E
C     I_STOK       storage level of pairs
C                                CANDIDATES impact    E/S
C     CAND_N       boites resultats nodes C     CAND_E       adresses des boites resultat elements
C                  MULTIMP*NSN MAX SIZE NOW ALLOWED FOR THE
C                  COUPLES NODES,ELT CANDIDATES
C     NOINT        INTERFACE USER NUMBER
C     TZINF        TAILLE ZONE INFLUENCE
C     MAXBOX       TAILLE MAX BUCKET
C     MINBOX       TAILLE MIN BUCKET
C=======================================================================
C
C
C    1- TEST ARRET = BOITE VIDE
C                    BOITE TROP PETITE 
C                    BOITE NE CONTENANT QU'ONE NODE C                    NO MORE MEMORY AVAILABLE
C
C-----------------------------------------------------------
C
C      IF(MEMX>ADD(2,1)+NB_EC)THEN
C        WRITE(ISTDO,*)' *******MEM MAX=',MEMX
C        MEMX=-1
C      ELSEIF(MEMX/=-1)THEN
C        MEMX=ADD(2,1)+NB_EC
C      ENDIF
C--------------------TEST ON EMPTY BOXES--------------
C
      IF(nb_ec==0.OR.nb_nc==0) THEN
C       write(6,*)" BOITE VIDE"
C       MUST COPY THE BOTTOMS OF THE STACKS INTO CORRESPONDING BAS_DE_PILE
C       BEFORE DESCENDING AGAIN INTO THE MIDDLE BRANCH
        CALL i7dstk(i_add,nb_nc,nb_ec,add,bpn,pn,bpe,pe)
        RETURN
      ENDIF
C
C-------------------TEST ON END OF BRANCH / MEMORY OVERFLOW------------
C
      dx = xyzm(4,i_add) - xyzm(1,i_add)
      dy = xyzm(5,i_add) - xyzm(2,i_add)
      dz = xyzm(6,i_add) - xyzm(3,i_add)
      dsup= max(dx,dy,dz)
C
      IF(add(2,1)+nb_ec>=maxsiz) THEN
C       NO MORE SPACE IN THE BOXED ELEMENTS STACK, BOXES TOO SMALL
        IF ( nb_n_b == numnod) THEN
C          WRITE(IOUT,*)'***ERROR INFINITE LOOP DETECTED '
C          WRITE(ISTDO,*)'***ERROR INFINITE LOOP DETECTED '
C          CALL ARRET(2)
           IF (istamping == 1) THEN
             CALL ancmsg(msgid=775,
     .                   msgtype=msgerror,
     .                   anmode=aninfo,
     .                   i1=id,
     .                   c1=titr)
           ELSE
             CALL ancmsg(msgid=685,
     .                   msgtype=msgerror,
     .                   anmode=aninfo,
     .                   i1=id,
     .                   c1=titr)
           ENDIF
        ENDIF
        i_mem = 1
        RETURN
      ENDIF
      IF(dsup<minbox.OR.
     .   nb_nc<=nb_n_b.AND.dsup<maxbox.OR.   
     .   nb_nc<=nb_n_b.AND.nb_ec==1) THEN
C
C       write(6,*)" NOUVELLE BOITE "
C       1- STORAGE OF THE CANDIDATE NODE(S) AND CORRESPONDING ELEMENTS
C       REMOVE THE USELESS ONES
        DO 20 i=1,nb_ec
            ne = bpe(i)
            n1=msr(irect(1,ne))
            n2=msr(irect(2,ne)) 
            n3=msr(irect(3,ne)) 
            n4=msr(irect(4,ne)) 
            DO 20 j=1,nb_nc
              nn=nsv(bpn(j))
              IF(nn/=n1.AND.nn/=n2.AND.nn/=n3.AND.nn/=n4) THEN
                j_stok = j_stok + 1
                prov_n(j_stok) = bpn(j)
                prov_e(j_stok) = ne
                IF(j_stok==nvsiz) THEN
                  lft = 1
                  llt = nvsiz
                  nft = 0
                  j_stok = 0
                  CALL i21cor3t(x    ,irect,nsv  ,prov_e ,prov_n,
     2                        gapv ,igap  ,gap   ,gap_s,gapmin ,
     3                        gapmax,xm0  ,depth ,drad ,ix1    ,
     4                        ix2   ,ix3   ,ix4   ,nsvg    ,x1     ,
     5                        x2    ,x3    ,x4    ,y1      ,y2     ,
     6                        y3    ,y4    ,z1    ,z2      ,z3     ,
     7                        z4    ,xi    ,yi    ,zi      ,dgapload)
                  CALL i7dst3(ix3,ix4,x1 ,x2 ,x3 ,
     1                  x4 ,y1 ,y2 ,y3 ,y4 ,
     2                  z1 ,z2 ,z3 ,z4 ,xi ,
     3                  yi ,zi ,x0 ,y0 ,z0 ,
     4                  nx1,ny1,nz1,nx2,ny2,
     5                  nz2,nx3,ny3,nz3,nx4,
     6                  ny4,nz4,p1 ,p2 ,p3 ,
     7                  p4 ,lb1,lb2,lb3,lb4,
     8                  lc1,lc2,lc3,lc4,llt)
                  CALL i7pen3(marge,gapv,n11,n21,n31,
     1                  pene ,nx1 ,ny1,nz1,nx2,
     2                  ny2  ,nz2 ,nx3,ny3,nz3,
     3                  nx4  ,ny4 ,nz4,p1 ,p2 ,
     4                  p3   ,p4,llt)
                  IF(i_stok+nvsiz<multimp*nsn) THEN
                    CALL i7cmp3(i_stok,cand_e  ,cand_n,1,pene,
     1                  prov_n,prov_e)
                  ELSE
                    i_bid = 0
                    CALL i7cmp3(i_bid,cand_e,cand_n,0,pene,
     1                  prov_n,prov_e)
                    IF(i_stok+i_bid<multimp*nsn) THEN
                      CALL i7cmp3(i_stok,cand_e,cand_n,1,pene,
     1                  prov_n,prov_e)
                    ELSE
                      i_mem = 2 ! TOO MANY POSSIBLE IMPACTS : return to upgrade Multimp
                      RETURN
                    ENDIF
                  ENDIF
                ENDIF
C               write(6,*)"Node candidat",BPN(J)
C               write(6,*)"Element  candidat",NE
              ENDIF  
   20   CONTINUE 
C       MUST COPY THE BOTTOMS OF THE STACKS INTO CORRESPONDING BAS_DE_PILE
C       BEFORE DESCENDING AGAIN INTO THE MIDDLE BRANCH
        CALL i7dstk(i_add,nb_nc,nb_ec,add,bpn,pn,bpe,pe)
        RETURN
      ENDIF
C
C-----------------------------------------------------------
C
C
C    2- SORTING PHASE ON THE MEDIANE ACCORDING TO THE LARGEST DIRECTION
C                    
C                   
C-----------------------------------------------------------
C
C
C    1- DETERMINER LA DIRECTION A DIVISER X,Y OU Z
C
      dir = 1
      IF(dy==dsup) THEN
        dir = 2
      ELSE IF(dz==dsup) THEN
        dir = 3
      ENDIF
      seuil =(xyzm(dir+3,i_add)+xyzm(dir,i_add))/2
C
C    2- DIVISER LES NODES EN TWO ZONES 
C
      nb_ncn= 0
      addnn= add(1,1)
      inf = 0
      sup = 0
      IF(igap==0)THEN
       DO i=1,nb_nc
        IF(x(dir,nsv(bpn(i)))<seuil) THEN
C         STORED IN THE BOTTOM OF THE BP STACK
          addnn = addnn + 1
          pn(addnn) = bpn(i)
          inf = 1
        ELSE
          nb_ncn = nb_ncn + 1
          bpn(nb_ncn) = bpn(i)
C         ON STOCKE EN ECRASANT PROGRESSIVEMENT BPN
          sup = 1
        ENDIF
       END DO
      ELSE
       gapsmx  = zero
       bgapsmx = zero
       DO i=1,nb_nc
        IF(x(dir,nsv(bpn(i)))<seuil) THEN
C         STORED IN THE BOTTOM OF THE BP STACK
          addnn = addnn + 1
          pn(addnn) = bpn(i)
          gapsmx = max(gapsmx,max(gap_s(bpn(i)),depth,drad))
          inf = 1
        ELSE
C         ON STOCKE EN ECRASANT PROGRESSIVEMENT BPN
          nb_ncn = nb_ncn + 1
          bpn(nb_ncn) = bpn(i)
          bgapsmx = max(bgapsmx,max(gap_s(bpn(i)),depth,drad))
          sup = 1
        ENDIF
       END DO
      END IF
CC
CC    3- DIVISER LES ELEMENTS 
CC
C      NB_ECN= 0
C      ADDNE= ADD(2,1)
C      SEUILI = SEUIL-TZINF
C      SEUILS = SEUIL+TZINF
C      DO 85 I=1,NB_EC
C        INF = 0
C        SUP = 0
C        DO 80 J=1,4
C          IP = IRECT(J,BPE(I))
C          IF(X(DIR,IP)<SEUILS) THEN
C            INF = 1
C            IF(SUP==1) GOTO 81
C          ENDIF
C          IF(X(DIR,IP)>=SEUILI) THEN
C            SUP = 1
C            IF(INF==1) GOTO 81
C          ENDIF
C   80   CONTINUE
C
C   81   CONTINUE
C        IF(INF==1) THEN
C         STORED IN THE BOTTOM OF THE BP STACK
C          ADDNE = ADDNE + 1
C          PE(ADDNE) = BPE(I)
C        ENDIF
C        IF(SUP==1) THEN
C         ON STOCKE EN ECRASANT PROGRESSIVEMENT BPE
C          NB_ECN = NB_ECN + 1
C          BPE(NB_ECN) = BPE(I)
C        ENDIF
C   85 CONTINUE
C
C    3- DIVIDE THE ELEMENTS
C
C 2 LIGNES PROV FOR TEST
C      INF = 1
C      SUP = 1
C
      nb_ecn= 0
      addne= add(2,1)
      IF(igap==0)THEN
       DO i=1,nb_ec
        xmx = max(xm0(dir,irect(1,bpe(i))),xm0(dir,irect(2,bpe(i))),
     .            xm0(dir,irect(3,bpe(i))),xm0(dir,irect(4,bpe(i))))
     .      + tzinf
        xmn = min(xm0(dir,irect(1,bpe(i))),xm0(dir,irect(2,bpe(i))),
     .            xm0(dir,irect(3,bpe(i))),xm0(dir,irect(4,bpe(i))))
     .      - tzinf
        IF(xmn<seuil.AND.inf==1) THEN
C         STORED IN THE BOTTOM OF THE BP STACK
          addne = addne + 1
          pe(addne) = bpe(i)
        ENDIF
        IF(xmx>=seuil.AND.sup==1) THEN
C         ON STOCKE EN ECRASANT PROGRESSIVEMENT BPE
          nb_ecn = nb_ecn + 1
          bpe(nb_ecn) = bpe(i)
        ENDIF
       ENDDO
      ELSE       
       DO i=1,nb_ec
        xmn = min(xm0(dir,irect(1,bpe(i))),xm0(dir,irect(2,bpe(i))),
     .            xm0(dir,irect(3,bpe(i))),xm0(dir,irect(4,bpe(i))))
     -       -max(min(max(gapsmx,gapmin),gapmax)+dgapload,depth,drad)
     -       -marge
        IF(xmn<seuil.AND.inf==1) THEN
C         STORED IN THE BOTTOM OF THE BP STACK
          addne = addne + 1
          pe(addne) = bpe(i)
        ENDIF
        xmx = max(xm0(dir,irect(1,bpe(i))),xm0(dir,irect(2,bpe(i))),
     .            xm0(dir,irect(3,bpe(i))),xm0(dir,irect(4,bpe(i))))
     +       +max(min(max(bgapsmx,gapmin),gapmax)+dgapload,depth,drad)
     +       +marge
        IF(xmx>=seuil.AND.sup==1) THEN
C         ON STOCKE EN ECRASANT PROGRESSIVEMENT BPE
          nb_ecn = nb_ecn + 1
          bpe(nb_ecn) = bpe(i)
        ENDIF
       ENDDO
      END IF
C
C    4- REMPLIR LES TABLEAUX D'ADRESSES
C
      add(1,2) = addnn
      add(2,2) = addne
C-----fill the mins of the next box and the max of the current one
C     (i.e. threshold is a max for the current)
C     We're going to go down and so we define a new box
C     fill the max of the new box
C     initialises in i7buc1 a 1.E30 comme ca on recupere
C     either XMAX or the max of the box
      xyzm(1,i_add+1) = xyzm(1,i_add) 
      xyzm(2,i_add+1) = xyzm(2,i_add)
      xyzm(3,i_add+1) = xyzm(3,i_add)
      xyzm(4,i_add+1) = xyzm(4,i_add)
      xyzm(5,i_add+1) = xyzm(5,i_add)
      xyzm(6,i_add+1) = xyzm(6,i_add)
      xyzm(dir,i_add+1) = seuil
      xyzm(dir+3,i_add) = seuil
C
      nb_nc = nb_ncn
      nb_ec = nb_ecn
C     increment the descent level before exiting
      i_add = i_add + 1
      IF(i_add>=1000) THEN
C       TROP NIVEAUX PILE ON VAS LES PRENDRE PLUS GRANDES...
        IF ( nb_n_b == numnod) THEN
C          WRITE(IOUT,*)'***COMPUTATION STOPPED WHILE INFINITELY LOOPING'
C          WRITE(ISTDO,*)'***COMPUTATION STOPPED WHILE INFINITELY LOOPING'
C          CALL ARRET(2)
           IF (istamping == 1) THEN
             CALL ancmsg(msgid=775,
     .                   msgtype=msgerror,
     .                   anmode=aninfo,
     .                   i1=id,
     .                   c1=titr)
           ELSE
             CALL ancmsg(msgid=685,
     .                   msgtype=msgerror,
     .                   anmode=aninfo,
     .                   i1=id,
     .                   c1=titr)
           ENDIF
        ENDIF
        i_mem = 1
        RETURN
      ENDIF
C
C     this return is only reached if ok = 0
      RETURN
      END