sesa20.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/.
!||====================================================================
!||    sesa20     ../engine/source/materials/mat/mat026/sesa20.F
!||--- called by ------------------------------------------------------
!||    mmain      ../engine/source/materials/mat_share/mmain.F90
!||--- calls      -----------------------------------------------------
!||    mintp_re   ../common_source/eos/mintp_re.F
!||    mintp_rt   ../common_source/eos/mintp_rt.F
!||====================================================================
      SUBROUTINE sesa20(
     1   PM,      EINT,    RHO,     TEMP,
     2   Z,       QOLD,    SESAME,  SOUND,
     3   XK,      VOLN,    MAT,     C,
     4   DV2,     POLD,    SSP,     RHO0,
     5   P01,     NEL)
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      "param_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER, INTENT(IN) :: NEL
      my_real
     .   PM(NPROPM,*), EINT(*), RHO(*), TEMP(*), Z(*), QOLD(*),
     .   SESAME(*), SOUND(*), XK(*), VOLN(MVSIZ), C(*), DV2(*), POLD(*), SSP(*),
     .   RHO0(*), P01(*)
      INTEGER MAT(*)
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER I, MX, NR, NT, IDR, IDT, IDP, IDE, IDQ
      my_real
     .   dpdm(mvsiz), 
     .   espe(mvsiz), e01(mvsiz), 
     .   dtde, dpdr, sph, rl, tl, zl, dxdr, zz, deltat, xlamb, 
     .   stefan, opac, rossel, xkr
      real*8 ne, na, atom
C-----------------------------------------------
      DATA na /6.0225e+23/
C---------------------------------
      DO i=1,nel
        zl = zero
        dxdr = zero
        mx = mat(i)
        dpdm(i) = four_over_3 * pm(22,mx)
        e01(i) =eint(i)-(pold(i)+qold(i))*dv2(i)
        e01(i) = max(zero,e01(i))
        espe(i)=e01(i)/max(em30,voln(i)*rho(i))
        nr  = nint(pm(25,mx))
        nt  = nint(pm(26,mx))
        idr = nint(pm(27,mx))
        idt = idr + nr
        idp = idt + nt
        ide = idp + nr * nt
 
        CALL mintp_re(sesame(idr),nr,sesame(idt),nt,sesame(ide),rho(i),temp(i),espe(i),dtde)
        dtde = max(zero,dtde)
        CALL mintp_rt(sesame(idr),nr,sesame(idt),nt,sesame(idp),rho(i),temp(i),p01(i),dpdr)
 
        dpdm(i)=dpdm(i)+rho0(i)*dpdr
        sph = rho(i)/max(em15,dtde)
        !------------------------------------
        !     IONISATION
        !------------------------------------
        nr  = nint(pm(28,mx))
        nt  = nint(pm(29,mx))
        idr = nint(pm(30,mx))
        idt = idr + nr
        idq = idt + nt
        rl = log10(rho(i))
        tl = log10(temp(i))
        CALL mintp_rt(sesame(idr),nr,sesame(idt),nt,sesame(idq),rl,tl,zl,dxdr)
        IF(zl <= -20)THEN
          !machine epsilon
          z(i)=em20
        ELSEIF(zl >= 20)THEN
          !machine infinity
          z(i)=ep20
        ELSE
          z(i)=exp(zl*log(ten)) !TEN**ZL
        ENDIF  
 
        !(CONDUNCTIVITE MAX POUR CALCUL DE DT)
        zz = max(em10,z(i))
        deltat= threep44 * zep26 * log(zz) / zz
        IF(deltat > zero)THEN
          deltat= one / ( one +  deltat)
        ELSE
          deltat= zero
        ENDIF
        atom  = pm(37,mx)
        ne    = rho(i)*na*zz/atom
        xlamb = pm(36,mx)*temp(i)**three_half/sqrt(ne)
        xlamb = max(one,xlamb)
        xlamb = max(em10, log(xlamb))
        xk(i) = zep4*deltat*pm(35,mx) * temp(i)**twop5 / (zz*xlamb)
        stefan = pm(51,mx)
        IF(stefan > zero.AND.temp(i) >= ep04)THEN
          nr     = nint(pm(48,mx))
          nt     = nint(pm(49,mx))
          idr    = nint(pm(50,mx))
          idt = idr + nr
          idq = idt + nt
          rl = log10(rho(i))
          tl = log10(temp(i))
          CALL mintp_rt(sesame(idr),nr,sesame(idt),nt,sesame(idq),rl,tl,opac,dxdr)
          opac = ten**opac
          rossel = one / ( rho(i) * opac)
          xkr = sixteen * stefan * temp(i)**3 * rossel * third
          xk(i) = xk(i) + xkr
        ENDIF
        xk(i) = min(xk(i),pm(51,mx))
      ENDDO
      
      !------------------------------------
      !     VITESSE DU SON
      !------------------------------------
      DO i=1,nel
        ssp(i)=sqrt(abs(dpdm(i))/rho0(i))
        sound(i) = ssp(i)
      ENDDO
C
      RETURN
      END