c*********************************************************************

	SUBROUTINE opa_compton(xchim,t,ro,kap,dkapt,dkapro,dkapx)

c routine public du module mod_opa
c Opacités free-free & Compton utilisées pour 10^7°K < T < 3.310^9°K

c Utilisation de la tabulation de Cox & Giuli §16.a, table 16.1
c prolongation gaunt=0 à 3.310^9°K

c entrées :
c        xchim(1)=X : comp. chim. par gramme
c        t : température K
c        ro : densité cgs

c sorties :
c        kappa : opacité gr / cm2)
c        dkapdt : kappa / d t
c        dkapdr : kappa / d densité
c        dkapdx : kappa / d xchim(1)

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

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

	USE mod_donnees, ONLY : amu, clight, echarg, eve, ihe4,kbol, langue,
	1 me, nchim, nucleo, pi, zi
	USE mod_kind
	USE mod_numerique, ONLY : bsp1dn, no_croiss
	USE mod_variables, ONLY : sortie

	IMPLICIT NONE

        REAL (kind=dp), INTENT(in), DIMENSION(nchim) :: xchim
	REAL (kind=dp), INTENT(in) :: ro, t
	REAL (kind=dp), INTENT(out) ::	kap, dkapt, dkapro, dkapx

	INTEGER, PARAMETER :: mkt=3, nkt=12
	REAL (kind=dp), DIMENSION(1,nkt) :: g
	REAL (kind=dp), SAVE, DIMENSION(nkt+mkt) :: at
	REAL (kind=dp),	ALLOCATABLE, DIMENSION(:) :: g1
	REAL (kind=dp), SAVE, DIMENSION(nkt) :: a
        REAL (kind=dp), SAVE, ALLOCATABLE, DIMENSION(:) :: zsa, z2sa
	REAL (kind=dp), DIMENSION(1) :: corr, dcorr

	REAL (kind=dp), SAVE :: cte1, cte2, cte3, cte4
	REAL (kind=dp) :: kap_ff, kt, xh, xh2, z

	INTEGER, SAVE :: knot, l=1

	LOGICAL, SAVE :: init=.TRUE.

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

2000	FORMAT(8es10.3)

c----------------------initialisations début------------------------
	IF(init)THEN
	 init=.FALSE.
	 cte1=8.d0*pi/3.d0*((echarg/clight)**2/me)**2/amu
	 IF(ALLOCATED(zi))THEN
	  cte2=cte1*zi(1)**2/nucleo(1)
	 ELSE
	  cte2=cte1/nucleo(1)
	 ENDIF
	 cte3=kbol/eve*1.d-3
	 cte4=3.76d22*2.d0	!facteur 2 car <Z/A>~1/2(1+X) de Cox & Guili 16.90

c initialisations, températures en Kev (1Kev ~ 10 T6)
	 ALLOCATE(g1(nkt))
	 a=(/ 1.d0, 2.d0, 4.d0, 6.d0, 9.d0, 14.d0, 20.d0, 30.d0, 50.d0,
	1 80.d0, 125.d0, 330.d0 /)
	 g1=(/ 1.d0, 0.95055d0,0.9044d0, 0.8626d0, 0.8087d0, 0.7279d0, 0.6525d0,
	1 0.5590d0, 0.4408d0, 0.3411d0, 0.2579d0, 0.00d0 /)
	 g=reshape(g1,shape(g)) ; DEALLOCATE(g1)
	 CALL bsp1dn(1,g,a,at,nkt,mkt,knot,.FALSE.,a(1),l,corr,dcorr)
	 IF(no_croiss)PRINT*,'Pb. at 1 in opa_compton'

c écritures
	 z=SUM(xchim(ihe4+1:nchim))
	 SELECT CASE(langue)
	 CASE('english')
	  WRITE(*,1001)t,z ; WRITE(2,1001)t,z
1001	  FORMAT(/,'Use of ff+Compton opacities, at least once, T=',es10.3,
	1' > 4d7 and/or Z=',es10.3,' > 0.1'/)
	 CASE DEFAULT
	  WRITE(*,1)t,z ; WRITE(2,1)t,z
1	  FORMAT(/,'Utilisation des opacités ff+Compton au moins 1 fois, T=',es10.3,
	1' > 4d7 et/ou Z=',es10.3,' > 0.1'/)
	 END SELECT

c Zi / Ai et Zi^2 / Ai
	 ALLOCATE(zsa(nchim), z2sa(nchim))
	 zsa=zi/nucleo ; z2sa=zi*zsa

	ENDIF
c----------------------initialisations fin------------------------

c < Z / A > et < Z^2 / A >
	xh=DOT_PRODUCT(xchim,zsa) ; xh2=DOT_PRODUCT(xchim,z2sa)

c diffusion Thomson
	kap=cte1*xh

c correction Compton
	kt=MAX(a(1),MIN(cte3*t,a(nkt)))		!Kt en Kev
	CALL bsp1dn(1,g,a,at,nkt,mkt,knot,.TRUE.,kt,l,corr,dcorr)
	IF(no_croiss)THEN
	 PRINT*,'ARRET, Pb. at 2 in opa_compton, T, xchim'
	 WRITE(*,2000)t ; WRITE(*,2000)xchim ; CALL sortie('opa_compton')
	ENDIF

c dérivée / T
	dkapt=kap*dcorr(1)*cte3

c diffusion Compton
	kap=kap*corr(1)

c dérivée / ro
	dkapro=0.d0

c dérivée / X
	dkapx=cte2*corr(1)

c opacité free-free selon Kippenhahan & Weigert 17.5-17.6, Cox & Guili 16.95
c 1+X et X+Y+B sont respectivement remplacés par 2 xh=< Z/A > et xh2= < Z^2 / A >
	kap_ff=cte4*xh*xh2*ro/SQRT(t**7)
	dkapro=dkapro+kap_ff/ro
	dkapt=dkapt-3.5d0*kap_ff/t
	dkapx=dkapx+kap_ff*(zsa(1)/xh+z2sa(1)/xh2)
	kap=kap+kap_ff

	RETURN

	END SUBROUTINE opa_compton