exponential_ function(x,c=log(10),k1=200,k2=100,tol=1e-4,M=1e6,pi1=2/3,pi2=1/3,beta1=2,beta2=0.5)
{

# Test for homogeneity in a two-component exponential scale mixture. 
# Results are given for the likelihood ratio test, the modified likelihood 
# ratio test (Chen, Chen, and Kalbfleisch), the D-Test based on unpenalized 
# maximum likelihood estimation with and without weighting functions, and 
# the D-Test based on penalized maximum likelihood estimation with and
# without weighting functions.
#
# Scenario: 
#  	H_0: There exists beta_0 such that
#		p_1 f(x | beta_1) + p_2 f(x | beta_2)
#	     	may be written as f(x | beta_0).
#	H_1: There is no such  beta_0.
#
# Required Input:
# 	x  is a vector of data.
#
# Optional Inputs:
#	c  is a penalty parameter governing the strength of 
#	        a  log(4 p_1 p_2)  term that is attached to
#		the log likelihood function in a penalized maximum
#		likelihood estimation scheme.
#  	k1  is the maximum number of iterations for computing
#		unpenalized maximum likelihood estimates using the 
#		expectation maximization algorithm.
#	k2  is the maximum number of iterations for computing
#		penalized maximum likelihood estimates using the
#		expectation maximization algorithm, with the
#		unpenalized maximum likelihood estimates taken as
#		the starting values.
#	tol  is the largest sum of relative changes in parameter
#		estimates for which the expectation maximization
#		algorithm will terminate short of the maximum
#		number of iterations.
#	M  defines the parameter space for  beta  as  [1/M, M].
#     pi1, pi2, beta1, beta2   are initial values for the expectation
#				   maximization algorithm used to obtain
#				   unpenalized maximum likelihood estimates.
#				   These initial values can be obtained
#				   by method of moments estimation or
#				   by visual inspection of an undersmoothed
#				   histogram of the data.
#
# Outputs:
# 	lrt   (Unpenalized) Likelihood Ratio Test:
#		critical values must be estimated via simulation
#	mlrt	Modified (Penalized) Likelihood Ratio Test 
#		developed by Chen, Chen, and Kalbfleisch:
#		for a large-sample level alpha test, compare to 
#		the upper 2*alpha quantile of a chi square distribution
#		on one degree of freedom
#	pmlrt   large-sample p-value for Modified Likelihood Ratio Test
#		('LARGE' = 0.30 or more)
#	dk	D-Test based on unpenalized maximum likelihood estimation:
#		critical values must be estimated via simulation;
#		if you are content with a simulation size of 10000,
#		you can access our simulation results through  exponentialc.S
#     dk.w1     with weighting function w_1(x) = x
#     dk.w2     with weighting function w_2(x) = x^2
#     dkpen   D-Test based on penalized maximum likelihood estimation:
#		for large-sample level alpha test, compare to the upper 
#		2*alpha quantile of a chi square distribution on one degree 
#		of freedom multiplied by (3/16) hat beta_0 / length(x)
#	pdkpen  large-sample p-value for D-Test based on 
#		penalized maximum likelihood estimation
#		('LARGE' = 0.30 or more)
#     dkpen.w1  and  pdkpen.w1   with weighting function w_1(x) = x:
#	 	compare to the upper 2*alpha quantile of a chi square 
#		distribution on one degree of freedom multiplied by
#		(3/32) / length(x)
#	dkpen.w2  and  pdkpen.w2   with weighting function w_2(x) = x^2:
#		compare to the upper 2*alpha quantile of a chi square
#		distribution on one degree of freedom multiplied by
#		(5/32 hat beta_0) / length(x)
#	nit1    number of expectation-maximization iterations to obtain
#		the unpenalized maximum likelihood estimates of the 
#		mixture parameters
#	beta1, beta2, pi1, pi2   unpenalized maximum likelihood estimates
#     nit2    number of expectation-maximization iterations to obtain
#		the penalized maximum likelihood estimates of the
#		mixture parameters
#	beta1pen, beta2pen, pi1pen, pi2pen   penalized maximum likelihood estimates
#
# Caveat:  
#     The  p-values  are based on asymptotic theory that, although presumed
#	applicable, has not been rigorously justified in the exponential scale 
#	case.  The principal difficulty is that one of the regularity conditions 
#	(uniform boundedness) is violated in the exponential scale case, unless 
#	the parameter space is severely restricted.


n_ length(x)

# n  is the "sample size"

beta0_ 1/mean(x)


beta1pen_ beta1
beta2pen_ beta2
pi1pen_ pi1
pi2pen_ pi2
lrt_ 0
mlrt_ 0
nit1_ 0
nit2_ 0 

storage.mode(x)_"double"

storage.mode(beta1)_"double"
storage.mode(beta2)_"double"
storage.mode(pi1)_"double"
storage.mode(pi2)_"double"

storage.mode(beta1pen)_"double"
storage.mode(beta2pen)_"double"
storage.mode(pi1pen)_"double"
storage.mode(pi2pen)_"double"

storage.mode(c)_"double"
storage.mode(M)_"double"
storage.mode(tol)_"double"

storage.mode(beta0)_"double"
storage.mode(lrt)_"double"
storage.mode(mlrt)_"double"


# Here we compute parameter estimates under both
# unpenalized and penalized maximum likelihood schemes.


w_.Fortran("exponential", 
		x=as.double(x),
	        n=as.integer(n), 
	        nit1=as.integer(nit1), 	        
		  nit2=as.integer(nit2), 
		  c=as.double(c),
		  M=as.double(M),
		  k1=as.integer(k1),
		  k2=as.integer(k2),
		  tol=as.double(tol),
		pi1=as.double(pi1), 
		pi2=as.double(pi2), 
		beta1=as.double(beta1), 
		beta2=as.double(beta2),
		pi1pen=as.double(pi1pen), 
		pi2pen=as.double(pi2pen), 
		beta1pen=as.double(beta1pen), 
		beta2pen=as.double(beta2pen),
		beta0=as.double(beta0),
		lrt=as.double(lrt),
		mlrt=as.double(mlrt))

pi1_ w$pi1
pi2_ w$pi2
beta1_ w$beta1
beta2_ w$beta2
lrt_ w$lrt
pi1pen_ w$pi1pen
pi2pen_ w$pi2pen
beta1pen_ w$beta1pen
beta2pen_ w$beta2pen
mlrt_ w$mlrt
nit1_ w$nit1
nit2_ w$nit2


# Here we compute the D-test statistics based on 
# the parameter estimates we have.

# unpenalized maximum likelihood

dk0_ beta0/2 + pi1*pi1*beta1/2 + pi2*pi2*beta2/2 
dk1_ 2*pi1*beta0*beta1/(beta0+beta1)+2*pi2*beta0*beta2/(beta0+beta2)
dk2_ 2*pi1*pi2*beta1*beta2/(beta1+beta2)

dk_ dk0 - dk1 + dk2


# weighting function w_1(x) = x

dknew0_ 1/4 + pi1*pi1/4 + pi2*pi2/4
dknew1_ 2*pi1*beta0*beta1/((beta0+beta1)^2)+2*pi2*beta0*beta2/((beta0+beta2)^2)
dknew2_ 2*pi1*pi2*beta1*beta2/((beta1+beta2)^2)

dk.w1_ dknew0 - dknew1 + dknew2


# weighting function w_2(x) = x^2

dk00_ -1*-1*beta0*beta0/(beta0+beta0)^3 
dk01_ -1*pi1*beta0*beta1/(beta0+beta1)^3
dk02_ -1*pi2*beta0*beta2/(beta0+beta2)^3

dk10_ pi1*-1*beta1*beta0/(beta1+beta0)^3 
dk11_ pi1*pi1*beta1*beta1/(beta1+beta1)^3
dk12_ pi1*pi2*beta1*beta2/(beta1+beta2)^3

dk20_ pi2*-1*beta2*beta0/(beta2+beta0)^3 
dk21_ pi2*pi1*beta2*beta1/(beta2+beta1)^3
dk22_ pi2*pi2*beta2*beta2/(beta2+beta2)^3

dk.w2_ 2*(dk00+dk01+dk02+dk10+dk11+dk12+dk20+dk21+dk22)


# penalized maximum likelihood

dk0_ beta0/2 + pi1pen*pi1pen*beta1pen/2 + pi2pen*pi2pen*beta2pen/2 
dk1_ 2*pi1pen*beta0*beta1pen/(beta0+beta1pen)+2*pi2pen*beta0*beta2pen/(beta0+beta2pen)
dk2_ 2*pi1pen*pi2pen*beta1pen*beta2pen/(beta1pen+beta2pen)

dkpen_ dk0 - dk1 + dk2


# weighting function w_1(x) = x

dknew0_ 1/4 + pi1pen*pi1pen/4 + pi2pen*pi2pen/4
dknew1_ 2*pi1pen*beta0*beta1pen/((beta0+beta1pen)^2)+2*pi2pen*beta0*beta2pen/((beta0+beta2pen)^2)
dknew2_ 2*pi1pen*pi2pen*beta1pen*beta2pen/((beta1pen+beta2pen)^2)

dkpen.w1_ dknew0 - dknew1 + dknew2


# weighting function w_2(x) = x^2

dk00_ -1*-1*beta0*beta0/(beta0+beta0)^3 
dk01_ -1*pi1pen*beta0*beta1pen/(beta0+beta1pen)^3
dk02_ -1*pi2pen*beta0*beta2pen/(beta0+beta2pen)^3

dk10_ pi1pen*-1*beta1pen*beta0/(beta1pen+beta0)^3 
dk11_ pi1pen*pi1pen*beta1pen*beta1pen/(beta1pen+beta1pen)^3
dk12_ pi1pen*pi2pen*beta1pen*beta2pen/(beta1pen+beta2pen)^3

dk20_ pi2pen*-1*beta2pen*beta0/(beta2pen+beta0)^3 
dk21_ pi2pen*pi1pen*beta2pen*beta1pen/(beta2pen+beta1pen)^3
dk22_ pi2pen*pi2pen*beta2pen*beta2pen/(beta2pen+beta2pen)^3

dkpen.w2_ 2*(dk00+dk01+dk02+dk10+dk11+dk12+dk20+dk21+dk22)



pmlrt_ 'LARGE'
if (mlrt > .2749959) pmlrt_ 0.5-0.5*pchisq(mlrt,1)

pdkpen_ 'LARGE'
sdkpen_ (16/3)*(1/beta0)*n*dkpen
if (sdkpen > .2749959) pdkpen_ 0.5-0.5*pchisq(sdkpen,1)

pdkpen.w1_ 'LARGE'
sdkpen.w1_ (32/3)*n*dkpen.w1
if (sdkpen.w1 > .2749959) pdkpen.w1_ 0.5-0.5*pchisq(sdkpen.w1,1)

pdkpen.w2_ 'LARGE'
sdkpen.w2_ (32/5)*beta0*n*dkpen.w2
if (sdkpen.w2 > .2749959) pdkpen.w2_ 0.5-0.5*pchisq(sdkpen.w2,1)

return(lrt,mlrt,pmlrt,dk,dk.w1,dk.w2,dkpen,pdkpen,dkpen.w1,pdkpen.w1,dkpen.w2,pdkpen.w2,nit1,beta1,beta2,pi1,pi2,nit2,beta1pen,beta2pen,pi1pen,pi2pen)
result
}