ie_bound.F File Reference

#include "implicit_f.inc"

Go to the source code of this file.

Functions/Subroutines
function	ie_bound (pext, pm, c0, c1, c2, c3, c4, c5, e0)

Function/Subroutine Documentation

ie_bound()

function ie_bound	(		pext,
			pm,
			c0,
			c1,
			c2,
			c3,
			c4,
			c5,
			e0 )

Definition at line 29 of file ie_bound.F.

C-----------------------------------------------
C   I m p l i c i t   T y p e s
C-----------------------------------------------
#include      "implicit_f.inc"
C-----------------------------------------------
C   I N P U T   O U T P U T   A r g u m e n t s 
C-----------------------------------------------
      my_real ie_bound, pext,c0,c1,c2,c3,c4,c5,e0,pm
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      my_real  pcrit, p0, mu, c45,tmp, mu0, x, e0_ref, p0_ref, c1_ref
C
C----------+------+----------------------------------------------------C
C VAR      | SIZE | DEFINITION                                         C
C----------+------+----------------------------------------------------C
C IE_BOUND |  1   | LOWER LIMIT FOR RELATIVE ENERGY (NEGATIVE OR NULL) C
C PCRIT    |  1   | CRITERIA ENSURING EXISTENCE OF AN ENERGY EXTREMUM  C
C MU0      |  1   | EXPANSION STATE WHERE TO COMPUTE LOWER LIMIT       C
C          |      | AT THIS POINT PRESSURE IS ZERO AND THEN dE IS ZERO C
C P0       |  1   | SIMPLIFICATION FOR 0-DEGREE POLYNOMIAL COEFFICIENT C
C----------+------+----------------------------------------------------C
 
C-----------------------------------------------
C   S o u r c e   L i n e s
C-----------------------------------------------
      p0=pext+c0
      mu0=-one
      pcrit=zero
      c45=max(c4,c5)
      ie_bound=-infinity !default : no minimum 
      e0_ref=zero
      p0_ref=zero
      c1_ref=zero
 
      !===================!
      !   NO MATERIAL     !
      !===================!
      IF( (c0 == zero) .AND. (c1 == zero) .AND. (c2 == zero) .AND. (c3 == zero) .AND. (c4 == zero) .AND. (c5 == zero))THEN
        ie_bound = zero
        RETURN      
      END IF      
      
      !===================!
      !    P(0) = 0       !
      !  =>Einf = E0      ! 
      ! dE=-(dP+Pext)dV=0 !           
      !===================!
      IF(abs(p0+c4*e0) < em20)THEN
        ie_bound = e0
        RETURN      
      END IF 
       
      !===================!
      !    PM+PEXT>0      !
      ! IE_BOUND DOES NOT !
      ! EXISTS            !
      !===================!
      IF ((pm+pext) > em20) THEN
      !In this case extremum does not exists. Pressure (P+Pext) will never be zero, then dE>0 is never zero and E has no extremum. limit E(->-1) = -INFINITY
        ie_bound=-infinity         
        RETURN
      END IF
      
      
      
 
      !Otherwise for following cases : {PM+PEXT>0}
  
 
      !===================!
      !    PERFECT GAS    !
      !===================!
      IF((c2 == zero) .AND. (c3 == zero) .AND. (p0 == c1) .AND. (c4 /= zero) .AND. (c4 == c5))THEN
       !this condition describes all formulation possible for a perfect gas, included relative pressure and relative energy.
        IF (c1 == 0)THEN
        !total energy
          ie_bound = zero
        ELSE
        !relative energy
          ie_bound = -p0/c4
        END IF
        RETURN
      END IF     
 
      !======================!
      ! INCOMPRESSIBLE GAS   !
      ! LINEAR ELASTIC SOLID !      
      !======================! 
      IF((c2 == zero) .AND. (c3 == zero) .AND. (c4 == zero) .AND. (c5 == zero) .AND. (c1 /= zero))THEN
        IF (p0 < c1)THEN
             x=-p0/c1  !C1-P0>0 then x  > -1, extremum is well defined
             ie_bound = -1/(1+x)*p0+c1/(1+x)+c1*log(1+x)+e0+p0-c1
        ELSE
             ie_bound = -infinity  !limit E(mu=-1) = -sgn(C1) * infinity, with C1>0        
           END IF        
           RETURN
      END IF                
 
      !======================!
      ! MIE-GRUNEISEN        !
      ! GAMMA 1st order      !      
      !======================! 
      !  extremum estimated using one from case gamma constant
 
      !======================!
      ! MIE-GRUNEISEN        !
      ! GAMMA CONSTANT       !      
      !======================! 
      IF((c1 > zero) .AND. (c4 >= c5) .AND. (c5 > zero))THEN
      !similar to E0=0
      !P=P0      +        C1x + C2x+C3x+(C4+C5x)(dE+E0)
      !P=P0+C4E0 + (C1+C5E0)x + C2x+C3x+(C4+C5x) dE    : solving dE then puis on ajoute E0
        e0_ref=e0
        c1_ref=c1
        p0_ref=p0
        IF(e0/=zero)THEN
          e0=zero
          p0=p0+c4*e0_ref
          c1=c1+c4*e0_ref
        END IF
      
      ! existing condition for C4/=1 et C4/=2
      IF((c4 /= 1) .AND. (c4 /= 2))THEN
        !diverging test (expansion)
        IF(c1-p0 < zero)THEN
          ie_bound=-infinity
          GOTO 100
        END IF
        !diverging test (compression)
        x=(-4*c2+2*c1+2*c4*p0+2*c4*c2-3*c4*c1-3*c4**2*p0+6*c3+c4**2*c1+c4**3*p0)/(2+c4**2-3*c4)
        IF(x < zero)THEN
          ie_bound=-infinity
          GOTO 100
        END IF
        !expansion extremum
        IF(p0 >= zero)THEN
          IF (c1 == p0)THEN
            ie_bound=-c1/c4   !limite en mu=-1
            GOTO 100
          ELSE !extremum in  ] -1 et 0 [
            x=(c4*p0+c1)/(c1-p0)  !necessarily positive
            x=exp(-log(x)/(1+c4))-1
            ie_bound=(-(1+x)**(-1-c4)*(c4*p0+c1+c1*x+c4*c1*x) / c4/(1+c4)+(c4*p0+c1)/c4/(1+c4))*(1+x)**c4
            GOTO 100
            END IF !(c1==p0)
          END IF !(P0>=ZERO)
        !extremum in compression
        IF(p0 < zero)THEN
        !extremum estimation in compression (C2,C3 needed) :
        ie_bound=-infinity ! ameliorer
          !find an estimator
          !fixed point method ?
          GOTO 100
        END IF
        
      ELSEIF((c4 == one) .OR. (c4 == two))THEN
        !diverging test in expansion
        IF(c1-p0 < zero)THEN
          ie_bound=-infinity
          GOTO 100
        END IF
        !diverging test in compression
          !diverge if dominant coefficient is negative
        IF (c3 < zero)THEN
          ie_bound = -infinity
          GOTO 100
        END IF
          !diverge if dominant coefficient is negative
        IF(c3 == zero)THEN
          IF (c2 < zero)THEN
            ie_bound=-infinity
            GOTO 100
          END IF
        END IF
        ! otherwise extremum does exist
        ! expansion if P0>0
        IF(p0 >= zero)THEN
           !C1**2-P0^2 is here necessarily positive because C1-P0>0 C1>0 P0>0
          x=-(c1+p0-sqrt(c1**2-p0**2))/(c1+p0)
          ie_bound=(1+x)*(half*c1+half*p0)-half*(c1+2*c1*x+p0)/(1+x)
          GOTO 100
        END IF
        ! compression
        IF(p0 < zero)THEN
        !extremum estimation in compression (C2,C3 needed) :
        ie_bound=-infinity ! to improve
          !find an estimator
          !fixed point method ?
          GOTO 100
        END IF
      END IF !condition for C4 (1 or 2)
      END IF !mie gruneisen gamma constant
      
  100 IF (e0_ref /= zero)THEN
        IF(ie_bound /= -infinity)ie_bound=ie_bound+e0_ref
        c1=c1_ref 
        e0=e0_ref
      END IF
      
      RETURN