gradienttest.cc File Reference

#include "mbconfig.h"
#include <ctime>
#include <cstdlib>
#include <iostream>
#include <iomanip>
#include <functional>
#include <typeinfo>
#include <cmath>
#include "clock_time.h"
#include "gradient.h"

Include dependency graph for gradienttest.cc:

Go to the source code of this file.

Macros
#define	M_PI   3.14159265358979
#define	GRADIENT_DEBUG   1

Functions
void	tic ()
void	tic (doublereal &dTime)
doublereal	toc ()
doublereal	random1 ()
template<typename T >
bool	bCompare (const T &a, const T &b, doublereal dTolRel=0.)
template<index_type N_SIZE>
bool	bCompare (const Gradient< N_SIZE > &a, const Gradient< N_SIZE > &b, doublereal dTolRel=0.)
template<int N_SIZE>
void	testRangeVector ()
template<int N_SIZE>
void	testMapVector ()
template<typename T >
void	func (const T &u, const T &v, const T &w, doublereal e, T &f)
template<typename T >
void	func_tmp (const T &u, const T &v, const T &w, doublereal e, T &f)
doublereal	sec (doublereal x)
void	funcDeriv (index_type N, doublereal u, const doublereal ud[], doublereal v, const doublereal vd[], doublereal w, const doublereal wd[], doublereal e, doublereal &f, doublereal fd[])
template<int N_SIZE>
void	testGradient (index_type N)
template<index_type N_SIZE>
void	testGradient2 ()
void	testGradientCopy (int N)
template<index_type N_SIZE>
void	testDifferentDofMaps (index_type N)
template<index_type N_SIZE>
void	testGradientLin (int N)
int	main (int argc, char *argv[])

Variables
int	NLoops = 1
int	NLoopsAss = 1
doublereal	dStartTime

Macro Definition Documentation

#define GRADIENT_DEBUG   1

Definition at line 71 of file gradienttest.cc.

Referenced by main().

#define M_PI   3.14159265358979

Definition at line 67 of file gradienttest.cc.

Referenced by PlaneHingeJoint::AfterConvergence(), PlaneRotationJoint::AfterConvergence(), Wheel4::AfterConvergence(), AxialRotationJoint::AfterConvergence(), PlanePinJoint::AfterConvergence(), CyclocopterPolimi::AfterConvergence(), ModLugreFriction::alpha(), ModLugreFriction::alphad_v(), ModLugreFriction::alphad_z(), TheodorsenAeroData::AssJac(), RotorTrimBase::AssRes(), DynamicInflowRotor::AssRes(), PetersHeRotor::AssRes(), CyclocopterPolimi::AssRes(), c81_aerod2(), c81_aerod2_u(), ScrewJointSh_c::ComputePitchAngle(), Control_valve::Control_valve(), CosineDriveCaller::CosineDriveCaller(), Motor::dGetPhiElectric(), PlaneHingeJoint::dGetPrivData(), PlaneRotationJoint::dGetPrivData(), AxialRotationJoint::dGetPrivData(), PlanePinJoint::dGetPrivData(), do_eig(), Dynamic_control_valve::Dynamic_control_valve(), Mat3x3::EigSym(), FourierSeriesDriveCaller::FourierSeriesDriveCaller(), getbladeparams(), GlauertRotor::GetInducedVelocity(), ManglerRotor::GetInducedVelocity(), CyclocopterPolimi::GetInducedVelocity(), MBDynParser::GetMatR2vec(), getrotorparams(), UniformRotor::Init(), GlauertRotor::Init(), ManglerRotor::Init(), DynamicInflowRotor::Init(), PetersHeRotor::Init(), main(), MatManip_test(), mp_ctg_t(), mp_tan_t(), Wheel2::Output(), GenericAerodynamicForce::Output(), Wheel4::Output(), Pressure_flow_control_valve::Pressure_flow_control_valve(), ReadGenericAerodynamicForce(), AerodynamicOutput::SetData(), SineDriveCaller::SineDriveCaller(), Unwrap(), and Wheel2::Wheel2().

Function Documentation

template<typename T >

bool bCompare	(	const T &	a,
		const T &	b,
		doublereal	dTolRel = 0.
	)

Definition at line 91 of file gradienttest.cc.

Referenced by func_tmp(), testGradient(), and testGradient2().

                                                               {
        doublereal dTolAbs = std::max<T>(1., std::max<T>(std::abs(a), std::abs(b))) * dTolRel;

        if (!(std::abs(a - b) <= dTolAbs)) {
                std::cerr << "a = " << a << " != b = " << b << std::endl;
                return false;
        }

        return true;
}

template<index_type N_SIZE>

bool bCompare	(	const Gradient< N_SIZE > &	a,
		const Gradient< N_SIZE > &	b,
		doublereal	dTolRel = 0.
	)

Definition at line 103 of file gradienttest.cc.

References grad::Gradient< N_SIZE >::dGetDerivativeLocal(), grad::Gradient< N_SIZE >::dGetValue(), grad::Gradient< N_SIZE >::iGetEndIndexLocal(), and grad::Gradient< N_SIZE >::iGetStartIndexLocal().

                                                                                             {
        doublereal dTolAbs = std::max(1., std::max(std::abs(a.dGetValue()), std::abs(b.dGetValue()))) * dTolRel;

        if (std::abs(a.dGetValue() - b.dGetValue()) > dTolAbs) {
                std::cerr << "a = " << a.dGetValue() << " != b = " << b.dGetValue() << std::endl;
                return false;
        }

        index_type iStart = std::min(a.iGetStartIndexLocal(), b.iGetStartIndexLocal());
        index_type iEnd = std::max(a.iGetEndIndexLocal(), b.iGetEndIndexLocal());

        for (index_type i = iStart; i < iEnd; ++i) {
                doublereal dTolAbs = std::max(1., std::max(std::abs(a.dGetDerivativeLocal(i)), std::abs(b.dGetDerivativeLocal(i)))) * dTolRel;

                if (std::abs(a.dGetDerivativeLocal(i) - b.dGetDerivativeLocal(i)) > dTolAbs) {
                        std::cerr << "ad(" << i << ") = " << a.dGetDerivativeLocal(i) << " != bd(" << i << ") = " << b.dGetDerivativeLocal(i) << std::endl;
                        return false;
                }
        }

        return true;
}

Here is the call graph for this function:

template<typename T >

void func	(	const T &	u,
		const T &	v,
		const T &	w,
		doublereal	e,
		T &	f
	)

Definition at line 309 of file gradienttest.cc.

References grad::cos(), grad::pow(), grad::sin(), and grad::tan().

Referenced by MathParser::StaticNameSpace::GetFunc(), MathParser::StaticNameSpace::IsFunc(), MBDynErrBase::MBDynErrBase(), ParseScalarFunction(), ReadCL1D(), ReadCL3D(), ReadCL6D(), ReadDC1D(), ReadDC3D(), ReadDC3x3D(), ReadDC6D(), ReadDC6x6D(), ReadDescription(), ReadDriveCallerData(), ReadDriveData(), ReadGustData(), ReadTimeStepData(), grad::LocalDofMap::Reset(), testGradient(), testGradient2(), ErrIndexOutOfRange::WriteMsg(), and MathParser::StaticNameSpace::~StaticNameSpace().

                                                                  {
    f = ((((3 * u + 2 * v) * (u - v) / (1 - w) * pow(u/w, v - 1) * sin(v) * cos(w) * (1 - tan(-w + v))) * e + 1. - 11. + 4.5 - 1.) * 3.5 / 2.8 + u - v) * w / u;
}

Here is the call graph for this function:

template<typename T >

void func_tmp	(	const T &	u,
		const T &	v,
		const T &	w,
		doublereal	e,
		T &	f
	)

Definition at line 314 of file gradienttest.cc.

References bCompare(), grad::cos(), GRADIENT_ASSERT, grad::pow(), grad::sin(), and grad::tan().

                                                                      {
        f = u;
        f *= 3;
        f += 2 * v;
        f *= (u - v);
        f /= (1 - w);
        f *= pow(u/w, v - 1);
        f *= sin(v);
        f *= cos(w);
        f *= (1 - tan(-w + v));
        f *= e;

#if GRADIENT_DEBUG > 0
        const T tmp1 =
#endif
        f++;

        GRADIENT_ASSERT(tmp1 == f - 1);

        f -= 11.;
        f += 4.5;

#if GRADIENT_DEBUG > 0
        const T tmp2 =
#endif
        --f;

        GRADIENT_ASSERT(tmp2 == f);

#if GRADIENT_DEBUG > 0
        const T tmp3 = f;
        f += f;
        GRADIENT_ASSERT(bCompare(f, 2 * tmp3));
        f -= f;
        GRADIENT_ASSERT(bCompare(f, T()));
        f = tmp3;
#endif

        f *= 3.5;
        f /= 2.8;
        f += u;
        f -= v;
        f *= w;
        f /= u;
}

Here is the call graph for this function:

void funcDeriv	(	index_type	N,
		doublereal	u,
		const doublereal	ud[],
		doublereal	v,
		const doublereal	vd[],
		doublereal	w,
		const doublereal	wd[],
		doublereal	e,
		doublereal &	f,
		doublereal	fd[]
	)

Definition at line 364 of file gradienttest.cc.

References grad::cos(), grad::log(), grad::pow(), sec(), grad::sin(), and grad::tan().

Referenced by testGradient(), and testGradient2().

                                                                                                                                                                                           {
        const doublereal tmp1 = sin(v);
        const doublereal tmp2 = cos(w);
        const doublereal tmp3 = u/w;
        const doublereal tmp6 = pow(tmp3, v-1);
        const doublereal tmp4 = tmp1*tmp6*tmp2;
        const doublereal tmp5 = (tan(w-v)+1);
        const doublereal tmp7 = ((2*v+3*u)*tmp4*tmp5);
        const doublereal tmp8 = pow(sec(w-v),2);
        const doublereal tmp9 = (u-v);
        const doublereal tmp10 = (1-w);
        const doublereal tmp11 = (v-1);
        const doublereal tmp12 = (2*v+3*u);
        const doublereal tmp13 = tmp9*tmp4*tmp5;
        const doublereal tmp14 = tmp12*tmp4*tmp5;
        const doublereal tmp15 = (tmp9*tmp12*tmp4*tmp8);
        const doublereal tmp16 = (tmp9*tmp11*tmp14);
        const doublereal tmp17 = tmp15/tmp10;
        const doublereal tmp18 = tmp6*tmp2*tmp5;
        const doublereal tmp19 = tmp9*tmp12*tmp1;
        const doublereal tmp20 = tmp7/tmp10;
        const doublereal tmp21 = 3.5 / 2.8;
        const doublereal tmp22 = e * tmp21;
        const doublereal tmp23 = ((((tmp9*tmp14)/tmp10) * e + 1. - 11. + 4.5 - 1.) * tmp21 + u - v);

    f = tmp23 * w / u;
    
    doublereal df_du = ((tmp16/(u*tmp10)+tmp20+(3*tmp13)/tmp10)*tmp22 + 1.) * w / u - tmp23 * w / (u * u);
    
    doublereal df_dv = (((tmp19*log(tmp3)*tmp18)/tmp10-tmp20+(2*tmp13)/tmp10+(tmp9*tmp12*cos(v)*tmp18)/tmp10-tmp17)*tmp22 - 1.) * w / u;
    
    doublereal df_dw = (-(tmp19*tmp6*sin(w)*tmp5)/tmp10-tmp16/(tmp10*w)+(tmp9*tmp14)/pow(1-w,2)+tmp17)*tmp22 * w / u + tmp23 / u;

    for (index_type i = 0; i < N; ++i) {
        fd[i] = df_du * ud[i] + df_dv * vd[i] + df_dw * wd[i];
    }
}

Here is the call graph for this function:

int main	(	int	argc,
		char *	argv[]
	)

Definition at line 1817 of file gradienttest.cc.

References GRADIENT_DEBUG, NLoops, and testGradientCopy().

                                 {
#ifdef HAVE_FEENABLEEXCEPT
        feenableexcept(FE_INVALID|FE_DIVBYZERO|FE_OVERFLOW);
#endif

        if (argc > 1) {
                NLoops = atol(argv[1]);
        }

        if (argc > 2) {
                NLoopsAss = atol(argv[2]);
        }

    testRangeVector<0>();
    testRangeVector<6>();
    testMapVector<0>();
    testMapVector<19>();

    for (index_type i = 1; i <= 12; ++i) {
        testGradient<0>(i);
    }

    testGradient<1>(1);
    testGradient<2>(2);
    testGradient<4>(4);
    testGradient<6>(6);
    testGradient<8>(8);
    testGradient<10>(10);
    testGradient<12>(12);
    testGradient2<0>();
    testGradient2<3>();

    testGradientCopy(0);
    testGradientCopy(1);
    testGradientCopy(2);
    testGradientCopy(4);
    testGradientCopy(6);
    testGradientCopy(8);
    testGradientCopy(10);
    testGradientCopy(15);
    testGradientLin<4>(4);
    testGradientLin<0>(4);

#ifdef HAVE_BLITZ
    gradAssTest::testAssembly();
#endif

    testDifferentDofMaps<3>(1);
    testDifferentDofMaps<4>(2);
    testDifferentDofMaps<5>(3);
    testDifferentDofMaps<6>(4);
    testDifferentDofMaps<7>(5);
    testDifferentDofMaps<8>(6);

    for (integer i = 1; i <= 10; ++i) {
        testDifferentDofMaps<0>(i);
        testDifferentDofMaps<12>(i);
        testDifferentDofMaps<14>(i);
        testDifferentDofMaps<16>(i);
    }

#ifndef NDEBUG
    std::cerr << "All tests passed" << std::endl;
#else
    std::cerr << "No tests have been done" << std::endl;
#endif

    std::cerr << "GRADIENT_DEBUG=" << GRADIENT_DEBUG << std::endl;

    return 0;
}

Here is the call graph for this function:

doublereal random1

(

)

Definition at line 86 of file gradienttest.cc.

Referenced by testGradient(), and testGradient2().

                     {
        return 2 * (doublereal(rand()) / RAND_MAX) - 1;
}

doublereal sec

(

doublereal

x

)

Definition at line 360 of file gradienttest.cc.

References grad::cos().

Referenced by funcDeriv().

                             {
    return 1./cos(x);
}

Here is the call graph for this function:

template<index_type N_SIZE>

void testDifferentDofMaps

(

index_type

N

)

Definition at line 730 of file gradienttest.cc.

References grad::Gradient< N_SIZE >::dGetDerivativeGlobal(), grad::Gradient< N_SIZE >::dGetValue(), and grad::MapVectorBase::GLOBAL.

                                        {
        {
                LocalDofMap dofMap1, dofMap2;

                Gradient<N_SIZE> g1(1000., MapVector<N_SIZE>(&dofMap1, 100, 100 + N, MapVectorBase::GLOBAL, 1.));
                Gradient<N_SIZE> g2(2000., MapVector<N_SIZE>(&dofMap2, 102, 102 + N, MapVectorBase::GLOBAL, 2.));

                std::cerr << "testDifferentDofMaps<" << N_SIZE << ">(" << N << "): operator+=()\n";
                std::cerr << "g1=" << g1 << std::endl;
                std::cerr << "g2=" << g2 << std::endl;

                g1 += g2;

                std::cerr << "g1=" << g1 << std::endl;

                assert(g1.dGetValue() == 3000.);

                for (index_type i = 100; i < 102 + N; ++i) {
                        if (i < 100 + N || i >= 102) {
                                const doublereal d = g1.dGetDerivativeGlobal(i);
                                if (i >= 102 && i < 100 + N) {
                                        assert(d == 3.);
                                } else if (i < 100 + N) {
                                        assert(d == 1.);
                                } else {
                                        assert(d == 2.);
                                }
                        }
                }
        }
        {
                LocalDofMap dofMap1, dofMap2;

                const doublereal u = 3, v = 4, du = 5., dv = 6.;
                Gradient<N_SIZE> g1(u, MapVector<N_SIZE>(&dofMap1, 100, 100 + N, MapVectorBase::GLOBAL, du));
                Gradient<N_SIZE> g2(v, MapVector<N_SIZE>(&dofMap2, 102, 102 + N, MapVectorBase::GLOBAL, dv));

                std::cerr << "testDifferentDofMaps<" << N_SIZE << ">(" << N << "): operator*=()\n";
                std::cerr << "g1=" << g1 << std::endl;
                std::cerr << "g2=" << g2 << std::endl;

                g1 *= g2;

                std::cerr << "g1=" << g1 << std::endl;

                assert(g1.dGetValue() == u * v);

                for (index_type i = 100; i < 102 + N; ++i) {
                        if (i < 100 + N || i >= 102) {
                                const doublereal d = g1.dGetDerivativeGlobal(i);
                                if (i >= 102 && i < 100 + N) {
                                        assert(d == du * v + u * dv);
                                } else if (i < 100 + N) {
                                        assert(d == du * v);
                                } else {
                                        assert(d == u * dv);
                                }
                        }
                }
        }
}

Here is the call graph for this function:

template<int N_SIZE>

void testGradient

(

index_type

N

)

Definition at line 403 of file gradienttest.cc.

References bCompare(), grad::Gradient< N_SIZE >::DerivativeResizeReset(), grad::Gradient< N_SIZE >::dGetDerivativeGlobal(), grad::Gradient< N_SIZE >::dGetValue(), dof, grad::fabs(), func(), funcDeriv(), grad::MapVectorBase::GLOBAL, NLoops, random1(), grad::Gradient< N_SIZE >::SetDerivativeGlobal(), grad::Gradient< N_SIZE >::SetValuePreserve(), grad::sqrt(), and tic().

Referenced by testMatVecProductGradient().

                                {
    assert(N_SIZE == 0 || N_SIZE == N);

    srand(0);

    LocalDofMap dof;
    Gradient<N_SIZE> u, v, w, f, f_tmp;
    doublereal f_d;
    
    u.SetValuePreserve(fabs(random1()) * 10);
    v.SetValuePreserve(random1() * 3);
    w.SetValuePreserve(random1() * 0.01);
    
    u.DerivativeResizeReset(&dof, 0, N, MapVectorBase::GLOBAL, 0.);
    v.DerivativeResizeReset(&dof, 0, N, MapVectorBase::GLOBAL, 0.);
    w.DerivativeResizeReset(&dof, 0, N, MapVectorBase::GLOBAL, 0.);
    
    for (index_type i = 0; i < N; ++i) {
        u.SetDerivativeGlobal(i, random1() * 0.01);
        v.SetDerivativeGlobal(i, random1() * 3);
        w.SetDerivativeGlobal(i, random1() * 10);
    }
    
    doublereal* ud_C = new doublereal[N];
    doublereal* vd_C = new doublereal[N];
    doublereal* wd_C = new doublereal[N];
    doublereal f_C;
    doublereal* fd_C = new doublereal[N];

    for (index_type i = 0; i < N; ++i) {
        ud_C[i] = u.dGetDerivativeGlobal(i);
        vd_C[i] = v.dGetDerivativeGlobal(i);
        wd_C[i] = w.dGetDerivativeGlobal(i);
    }

    doublereal t_AD = 0., t_AD_tmp = 0., t_C = 0., t_d = 0.;

    for (int i = 0; i < NLoops; ++i) {
        const doublereal e = random1();
        doublereal start, stop;

        tic(start);
        func(u, v, w, e, f);
        tic(stop);

        t_AD += stop - start;

        tic(start);
        func(u.dGetValue(), v.dGetValue(), w.dGetValue(), e, f_d);
        tic(stop);

        t_d += stop - start;

        tic(start);
        func(u, v, w, e, f_tmp);
        tic(stop);

        t_AD_tmp += stop - start;

                tic(start);
                funcDeriv(N,
                                  u.dGetValue(),
                                  ud_C,
                                  v.dGetValue(),
                                  vd_C,
                                  w.dGetValue(),
                                  wd_C,
                                  e,
                                  f_C,
                                  fd_C);
                tic(stop);

                t_C += stop - start;

                assert(bCompare(f, f_tmp, std::numeric_limits<doublereal>::epsilon()));
            assert(bCompare(f.dGetValue(), f_C, sqrt(std::numeric_limits<doublereal>::epsilon())));
            assert(bCompare(f.dGetValue(), f_d, std::numeric_limits<doublereal>::epsilon()));

            for (index_type j = 0; j < N; ++j) {
                doublereal err = std::abs(f.dGetDerivativeGlobal(j) / fd_C[j] - 1.);
                assert(err < sqrt(std::numeric_limits<scalar_deriv_type>::epsilon()));
            }
    }
    
    std::cerr << "testGradient<" << N_SIZE << ">(" << N << ")" << std::endl;

    std::cout << "dof map:" << std::endl << dof << std::endl;
    std::cout << "u=" << u << std::endl;
    std::cout << "v=" << v << std::endl;
    std::cout << "w=" << w << std::endl;
    std::cout << "f=" << f << std::endl;
    std::cout << "f_tmp=" << f_tmp << std::endl;

    std::cout << "f_C=" << f_C << std::endl;
    std::cout << "fd_C=" << std::endl;

    for (int i = 0; i < N; ++i) {
        std::cout << " " << fd_C[i] << std::endl;
    }

    std::cout << std::endl;

    std::cerr << "t_AD=" << t_AD << "s" << std::endl;
    std::cerr << "t_d=" << t_d << "s" << std::endl;
    std::cerr << "t_AD_tmp=" << t_AD_tmp << "s" << std::endl;
    std::cerr << "t_C=" << t_C << "s" << std::endl;
    std::cerr << "overhead:" << t_AD / t_C << std::endl;
    std::cerr << "overhead tmp:" << t_AD_tmp / t_C << std::endl;

    delete [] ud_C;
    delete [] vd_C;
    delete [] wd_C;
    delete [] fd_C;
}

Here is the call graph for this function:

template<index_type N_SIZE>

void testGradient2

(

)

Definition at line 519 of file gradienttest.cc.

References grad::Alias(), bCompare(), grad::Gradient< N_SIZE >::DerivativeResizeReset(), grad::Gradient< N_SIZE >::dGetDerivativeGlobal(), grad::Gradient< N_SIZE >::dGetValue(), dof, grad::fabs(), func(), funcDeriv(), grad::MapVectorBase::GLOBAL, GRADIENT_ASSERT, NLoops, random1(), grad::Gradient< N_SIZE >::SetValuePreserve(), grad::sqrt(), and tic().

                     {
        assert(N_SIZE == 0 || N_SIZE == 3);

        const doublereal dTol = sqrt(std::numeric_limits<doublereal>::epsilon());

    srand(0);

    LocalDofMap dof;
    Gradient<N_SIZE> u, v, w, f;

    u.SetValuePreserve(fabs(random1()) * 10);
    v.SetValuePreserve(fabs(random1()) * 3);
    w.SetValuePreserve(fabs(random1()) * 0.01);

    u.DerivativeResizeReset(&dof, 0, MapVectorBase::GLOBAL, random1() * 0.01);
    v.DerivativeResizeReset(&dof, 1, MapVectorBase::GLOBAL, random1() * 3);
    w.DerivativeResizeReset(&dof, 2, MapVectorBase::GLOBAL, random1() * 10);

    doublereal ud_C[3];
    doublereal vd_C[3];
    doublereal wd_C[3];
    doublereal f_C;
    doublereal fd_C[3];

    for (index_type i = 0; i < 3; ++i) {
        ud_C[i] = u.dGetDerivativeGlobal(i);
        vd_C[i] = v.dGetDerivativeGlobal(i);
        wd_C[i] = w.dGetDerivativeGlobal(i);
    }

    doublereal t_AD = 0., t_C = 0.;

    for (int i = 0; i < NLoops; ++i) {
        const doublereal e = random1();
        doublereal start, stop;

        tic(start);
        func(u, v, w, e, f);
        tic(stop);

        t_AD += stop - start;

                tic(start);
                funcDeriv(3,
                                  u.dGetValue(),
                                  ud_C,
                                  v.dGetValue(),
                                  vd_C,
                                  w.dGetValue(),
                                  wd_C,
                                  e,
                                  f_C,
                                  fd_C);
                tic(stop);

                t_C += stop - start;

            assert(bCompare(f.dGetValue(), f_C, dTol));
            //assert(bCompare(Y[0], f_C, dTol));

            for (index_type j = 0; j < 3; ++j) {
                assert(bCompare<scalar_deriv_type>(f.dGetDerivativeGlobal(j), fd_C[j], sqrt(std::numeric_limits<scalar_deriv_type>::epsilon())));
            //  assert(bCompare(jac[j], fd_C[j], dTol));
            }
    }

    std::cerr << "testGradient2<" << N_SIZE << ">(3)" << std::endl;

    std::cout << "dof map:" << std::endl << dof << std::endl;
    std::cout << "u=" << u << std::endl;
    std::cout << "v=" << v << std::endl;
    std::cout << "w=" << w << std::endl;
    std::cout << "f=" << f << std::endl;

    std::cout << "f_C=" << f_C << std::endl;
    std::cout << "fd_C=" << std::endl;

    for (int i = 0; i < 3; ++i) {
        std::cout << " " << fd_C[i] << std::endl;
    }

    std::cout << std::endl;

    std::cerr << "t_AD=" << t_AD << "s" << std::endl;
    std::cerr << "t_C=" << t_C << "s" << std::endl;
    std::cerr << "overhead:" << t_AD / t_C << std::endl;

    Gradient<N_SIZE> w2 = 3 * u + 2 * v + 0.5 * w;

    w = 3 * u + 2 * v + 0.5 * Alias(w);

    GRADIENT_ASSERT(bCompare(w, w2));

    std::cout << "w=" << w << std::endl;
    std::cout << "w2=" << w2 << std::endl;
}

Here is the call graph for this function:

void testGradientCopy

(

int

N

)

Definition at line 646 of file gradienttest.cc.

References grad::MapVectorBase::GLOBAL, NLoops, tic(), and toc().

Referenced by main().

                             {
        LocalDofMap map1, map2;
        bool bFirst = true;

        std::cerr << "testGradientCopy(" << N << ")\n";

        tic();

        for (int i = 0; i < NLoops; ++i) {
                MapVector<0> ud(&map1, 100, 102, MapVectorBase::GLOBAL, 1.);

                for (int i = 0; i < N; ++i) {
                        MapVector<0> dummy(&map1, 102 + i, MapVectorBase::GLOBAL, 1.5);
                }

                MapVector<0> vd(&map1, 200, 203, MapVectorBase::GLOBAL, 2.);

                for (int i = 0; i < N; ++i) {
                        MapVector<0> dummy(&map1, 203 + i, MapVectorBase::GLOBAL, 2.5);
                }

                MapVector<0> wd(&map1, 300, 302, MapVectorBase::GLOBAL, 3.);

                for (int i = 0; i < N; ++i) {
                        MapVector<0> dummy(&map1, 302 + i, MapVectorBase::GLOBAL, 3.5);
                }

                Gradient<0> u(100., ud);
                Gradient<0> v(200., vd);
                Gradient<0> w(300., wd);

                assert(u != 100.001);
                assert(u >= 100.);
                assert(u <= 100.);
                assert(u == 100.);
                assert(u > 99.999);
                assert(u < 100.001);
                assert(100. == u);
                assert(100.001 != u);
                assert(100. <= u);
                assert(100. >= u);
                assert(99.999 < u);
                assert(100.001 > u);

                assert(v != 200.001);
                assert(v >= 200.);
                assert(v <= 200.);
                assert(v == 200.);
                assert(v > 199.999);
                assert(v < 200.001);
                assert(200. == v);
                assert(200.001 != v);
                assert(200. <= v);
                assert(200. >= v);
                assert(199.999 < v);
                assert(200.001 > v);

                Gradient<0> x = v + w;

                if (bFirst) {
                        std::cout << std::endl;
                        std::cout << "map1=" << std::endl << map1 << std::endl;
                        std::cout << "u=" << u << std::endl;
                        std::cout << "v=" << v << std::endl;
                        std::cout << "w=" << w << std::endl;
                        std::cout << "x=" << x << std::endl;
                        std::cout << std::endl;
                }

                Gradient<5> y(x, &map2);

                if (bFirst) {
                        std::cout << "map2=" << std::endl << map2 << std::endl;
                        std::cout << "y=" << y << std::endl;
                        std::cout << std::endl;
                }

                bFirst = false;
        }

        std::cerr << "t=" << toc() << "s" << std::endl;
}

Here is the call graph for this function:

template<index_type N_SIZE>

void testGradientLin

(

int

N

)

Definition at line 793 of file gradienttest.cc.

References grad::Gradient< N_SIZE >::bIsEqual(), grad::Gradient< N_SIZE >::DerivativeResizeReset(), grad::Gradient< N_SIZE >::dGetDerivativeGlobal(), grad::Gradient< N_SIZE >::dGetValue(), grad::MapVectorBase::GLOBAL, and grad::Gradient< N_SIZE >::SetValuePreserve().

                            {
        Gradient<N_SIZE> A11, A12, A13, A14;
        Gradient<N_SIZE> x1, x2, x3, x4;
        LocalDofMap oDof;

        A11.SetValuePreserve(1);
        A12.SetValuePreserve(2);
        A13.SetValuePreserve(3);
        A14.SetValuePreserve(4);

        A11.DerivativeResizeReset(&oDof, 0, MapVectorBase::GLOBAL, 1);
        A12.DerivativeResizeReset(&oDof, 1, MapVectorBase::GLOBAL, 1);
        A13.DerivativeResizeReset(&oDof, 2, MapVectorBase::GLOBAL, 1);
        A14.DerivativeResizeReset(&oDof, 3, MapVectorBase::GLOBAL, 1);

        std::cout << "A11=" << A11 << std::endl;
        std::cout << "A12=" << A12 << std::endl;
        std::cout << "A13=" << A13 << std::endl;
        std::cout << "A14=" << A14 << std::endl;

        x1.SetValuePreserve(1);
        x2.SetValuePreserve(1);
        x3.SetValuePreserve(1);
        x4.SetValuePreserve(1);

        x1.DerivativeResizeReset(&oDof, 0, MapVectorBase::GLOBAL, 1);
        x2.DerivativeResizeReset(&oDof, 1, MapVectorBase::GLOBAL, 1);
        x3.DerivativeResizeReset(&oDof, 2, MapVectorBase::GLOBAL, 1);
        x4.DerivativeResizeReset(&oDof, 3, MapVectorBase::GLOBAL, 1);

        std::cout << "x1=" << x1 << std::endl;
        std::cout << "x2=" << x2 << std::endl;
        std::cout << "x3=" << x3 << std::endl;
        std::cout << "x4=" << x4 << std::endl;

        Gradient<N_SIZE> y1 = A11 * x1;
        Gradient<N_SIZE> y2 = A12 * x2;
        Gradient<N_SIZE> y3 = A13 * x3;
        Gradient<N_SIZE> y4 = A14 * x4;

        std::cout << "y1=" << y1 << std::endl;
        std::cout << "y2=" << y2 << std::endl;
        std::cout << "y3=" << y3 << std::endl;
        std::cout << "y4=" << y4 << std::endl;

        Gradient<N_SIZE> y_2 = y1 + y2;
        Gradient<N_SIZE> y_3 = y_2 + y3;
        Gradient<N_SIZE> y_4 = y_3 + y4;

        Gradient<N_SIZE> y_1 = y1 + y2 + y3 + y4;
        Gradient<N_SIZE> y = A11 * x1 + A12 * x2 + A13 * x3 + A14 * x4;

        std::cout << "y=" << y << std::endl;
        std::cout << "y_1=" << y_1 << std::endl;
        std::cout << "y_2=" << y_2 << std::endl;
        std::cout << "y_3=" << y_3 << std::endl;
        std::cout << "y_4=" << y_4 << std::endl;

        assert(y.bIsEqual(y_1));
        assert(y.bIsEqual(y_4));
        assert(y.dGetValue() == 10);
        assert(y.dGetDerivativeGlobal(0) == 2);
        assert(y.dGetDerivativeGlobal(1) == 3);
        assert(y.dGetDerivativeGlobal(2) == 4);
        assert(y.dGetDerivativeGlobal(3) == 5);
}

Here is the call graph for this function:

template<int N_SIZE>

void testMapVector

(

)

Definition at line 267 of file gradienttest.cc.

References grad::MapVector< N_SIZE >::dGetLocalVector(), dof, grad::MapVector< N_SIZE >::iGetEndIndexLocal(), grad::LocalDofMap::iGetGlobalDof(), grad::MapVector< N_SIZE >::iGetStartIndexLocal(), and grad::LocalDofMap::Size().

                     {
    assert(N_SIZE == 0 || N_SIZE == 19);
    LocalDofMap dof;
    
    const int iCount = 4;
    
    MapVector<N_SIZE> v[iCount] = {   MapVector<N_SIZE>(&dof, 501, 506+1, MapVector<N_SIZE>::GLOBAL, 1.),
                                      MapVector<N_SIZE>(&dof, 601, 606+1, MapVector<N_SIZE>::GLOBAL, 2.),
                                      MapVector<N_SIZE>(&dof, 701, 706+1, MapVector<N_SIZE>::GLOBAL, 3.),
                                      MapVector<N_SIZE>(&dof, 801, 801+1, MapVector<N_SIZE>::GLOBAL, 4.)   };

    MapVector<N_SIZE> v1(&dof, 0, 1+1, MapVector<N_SIZE>::LOCAL, 5.);
    MapVector<N_SIZE> v2(&dof, 1, 2+1, MapVector<N_SIZE>::LOCAL, 6.);                           

    for (index_type j = 0; j < dof.Size(); ++j) {
        for (int i = 0; i < iCount; ++i) {
            if (j >= v[i].iGetStartIndexLocal() && j < v[i].iGetEndIndexLocal()) {
                doublereal dVal = v[i].dGetLocalVector(j);

                std::cout << "v[" << i << "][" << j << "->" << dof.iGetGlobalDof(j) << "]=" << dVal << std::endl;
            }
        }
    }

    for (index_type j = 0; j < dof.Size(); ++j) {
        if (j >= v1.iGetStartIndexLocal() && j < v1.iGetEndIndexLocal()) {
            doublereal dVal = v1.dGetLocalVector(j);

            std::cout << "v1[" << j << "->" << dof.iGetGlobalDof(j) << "]=" << dVal << std::endl;
        }
    }

    for (index_type j = 0; j < dof.Size(); ++j) {
        if (j >= v2.iGetStartIndexLocal() && j < v2.iGetEndIndexLocal()) {
            doublereal dVal = v2.dGetLocalVector(j);

            std::cout << "v2[" << j << "->" << dof.iGetGlobalDof(j) << "]=" << dVal << std::endl;
        }
    }
}

Here is the call graph for this function:

template<int N_SIZE>

void testRangeVector

(

)

Definition at line 127 of file gradienttest.cc.

References grad::RangeVector< T, N_SIZE >::GetValue(), grad::RangeVector< T, N_SIZE >::GetVectorValue(), grad::RangeVector< T, N_SIZE >::iGetEndIndex(), grad::RangeVector< T, N_SIZE >::iGetEndIndexVector(), grad::RangeVector< T, N_SIZE >::iGetStartIndex(), grad::RangeVector< T, N_SIZE >::iGetStartIndexVector(), grad::RangeVector< T, N_SIZE >::Reserve(), grad::RangeVector< T, N_SIZE >::Reset(), grad::RangeVector< T, N_SIZE >::ResizePreserve(), grad::RangeVector< T, N_SIZE >::ResizeReset(), grad::RangeVector< T, N_SIZE >::SetValue(), and grad::RangeVector< T, N_SIZE >::SetVectorValue().

                       {
        std::cout << "testRangeVector<" << N_SIZE << ">();\n";
    assert(N_SIZE == 0 || N_SIZE >= 6);
    index_type start = 3;
    index_type end = 6;
    RangeVector<doublereal, N_SIZE> v(start, end, 0.);
    RangeVector<float, 2 * N_SIZE> vf(v), vf2;
    RangeVector<doublereal, N_SIZE> v2(v);

    int k = 0;
    
    for (index_type i = start; i < end; ++i) {
        v.SetValue(i, ++k);
    }
    
    v2 = v;
    vf = v;
    vf.Reset();
    vf.Reserve(N_SIZE == 0 ? vf.iGetEndIndex() * 2 : 2 * N_SIZE);
    vf = v2;
    const RangeVector<float, N_SIZE> v3(v2);
    v2.Reset();
    v2 = v3;
    v2.Reserve(N_SIZE == 0 ? v2.iGetEndIndex() * 2 : N_SIZE);

    std::cout << "v(" << v.iGetStartIndex() << ":" << v.iGetEndIndex() - 1 << ")\n";

    for (index_type i = v.iGetStartIndex(); i < v.iGetEndIndex(); ++i) {
        std::cout << "v[" << i << "]=" << v.GetValue(i) << std::endl;
    }

    for (index_type i = vf.iGetStartIndex(); i < vf.iGetEndIndex(); ++i) {
        std::cout << "vf[" << i << "]=" << vf.GetValue(i) << std::endl;
    }

    for (index_type i = v2.iGetStartIndex(); i < v2.iGetEndIndex(); ++i) {
        std::cout << "v2[" << i << "]=" << v2.GetValue(i) << std::endl;
    }

    for (index_type i = v3.iGetStartIndex(); i < v3.iGetEndIndex(); ++i) {
        std::cout << "v3[" << i << "]=" << v3.GetValue(i) << std::endl;
    }

    std::cout << std::endl;

    for (index_type i = 0; i < v.iGetEndIndex(); ++i)
    {
        assert(v.GetValue(i) == vf.GetValue(i));
        assert(v2.GetValue(i) == v.GetValue(i));
    }

    v2.ResizeReset(v.iGetStartIndex(), v.iGetEndIndex(), 0.);

    for (index_type i = v.iGetStartIndexVector(); i < v.iGetEndIndexVector(); ++i)
    {
        v2.SetVectorValue(i, v.GetVectorValue(i));
    }

    for (index_type i = v.iGetStartIndex(); i < v.iGetEndIndex(); ++i)
    {
        assert(v2.GetValue(i) == v.GetValue(i));
    }

    vf.ResizeReset(0, 12, 0.0f);

    for (index_type i = vf.iGetStartIndex(); i < vf.iGetEndIndex(); ++i)
    {
        vf.SetValue(i, float(i + 1));
    }

    vf2.ResizeReset(vf.iGetStartIndex(), vf.iGetEndIndex(), 0.0f);

    for (index_type i = vf.iGetStartIndexVector(); i < vf.iGetEndIndexVector(); ++i)
    {
        vf2.SetVectorValue(i, vf.GetVectorValue(i));
    }

    for (index_type i = vf.iGetStartIndex(); i < vf.iGetEndIndex(); ++i)
    {
        assert(vf.GetValue(i) == vf2.GetValue(i));
    }

    for (index_type i = vf2.iGetStartIndex(); i < vf2.iGetEndIndex(); ++i) {
        std::cout << "vf2[" << i << "]=" << vf2.GetValue(i) << std::endl;
    }

    std::cout << std::endl;

    v.ResizePreserve(1, 6);

    std::cout << "v(" << v.iGetStartIndex() << ":" << v.iGetEndIndex() - 1 << ")\n";

    for (index_type i = v.iGetStartIndex(); i < v.iGetEndIndex(); ++i) {
        std::cout << "v[" << i << "]=" << v.GetValue(i) << std::endl;
    }

    std::cout << std::endl;

    v.ResizePreserve(1, 5);

    std::cout << "v(" << v.iGetStartIndex() << ":" << v.iGetEndIndex() - 1 << ")\n";

    for (index_type i = v.iGetStartIndex(); i < v.iGetEndIndex(); ++i) {
        std::cout << "v[" << i << "]=" << v.GetValue(i) << std::endl;
    }

    std::cout << std::endl;

    v.ResizePreserve(0, 6);
    
    std::cout << "v(" << v.iGetStartIndex() << ":" << v.iGetEndIndex() - 1 << ")\n";

    for (index_type i = v.iGetStartIndex(); i < v.iGetEndIndex(); ++i) {
        std::cout << "v[" << i << "]=" << v.GetValue(i) << std::endl;
    }

    std::cout << std::endl;

    v.ResizePreserve(4, 5);

    std::cout << "v(" << v.iGetStartIndex() << ":" << v.iGetEndIndex() - 1 << ")\n";

    for (index_type i = v.iGetStartIndex(); i < v.iGetEndIndex(); ++i) {
        std::cout << "v[" << i << "]=" << v.GetValue(i) << std::endl;
    }

    std::cout << std::endl;

    v.ResizePreserve(0, 6);

    std::cout << "v(" << v.iGetStartIndex() << ":" << v.iGetEndIndex() - 1 << ")\n";

    for (index_type i = v.iGetStartIndex(); i < v.iGetEndIndex(); ++i) {
        std::cout << "v[" << i << "]=" << v.GetValue(i) << std::endl;
    }

    std::cout << std::endl;
}

Here is the call graph for this function:

void tic

(

)

Definition at line 866 of file gradienttest.cc.

References dStartTime, and mbdyn_clock_time().

Referenced by testGradient(), testGradient2(), and testGradientCopy().

           {
    dStartTime = mbdyn_clock_time();
}

Here is the call graph for this function:

void tic

(

doublereal &

dTime

)

Definition at line 862 of file gradienttest.cc.

References mbdyn_clock_time().

                            {
        dTime = mbdyn_clock_time();
}

Here is the call graph for this function:

doublereal toc

(

)

Definition at line 870 of file gradienttest.cc.

References dStartTime, and mbdyn_clock_time().

Referenced by testGradientCopy().

                 {
        return mbdyn_clock_time() - dStartTime;
}

Here is the call graph for this function:

Variable Documentation

doublereal dStartTime

Definition at line 860 of file gradienttest.cc.

Referenced by tic(), and toc().

int NLoops = 1

Definition at line 79 of file gradienttest.cc.

Referenced by main(), testGradient(), testGradient2(), and testGradientCopy().

int NLoopsAss = 1

Definition at line 80 of file gradienttest.cc.