inivol_set.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/.
!||====================================================================
!||    inivol_set           ../starter/source/initial_conditions/inivol/inivol_set.F
!||--- called by ------------------------------------------------------
!||    init_inivol          ../starter/source/initial_conditions/inivol/init_inivol.F90
!||--- uses       -----------------------------------------------------
!||    message_mod          ../starter/share/message_module/message_mod.F
!||====================================================================
      SUBROUTINE inivol_set(UVAR     , NUVAR , NEL     , KVOL         , MLW                         ,
     .                      ELBUF_TAB, NG    , NBSUBMAT, MULTI_FVM    ,
     .                      IDP      , IPART , NFT     , IMID,
     .                      MAT_PARAM)
C-----------------------------------------------
C   M o d u l e s
C-----------------------------------------------
      USE elbufdef_mod
      USE multi_fvm_mod
      USE elbufdef_mod
      USE message_mod
      USE constant_mod, ONLY : nine_over_10
      USE multimat_param_mod , ONLY : m51_n0phas, m51_nvphas
      USE matparam_def_mod  , ONLY : matparam_struct_
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      "com04_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER,INTENT(IN)                 :: IMID
      INTEGER,INTENT(IN)                 :: NEL, NUVAR, MLW, NG, NBSUBMAT,IDP,IPART(*),NFT
      my_real,INTENT(INOUT)              :: KVOL(NBSUBMAT,NEL)
      my_real,INTENT(INOUT)              :: uvar(nel,nuvar)
      TYPE(multi_fvm_struct),INTENT(IN)  :: MULTI_FVM
      TYPE(elbuf_struct_), TARGET, DIMENSION(NGROUP),INTENT(INOUT) :: ELBUF_TAB
      TYPE (MATPARAM_STRUCT_) ,DIMENSION(NUMMAT) ,INTENT(IN) :: MAT_PARAM
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER                 :: I, IMAT, KK
      TYPE(G_BUFEL_) ,POINTER :: GBUF
      TYPE(L_BUFEL_) ,POINTER :: LBUF
      my_real                 :: p, p1,p2,p3,p4,sumvf,ratio
      my_real                 :: vfrac0(nbsubmat)
      INTEGER                 :: default_SUBMAT_id(NEL)
C-----------------------------------------------
C   S o u r c e  L i n e s
C-----------------------------------------------
 
      ! DETERMINE DEFAULT SUBMAT ID
      vfrac0(1:nbsubmat) = mat_param(imid)%MULTIMAT%VFRAC(1:nbsubmat)
      default_submat_id(1:nel) = max(1, maxloc(vfrac0,1))
 
      !---CHECK UNOCCUPIPED SUBVOLUME AND FILL WITH DEFAULT PHASE IF SUM(vfrac>1) display an error message.
      DO i=1,nel
       IF(ipart(i+nft) /= idp)cycle
        sumvf = sum(kvol(1:nbsubmat,i))
        imat = default_submat_id(i)
        IF(sumvf > one+em06)THEN
          !not expected, user input must be checked : error message
          sumvf = sum(kvol(1:nbsubmat,i))
          sumvf = sumvf - kvol(imat,i)
          IF(sumvf <= one .AND. sumvf > zero)THEN  ! sumvf is here calculated in range 2:nbsumat
            kvol(imat,i)=one-sumvf
          ELSE
            sumvf = sum(kvol(1:nbsubmat,i))
            ratio=one/sumvf
            kvol(1:nbsubmat,i)=ratio*kvol(1:nbsubmat,i)
          ENDIF
        ELSEIF(sumvf < one-em06)THEN
          !fill unoccupied subvolume with phase-1
          kvol(imat,i) = kvol(imat,i) + one-sumvf
        ELSEIF(sumvf >= one-em06 .AND. sumvf <= one+em06)THEN
          !get rid of precision issue so that sumvf is exactly 1.000000
          ratio=one/sumvf
          kvol(1:nbsubmat,i)=ratio*kvol(1:nbsubmat,i)
        ENDIF
      ENDDO
 
      !---set initial volumetric fraction result (KVOL) in material buffer (+ posteriori checks)
      IF(mlw == 51)THEN
        DO imat=1,4
          kk = m51_n0phas + (imat-1)*m51_nvphas
          DO i=1,nel
            IF(ipart(i+nft) /= idp)cycle
            ratio = min(max(zero,kvol(imat,i)), one)
            uvar(i,1+kk)   = ratio
            uvar(i,23+kk)  = ratio
          ENDDO
        ENDDO
        DO i=1,nel
          IF(ipart(i+nft) /= idp)cycle
          kk  = m51_n0phas + (1-1)*m51_nvphas
          p1  = uvar(i,18+kk)
          kk  = m51_n0phas + (2-1)*m51_nvphas
          p2  = uvar(i,18+kk)
          kk  = m51_n0phas + (3-1)*m51_nvphas
          p3  = uvar(i,18+kk)
          kk  = m51_n0phas + (4-1)*m51_nvphas
          p4  = uvar(i,18+kk)
          sumvf=sum(kvol(1:nbsubmat,i))
          p   = kvol(1,i)*p1 + kvol(2,i)*p2 + kvol(3,i)*p3 + kvol(4,i)*p4
          uvar(i,4) = p
        ENDDO
      ELSEIF(mlw == 37)THEN
        DO i=1,nel
          IF(ipart(i+nft) /= idp)cycle
          uvar(i,4)  = kvol(1,i)
          uvar(i,5)  = kvol(2,i)
        ENDDO
      ELSEIF(mlw == 151)THEN
        gbuf => elbuf_tab(ng)%GBUF
        DO imat=1,multi_fvm%NBMAT
          lbuf => elbuf_tab(ng)%BUFLY(imat)%LBUF(1,1,1)
          DO i=1,nel
            IF(ipart(i+nft) /= idp)cycle
            lbuf%VOL(i) = kvol(imat,i) * gbuf%VOL(i)
          ENDDO
        ENDDO
      ENDIF
 
 
      RETURN
      END