c	INCLUDE 'mod_etat.f'
	
c********************************************************************

	PROGRAM test_etat

c	Auteur: P. Morel, Département J.D. Cassini, O.C.A.

c	EFF ne donne pas toujours des valeurs exactes pour
c	drox, dux, drotx et dutx car X et Y sont reliés par X+Y+Z=1
c	ce dont on ne semble pas tenir compte 

c-------------------------------------------------------------------------

	USE mod_donnees, ONLY: f_eos, ihe4, lit_nl, nbelem, nchim,
	1 nom_chemin, nom_elem, nom_etat, nom_fich2, x0, y0, z0
	USE mod_etat, ONLY: etat
	USE mod_kind
	USE mod_nuc, ONLY: nuc
			
	IMPLICIT NONE
	
	REAL (kind=dp), ALLOCATABLE, DIMENSION(:,:) :: jac
	REAL (kind=dp), ALLOCATABLE, DIMENSION(:) :: xchim, dxchim, ex
	REAL (kind=dp), DIMENSION(5) :: epsilo

	REAL (kind=dp), PARAMETER :: dd=1.d-6, unpdd=1.d0+dd
	REAL (kind=dp) :: be7, b8,et, ero, f17, hh, n13, o15, w_rot

	REAL (kind=dp) :: p, t, stor0, stor, dstor,
	1 ro, drop, drot, drox, u, dup, dut, dux,
	2 delta, deltap, deltat, deltax, cp,
	3 dcpp, dcpt, dcpx, gradad, dgradadp, dgradadt,
	4 dgradadx, alfa, beta, gamma1
	REAL (kind=dp) :: ro0, drop0, drot0, drox0, u0, dup0, dut0, dux0,
	2 delta0, deltap0, deltat0, deltax0, cp0,
	3 dcpp0, dcpt0, dcpx0, gradad0, dgradadp0, dgradadt0,
	4 dgradadx0, alfa0, beta0, gamma10	
		
	INTEGER :: ktest=1
		
	LOGICAL :: deriv
	
c-----------------------------------------------------------------------
	
2000	FORMAT(8es10.3)

c	lecture de test.don

	nom_fich2='test' ; CALL lit_nl(w_rot)
	
c appel d'initialisation pour tabulation des réactions nucléaires
c allocations fictives
	ALLOCATE(xchim(0),dxchim(0),jac(0,0),ex(0))
	CALL nuc(1.5d+07,1.5d+02,xchim,dxchim,jac,.FALSE.,0,
	1 epsilo,et,ero,ex,hh,be7,b8,n13,o15,f17)

c détermination des abondances initiales
	DEALLOCATE(xchim,dxchim,ex)
	ALLOCATE(xchim(nchim),dxchim(nchim),ex(nchim))
	CALL nuc(1.5d+07,1.5d+02,xchim,dxchim,jac,.FALSE.,1,
	1 epsilo,et,ero,ex,hh,be7,b8,n13,o15,f17)
	ktest=1
c	ktest=4
	B2: DO
	 SELECT CASE(ktest)
	 CASE(1)
	  p=1.754D+05 ; t=6.554D+03
	 CASE(2)
	  t=10.d0**(3.83644); p=10.d0**(5.34627)
	 CASE(3)
	  t=3.5d6 ; p=828473414360767.4
	 CASE(4)
	  t=6.554D+03 ; p=1.754D+05
	 CASE(5)
	  t=6.184E+06 ; p=2.711E+14
	 CASE(6)
	  t=4.050E+04 ; p=1.567E+05

c	 CASE()
c	  t= ; p=1.567E+05
	 CASE DEFAULT
	  EXIT B2
	 END SELECT
	 ktest=ktest+1	 
	 deriv=.TRUE.	
	 CALL etat(p,t,xchim,deriv,
	1 ro0,drop0,drot0,drox0,u0,dup0,dut0,dux0,
	2 delta0,deltap0,deltat0,deltax0,cp0,dcpp0,dcpt0,dcpx0,
	3 gradad0,dgradadp0,dgradadt0,dgradadx0,alfa0,beta0,gamma10)
	 WRITE(*,*)'p, t, X'
	 WRITE(*,2000)p,t,xchim(1)
	 WRITE(*,*)'ro, u, delta, cp, gradad, alfa, beta, gamma1'	 
	 WRITE(*,2000)ro0,u0,delta0,cp0,gradad0,alfa0,beta0,gamma10

	 deriv=.FALSE. ; WRITE(*,*)	
	 WRITE(*,*)
	 WRITE(*,3)p
3	 FORMAT('dérivées par rapport à p=',es10.3)
	 stor0=p ; stor=stor0*unpdd
	 IF(stor == 0.d0)stor=dd
	 dstor=stor-stor0
	 p=stor	
	 CALL etat(p,t,xchim,deriv,
	1 ro,drop,drot,drox,u,dup,dut,dux,
	2 delta,deltap,deltat,deltax,cp,dcpp,dcpt,dcpx,
	3 gradad,dgradadp,dgradadt,dgradadx,alfa,beta,gamma1)
	 p=stor0
	 WRITE(*,*)'drop,deltap,dgradadp'
	 WRITE(*,2000)drop0,(ro-ro0)/dstor,deltap0,(delta-delta0)/dstor,
	1 dgradadp0,(gradad-gradad0)/dstor
	 WRITE(*,*)'dup,dcpp'
	 WRITE(*,2000)dup0,(u-u0)/dstor,dcpp0,(cp-cp0)/dstor
	 PAUSE'dérivées par rapport à P'
	 WRITE(*,4)t
4	 FORMAT('dérivées par rapport à T=',es10.3)
	 stor0=t ; stor=stor0*unpdd
	 IF(stor == 0.d0)stor=dd
	 dstor=stor-stor0
	 t=stor
	 CALL etat(p,t,xchim,deriv,
	1 ro,drop,drot,drox,u,dup,dut,dux,
	2 delta,deltap,deltat,deltax,cp,dcpp,dcpt,dcpx,
	3 gradad,dgradadp,dgradadt,dgradadx,alfa,beta,gamma1)
	 t=stor0
	 WRITE(*,*)'drot,deltat,dgradadt'
	 WRITE(*,2000)drot0,(ro-ro0)/dstor,deltat0,(delta-delta0)/dstor,
	1 dgradadt0,(gradad-gradad0)/dstor
	 WRITE(*,*)'dut,dcpt'
	 WRITE(*,2000)dut0,(u-u0)/dstor,dcpt0,(cp-cp0)/dstor
	 PAUSE'dérivées par rapport à T'
	 
	 WRITE(*,5)xchim(1)
5	 FORMAT('dérivées par rapport à X=',es10.3)
	 stor0=xchim(1)
	 stor=stor0*unpdd
	 IF(stor == 0.d0)stor=dd
	 dstor=stor-stor0	 
	 xchim(1)=stor
	 CALL etat(p,t,xchim,deriv,
	1 ro,drop,drot,drox,u,dup,dut,dux,
	2 delta,deltap,deltat,deltax,cp,dcpp,dcpt,dcpx,
	3 gradad,dgradadp,dgradadt,dgradadx,alfa,beta,gamma1)
	 xchim(1)=stor
	 WRITE(*,*)'drox0,deltax,dgradadx'
	 WRITE(*,2000)drox0,(ro-ro0)/dstor,deltax0,(delta-delta0)/dstor,
	1 dgradadx0,(gradad-gradad0)/dstor
	 WRITE(*,*)'dux0,dcpx'
	 WRITE(*,2000)dux0,(u-u0)/dstor,dcpx0,(cp-cp0)/dstor
	 PAUSE'dérivées par rapport à X'
		
	ENDDO B2
	
	STOP
	
	END PROGRAM test_etat