C++ CppAD::log方法代码示例

bool PowInt(void)
{	bool ok = true;

	using CppAD::pow;
	using CppAD::exp;
	using CppAD::log;
	using namespace CppAD;


	// independent variable vector, indices, values, and declaration
	CPPAD_TESTVECTOR(AD<double>) U(1);
	U[0]     = 2.;
	Independent(U);

	// dependent variable vector and indices
	CPPAD_TESTVECTOR(AD<double>) Z(2);

	// dependent variable values
	Z[0]         = pow(U[0], 5);     // x = u^5
	Z[1]         = pow(U[0], -5);    // y = u^{-5}

	// create f: U -> Z and vectors used for derivative calculations
	ADFun<double> f(U, Z);
	CPPAD_TESTVECTOR(double) v( f.Domain() );
	CPPAD_TESTVECTOR(double) w( f.Range() );

	/*
	x_u = 5 * u^4
	y_u = - 5 * u^{-6}
	*/

	// check function values values
	double u = Value(U[0]);
	ok &= NearEqual(Z[0] , exp( log(u) * 5.),              1e-10 , 1e-10);
	ok &= NearEqual(Z[1] , exp( - log(u) * 5.),            1e-10 , 1e-10);

	// forward computation of partials 
	v[0] = 1.;
	w = f.Forward(1, v);
	ok &= NearEqual(w[0] , 5. * exp( log(u) * 4.),         1e-10 , 1e-10);
	ok &= NearEqual(w[1] , - 5. * exp( - log(u) * 6.),     1e-10 , 1e-10);

	return ok;
}

bool log10(void)
{	bool ok = true;

	using CppAD::log10;
	using CppAD::log;
	using namespace CppAD;

	// independent variable vector, indices, values, and declaration
	CPPAD_TESTVECTOR(AD<double>) U(1);
	size_t s = 0;
	U[s]     = 10.;
	Independent(U);

	// dependent variable vector, indices, and values
	CPPAD_TESTVECTOR(AD<double>) Z(2);
	size_t x = 0;
	size_t y = 1;
	Z[x]     = log10(U[s]);
	Z[y]     = log10(Z[x]);

	// define f : U -> Z and vectors for derivative calculations
	ADFun<double> f(U, Z);
	CPPAD_TESTVECTOR(double) v( f.Domain() );
	CPPAD_TESTVECTOR(double) w( f.Range() );

	// check values
	ok &= NearEqual(Z[x] , 1.,  1e-10 , 1e-10);
	ok &= NearEqual(Z[y] , 0.,  1e-10 , 1e-10);

	// forward computation of partials w.r.t. s
	double l10 = log(10.);
	v[s] = 1.;
	w    = f.Forward(1, v);
	ok &= NearEqual(w[x], 1./(U[s]*l10)         , 1e-10 , 1e-10); // dx/ds
	ok &= NearEqual(w[y], 1./(U[s]*Z[x]*l10*l10), 1e-10 , 1e-10); // dy/ds

	// reverse computation of partials of y
	w[x] = 0.;
	w[y] = 1.;
	v    = f.Reverse(1,w);
	ok &= NearEqual(v[s], 1./(U[s]*Z[x]*l10*l10), 1e-10 , 1e-10); // dy/ds

	return ok;
}

bool VecUnary(void)
{	
	using namespace CppAD;
	using CppAD::abs;
	using CppAD::sin;
	using CppAD::atan;
	using CppAD::cos;
	using CppAD::exp;
	using CppAD::log;
	using CppAD::sqrt;

	bool ok  = true;
	size_t n = 8;
	size_t i;

	CPPAD_TEST_VECTOR< AD<double> > X(n);
	VecAD<double>             Y(n);
	CPPAD_TEST_VECTOR< AD<double> > Z(n);


	for(i = 0; i < n; i++)
		X[i] = int(i);  // some compilers require the int here
	Independent(X);

	AD<double> j;

	j    = 0.;
	Y[j] = X[0]; 
	Z[0] = -Y[j];

	j    = 1.;
	Y[j] = X[1]; 
	Z[1] = sin( Y[j] );

	j    = 2.;
	Y[j] = X[2]; 
	Z[2] = abs( Y[j] );

	j    = 3.;
	Y[j] = X[3]; 
	Z[3] = atan( Y[j] );

	j    = 4.;
	Y[j] = X[4]; 
	Z[4] = cos( Y[j] );

	j    = 5.;
	Y[j] = X[5]; 
	Z[5] = exp( Y[j] );

	j    = 6.;
	Y[j] = X[6]; 
	Z[6] = log( Y[j] );

	j    = 7.;
	Y[j] = X[7]; 
	Z[7] = sqrt( Y[j] );

	
	ADFun<double> f(X, Z);
	CPPAD_TEST_VECTOR<double> x(n);
	CPPAD_TEST_VECTOR<double> z(n);

	for(i = 0; i < n; i++)
		x[i] = 2. / double(i + 1);
	x[7] = abs( x[7] );

	z    = f.Forward(0, x);

	ok  &= NearEqual(z[0],      - x[0],  1e-10, 1e-10);
	ok  &= NearEqual(z[1], sin( x[1] ),  1e-10, 1e-10);
	ok  &= NearEqual(z[2], abs( x[2] ),  1e-10, 1e-10);
	ok  &= NearEqual(z[3], atan(x[3] ),  1e-10, 1e-10);
	ok  &= NearEqual(z[4], cos( x[4] ),  1e-10, 1e-10);
	ok  &= NearEqual(z[5], exp( x[5] ),  1e-10, 1e-10);
	ok  &= NearEqual(z[6], log( x[6] ),  1e-10, 1e-10);
	ok  &= NearEqual(z[7], sqrt(x[7] ),  1e-10, 1e-10);

	return ok;
}

bool CondExp(void)
{	bool ok = true;

	using CppAD::isnan;
	using CppAD::AD;
	using CppAD::NearEqual;
	using CppAD::log; 
	using CppAD::abs;
	double eps  = 100. * CppAD::numeric_limits<double>::epsilon();
	double fmax = std::numeric_limits<double>::max();

	// domain space vector
	size_t n = 5;
	CPPAD_TESTVECTOR(AD<double>) X(n);
	size_t j;
	for(j = 0; j < n; j++)
		X[j] = 1.;

	// declare independent variables and start tape recording
	CppAD::Independent(X);

	AD<double> Sum  = 0.;
	AD<double> Zero = 0.;
	for(j = 0; j < n; j++)
	{	// if x_j > 0, add x_j * log( x_j ) to the sum
		Sum += CppAD::CondExpGt(X[j], Zero, X[j] * log(X[j]),    Zero);
	}

	// range space vector 
	size_t m = 1;
	CPPAD_TESTVECTOR(AD<double>) Y(m);
	Y[0] = Sum;

	// create f: X -> Y and stop tape recording
	CppAD::ADFun<double> f(X, Y);

	// vectors for arguments to the function object f
	CPPAD_TESTVECTOR(double) x(n);   // argument values
	CPPAD_TESTVECTOR(double) y(m);   // function values 
	CPPAD_TESTVECTOR(double) w(m);   // function weights 
	CPPAD_TESTVECTOR(double) dw(n);  // derivative of weighted function

	// a case where x[j] > 0 for all j
	double check  = 0.;
	for(j = 0; j < n; j++)
	{	x[j]   = double(j + 1); 
		check += x[j] * log( x[j] );
	}

	// function value 
	y  = f.Forward(0, x);
	ok &= NearEqual(y[0], check, eps, eps);

	// compute derivative of y[0]
	w[0] = 1.;
	dw   = f.Reverse(1, w);
	for(j = 0; j < n; j++)
		ok &= NearEqual(dw[j], log(x[j]) + 1., eps, eps); 

	// a case where x[3] is equal to zero
	check -= x[3] * log( x[3] );
	x[3]   = 0.;

	// function value 
	y   = f.Forward(0, x);
	ok &= NearEqual(y[0], check, eps, eps);

	// check derivative of y[0]
	f.check_for_nan(false);
	w[0] = 1.;
	dw   = f.Reverse(1, w);
	for(j = 0; j < n; j++)
	{	if( x[j] > 0 )
			ok &= NearEqual(dw[j], log(x[j]) + 1., eps, eps); 
		else
		{	// In this case computing dw[j] is computed using 
			//      log(x[j]) + x[j] / x[j] 
			// which has limit minus infinity but computes as nan.
			ok &= ( isnan( dw[j] ) || dw[j] <= -fmax );
		}
	}
	
	return ok;
}