i9wal2.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/.
!||====================================================================
!||    i9wal2            ../engine/source/interfaces/int09/i9wal2.F
!||--- called by ------------------------------------------------------
!||    i9wale            ../engine/source/interfaces/int09/i9wale.F
!||--- calls      -----------------------------------------------------
!||    bcs2              ../engine/source/constraints/general/bcs/bcs2.f
!||    i9grd2            ../engine/source/interfaces/int09/i9grd2.F
!||    spmd_glob_dsum9   ../engine/source/mpi/interfaces/spmd_th.F
!||    spmd_glob_isum9   ../engine/source/mpi/interfaces/spmd_th.F
!||    spmd_ibcast       ../engine/source/mpi/generic/spmd_ibcast.F
!||    spmd_rbcast       ../engine/source/mpi/generic/spmd_rbcast.F
!||--- uses       -----------------------------------------------------
!||    elbufdef_mod      ../common_source/modules/mat_elem/elbufdef_mod.F90
!||    element_mod       ../common_source/modules/elements/element_mod.F90
!||====================================================================
      SUBROUTINE i9wal2(X,V   ,W         ,A         ,CRST      ,
     2   NSV       ,ILOC      ,IRTL      ,ICODE     ,ISKEW     ,
     3   SKEW      ,MSR       ,LMSR      ,NSEG      ,IRECTS    ,
     4   IRECT     ,UPW       ,IXQ       ,ELBUF_TAB ,
     5   IPARG     ,PM        ,NALE      ,EE        ,IELES     ,
     6   IELEM     ,TSTIF     ,INTTH     ,IEULT     ,STENS     ,
     7   ISIZES    ,ISIZEM, NRTS, NRTM,NSN,NMN) 
C-----------------------------------------------
C   M o d u l e s
C-----------------------------------------------
      USE elbufdef_mod      
      use element_mod , only : nixq      
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      "com01_c.inc"
#include      "com08_c.inc"
#include      "scr08_a_c.inc"
#include      "param_c.inc"
#include      "task_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER, INTENT(IN) :: NRTS, NRTM,NSN,NMN
      INTEGER IRECT(4,*), NSV(NSN), ILOC(*), IRTL(NSN), ICODE(*), ISKEW(*),
     .   MSR(*), IRECTS(4,*), LMSR(*), NSEG(*),IXQ(NIXQ,*),
     .   IPARG(NPARG,*), IELES(*), IELEM(*),NALE(*) ,
     .   INTTH, IEULT, ISIZES, ISIZEM
C     REAL
      my_real
     .   upw, tstif,ttt, stens,
     .   x(3,*), v(3,*), w(3,*), a(3,*), crst(2,*), skew(lskew,*),
     .   pm(npropm,*),ee(*)
      TYPE (ELBUF_STRUCT_), TARGET, DIMENSION(NGROUP) :: ELBUF_TAB
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER II, N, L, JJ, NN, LL, LL1, LL2, LG, KK, KKK, K1, 
     . I1, IERR, IGROU, IELN,
     .   ixx(4), iperm(2), jperm(2),
     .   itemp(2), is, im, ilen,
     .   tags(isizes),tagm(isizem), lists(isizes),listm(isizem),
     .   icomerr(isizem+isizes),icomngr(isizem+isizes),
     .   icomnel(isizem+isizes)
C     REAL
      my_real h(2),
     . vmy, vmz, vy, vz, vv, ny, nz, vt,
     . nny, nnz, fac, p, ty, tz,
     .   efric, vols, volm, ts, tm ,tstifm, tstifs, dvn,
     .   tstift, phi, areas, aream, vn, wn, stensy, stensz,
     .   comarea(isizem+isizes),comstf(isizem+isizes),
     .   comt(isizem+isizes),comvol(isizem+isizes),combuf(isizem+isizes)
 
c
      DATA iperm/ 2, 1/,jperm/ 1, -1/
C-----------------------------------------------
C
C preparation phase for spmd
C
      IF(intth/=zero) THEN
       IF(ispmd==0) THEN
        DO ii = 1, nrts
          tags(ii) = 0
        ENDDO
        DO ii = 1, nrtm
          tagm(ii) = 0
        ENDDO
        is = 0
        im = 0
        DO ii = 1, nsn
          l = irtl(ii)
          IF(iloc(ii)>0.AND.nmn>0)THEN
           IF(tagm(l)==0)THEN
             im = im + 1
             listm(im) = l
             tagm(l) = im
           END IF
           ll1=nseg(ii)
           ll2=nseg(ii+1)-1
           DO ll=ll1,ll2
             lg  = lmsr(ll)
             IF(tags(lg)==0) THEN
               is = is + 1
               lists(is) = lg
               tags(lg) = is
             ENDIF
           ENDDO
          ENDIF
        ENDDO
C
C Elem lists compacting
C
        itemp(1) = is
        itemp(2) = im
       ENDIF
C
C Sending second/mains facet list in contact
C
       IF(nspmd > 1) THEN
         CALL spmd_ibcast(itemp,itemp,1,2,0,2)
         is = itemp(1)
         im = itemp(2)
         ilen = im+is
         CALL spmd_ibcast(listm,listm,1,im,0,2)
         CALL spmd_ibcast(lists,lists,1,is,0,2)
       END IF
       DO ii = 1, im
         l = listm(ii)
         ix(1) = msr(irect(1,l))
         ix(2) = msr(irect(2,l))
         IF(ielem(l)>0) THEN
          CALL i9grd2(
     1       ierr ,comarea(ii),comstf(ii),comt(ii),comvol(ii),
     2       ielem(l)   ,x  ,ixq(1,ielem(l)),ix        ,
     3       iparg,pm         ,elbuf_tab , igrou   ,ieln    )
           icomerr(ii) = ierr
           icomngr(ii) = igrou
           icomnel(ii) = ieln
         ELSE
           comarea(ii) = zero
           comstf(ii)  = zero
           comt(ii)    = zero
           comvol(ii)  = zero
           icomerr(ii) = 0
           icomngr(ii) = 0
           icomnel(ii) = 0
         ENDIF
         combuf(ii) = zero
       ENDDO
C
       DO ii = 1, is
         l = lists(ii)
         ixx(1)=nsv(irects(1,l))
         ixx(2)=nsv(irects(2,l))
         IF(ieles(l)>0) THEN
          CALL i9grd2(
     1     ierr ,comarea(im+ii),comstf(im+ii),comt(im+ii),comvol(im+ii),
     2     ieles(l)      ,x    ,ixq(1,ieles(l))    ,ixx                ,
     3     iparg,pm            ,elbuf_tab ,igrou   ,ieln               )
           icomerr(im+ii) = ierr
           icomngr(im+ii) = igrou
           icomnel(im+ii) = ieln
         ELSE
           comarea(im+ii) = zero
           comstf(im+ii)  = zero
           comt(im+ii)    = zero
           comvol(im+ii)  = zero
           icomerr(im+ii) = 0
           icomngr(im+ii) = 0
           icomnel(im+ii) = 0
         ENDIF
         combuf(im+ii) = zero
       ENDDO
C
       IF (nspmd > 1) THEN
C
C gather values
C
                 CALL spmd_glob_dsum9(comarea,ilen)
                 CALL spmd_glob_dsum9(comstf,ilen)
                 CALL spmd_glob_dsum9(comt,ilen)
                 CALL spmd_glob_dsum9(comvol,ilen)
                 CALL spmd_glob_isum9(icomerr,ilen)
                 CALL spmd_glob_isum9(icomngr,ilen)
                 CALL spmd_glob_isum9(icomnel,ilen)
C   node part on P0 only
              IF(ispmd/=0) GOTO 900
       END IF
C   interface processed by p0
      ELSE
       IF(ispmd/=0) RETURN
      ENDIF
C
      DO 800 ii=1,nsn
      ll1=nseg(ii)
      ll2=nseg(ii+1)-1
      n=nsv(ii)
      IF(iloc(ii)>0.AND.nmn>0)THEN
C---------------------------------
C       CONTACT
C---------------------------------
        l=irtl(ii)
        DO 10 jj=1,2
        nn=irect(jj,l)
 10     ix(jj)=msr(nn)
C
        h(1) = half*(one - crst(1,ii))
        h(2) = half*(one + crst(1,ii))
C---------------------------------
C       mesh velocity
C---------------------------------
        vmy=zero
        vmz=zero
C
        DO jj=1,2
         vmy=vmy+w(2,ix(jj))*h(jj)
         vmz=vmz+w(3,ix(jj))*h(jj)
        ENDDO
C
        w(2,n)=vmy
        w(3,n)=vmz
C---------------------------------
C       PONT THERMIQUE
C---------------------------------
        IF(intth/=zero)THEN
           kk = tagm(l)
           efric = half * ee(ii) / (ll2-ll1+1)
             ierr = icomerr(kk)
             aream = comarea(kk)
             tstifm = comstf(kk)
             tm = comt(kk)
             volm = comvol(kk)
           IF(ierr==0) THEN
             DO ll = ll1,ll2
               lg = lmsr(ll)
               jj = tags(lg) + im
               ierr = icomerr(jj)
               areas = comarea(jj)
               tstifs = comstf(jj)
               ts = comt(jj)
               vols = comvol(jj)
               IF(ierr==0) THEN
                 tstift = tstifm + tstifs + tstif
                 phi = areas * dt1 * (tm-ts) / tstift
                 combuf(jj) = combuf(jj)
     +                      + (efric+phi)/vols
                 combuf(kk) = combuf(kk)
     +                      + (efric-phi)/volm
               ENDIF
             ENDDO
           ENDIF
c         ENDIF
        ENDIF
C
      ELSEIF(iloc(ii)<0.OR.nmn==0)THEN
C---------------------------------------
C       PAS DE CONTACT => SURFACE LIBRE
C---------------------------------------
        iloc(ii) = -iloc(ii)
C
        vy = v(2,n) - w(2,n)
        vz = v(3,n) - w(3,n)
        vv = max(em30,sqrt(vy**2+vz**2))
        nny = zero
        nnz = zero
C------------------------------------------------------
C       loop over facets connected to node ii
C------------------------------------------------------
        DO 300 ll=ll1,ll2
         lg=lmsr(ll)
         DO 200 kkk=1,2
          kk=kkk
 200     IF(irects(kk,lg)==ii) GO TO 250
 250     CONTINUE
C------------------------------------------------------
C         calculation of the normal with upwind on the upstream side
C------------------------------------------------------
          k1 = iperm(kk)
          i1 = nsv(irects(k1,lg))
          ty = x(2,i1) - x(2,n)
          tz = x(3,i1) - x(3,n)
          ttt = max(em30,sqrt(ty**2+tz**2))
C          VT = VY*TY + VZ*TZ
          vt = v(2,n)*ty + v(3,n)*tz
          p = onep0001 - upw*(half + sign(half,vt))
          ny = tz
          nz =-ty
C          FAC = P / MAX(EM30,SQRT(NY**2+NZ**2))
          fac = p * jperm(kk)
          nny = nny + ny*fac
          nnz = nnz + nz*fac
C-------------------------------------
C       surface tension
C-------------------------------------
          stensy = stens * ty / ttt
          stensz = stens * tz / ttt
          a(2,n)  = a(2,n)  + stensy
          a(3,n)  = a(3,n)  + stensz
 300    CONTINUE
        fac = max(em30,sqrt(nny**2+nnz**2))
        nny = nny/fac
        nnz = nnz/fac
C---------------------------------
C       grid boundary conditions
C---------------------------------
        IF(icode(n)/=0)THEN
C---------------------------------
C       W LAGRANGIEN ACCORDING TO N
C---------------------------------
          dvn = vy * nny + vz * nnz
          w(2,n) = w(2,n) + dvn * nny
          w(3,n) = w(3,n) + dvn * nnz
          CALL bcs2(w(1,n),skew(1,iskew(n)),iskew(n),icode(n))
          vn = v(2,n)*nny + v(3,n)*nnz
          wn = w(2,n)*nny + w(3,n)*nnz
C-------------------------------------
C         W LAGRANGIEN ACCORDING TO N + BCS
C-------------------------------------
          IF(abs(wn)>em30)THEN
            fac = vn / wn
            w(2,n) = w(2,n) * fac
            w(3,n) = w(3,n) * fac
          ENDIF
        ELSEIF(ieult/=0)THEN
C-------------------------------------
C         W LAGRANGIEN ACCORDING TO N
C         W EULERIEN ACCORDING TO T
C---------------------------------
          vn = v(2,n) * nny + v(3,n) * nnz
          w(2,n) = vn * nny
          w(3,n) = vn * nnz
        ELSE
C-------------------------------------
C         W LAGRANGIEN ACCORDING TO N
C         LIBRE ACCORDING TO T
C---------------------------------
          dvn = vy * nny + vz * nnz
          w(2,n) = w(2,n) + dvn * nny
          w(3,n) = w(3,n) + dvn * nnz
        ENDIF
      ENDIF
C
  800 CONTINUE
C
C finalization phase for spmd
C
  900 CONTINUE
      IF(intth/=zero) THEN
       IF(nspmd > 1) THEN
C
C Envoi buffer elements updates
C
         CALL spmd_rbcast(combuf,combuf,1,ilen,0,2)
       END IF
C
C Local Elbuf update
C
       DO ii = 1, im
         l = listm(ii)
         IF(ielem(l)>0) THEN
           igrou = icomngr(ii)
           ieln  = icomnel(ii)
           elbuf_tab(igrou)%GBUF%EINT(ieln) = 
     .     elbuf_tab(igrou)%GBUF%EINT(ieln) + combuf(ii)
         ENDIF
       ENDDO
C
       DO ii = 1, is
         l =  lists(ii)
         IF(ieles(l)>0) THEN
           igrou = icomngr(im+ii)
           ieln  = icomnel(im+ii)
           elbuf_tab(igrou)%GBUF%EINT(ieln) = 
     .     elbuf_tab(igrou)%GBUF%EINT(ieln) + combuf(im+ii)
         ENDIF
       ENDDO
      ENDIF
C
      RETURN
      END