source: trunk/ADOL-C/src/adouble.h @ 42

Last change on this file since 42 was 42, checked in by awalther, 11 years ago

set svn keywords property

  • Property svn:keywords set to Author Date Id Revision
File size: 42.8 KB
Line 
1/* ---------------------------------------------------------------------------
2 ADOL-C -- Automatic Differentiation by Overloading in C++
3
4 Revision: $Id: adouble.h 42 2009-07-15 18:37:17Z awalther $
5 Contents: adouble.h contains the basis for the class of adouble
6           included here are all the possible functions defined on
7           the adouble class.  Notice that, as opposed to ealier versions,
8           both the class adub and the class adouble are derived from a base
9           class (badouble).  See below for further explanation.
10
11 Copyright (c) Andrea Walther, Andreas Griewank, Andreas Kowarz,
12               Hristo Mitev, Sebastian Schlenkrich, Jean Utke, Olaf Vogel
13 
14 This file is part of ADOL-C. This software is provided as open source.
15 Any use, reproduction, or distribution of the software constitutes
16 recipient's acceptance of the terms of the accompanying license file.
17 
18---------------------------------------------------------------------------*/
19
20#if !defined(ADOLC_ADOUBLE_H)
21#define ADOLC_ADOUBLE_H 1
22
23/****************************************************************************/
24/*                                                         THIS FILE IS C++ */
25#ifdef __cplusplus
26
27#include <cstdio>
28#include <cstdlib>
29#include <iostream>
30#include <cmath>
31using std::cout;
32using std::cin;
33using std::cerr;
34using std::ostream;
35using std::istream;
36
37#include <common.h>
38
39/* NOTICE: There are automatic includes at the end of this file! */
40
41#undef TAPELESS
42#undef SAFE
43#if defined(ADOLC_TAPELESS)
44#  define TAPELESS
45#  undef SAFE
46#endif
47
48#if defined(SAFE_ADOLC_TAPELESS)
49#  define TAPELESS
50#  define SAFE
51#endif
52
53#if !defined(TAPELESS)
54
55/****************************************************************************/
56/*                                             FORWARD DECLARATIONS (TAPES) */
57
58/*--------------------------------------------------------------------------*/
59class adouble;
60class adub;
61class badouble;
62class adubv;
63/* class doublev;  that's history */
64
65/*--------------------------------------------------------------------------*/
66void ADOLC_DLL_EXPORT condassign( double &res, const double &cond,
67                                  const double &arg1, const double &arg2 );
68void ADOLC_DLL_EXPORT condassign( double &res, const double &cond,
69                                  const double &arg );
70
71double ADOLC_DLL_EXPORT fmin( const double &x, const double &y );
72double ADOLC_DLL_EXPORT fmax( const double &x, const double &y );
73
74
75/****************************************************************************/
76/*                                                           CLASS BADOUBLE */
77
78/**
79   The class badouble contains the basic definitions for
80   the arithmetic operations, comparisons, etc.
81   This is a basic class from which the adub and adouble are
82   derived.  Notice that the constructors/destructors for
83   the class badouble are of the trivial variety.  This is the
84   main difference among badoubles, adubs, and adoubles.
85*/
86class ADOLC_DLL_EXPORT badouble {
87    friend ADOLC_DLL_EXPORT class badoublev;
88protected:
89    locint location;
90    badouble( void ) {};
91    // must be public when using gcc >= 3.4 ( problems with value() )
92    // (see GCC 3.4 Release Series - Changes, New Features, and Fixes)
93    //
94    // badouble( const badouble& a ) {location = a.location;};
95    badouble( locint lo ) {
96        location = lo;
97    };
98
99public:
100    /*--------------------------------------------------------------------------*/
101    badouble( const badouble& a ) {
102        location = a.location;
103    }
104    ;           /* ctor */
105
106    inline locint loc( void ) const;                         /* Helpful stuff */
107
108    /*------------------------------------------------------------------------*/
109    badouble& operator >>= ( double& );                        /* Assignments */
110    badouble& operator <<= ( double );
111    void declareIndependent ();
112    void declareDependent ();
113    badouble& operator = ( double );
114    badouble& operator = ( const badouble& );
115    badouble& operator = ( const adub& );
116    double getValue() const;
117    inline double value() {
118        return getValue();
119    }
120    void setValue ( const double );
121    /* badouble& operator = ( const adouble& );
122       !!! olvo 991210: was the same as badouble-assignment */
123
124    /*--------------------------------------------------------------------------*/
125    friend ADOLC_DLL_EXPORT std::ostream& operator << ( std::ostream&, const badouble& );  /* IO friends */
126    friend ADOLC_DLL_EXPORT std::istream& operator >> ( std::istream&, const badouble& );
127
128    /*------------------------------------------------------------------------*/
129    badouble& operator += ( double );               /* Operation + Assignment */
130    badouble& operator += ( const badouble& );
131    badouble& operator -= ( double y );
132    badouble& operator -= ( const badouble& );
133    badouble& operator *= ( double );
134    badouble& operator *= ( const badouble& );
135    badouble& operator /= ( double );
136    badouble& operator /= ( const badouble& );
137    /* olvo 991122 n2l: new special op_codes */
138    badouble& operator += ( const adub& );
139    badouble& operator -= ( const adub& );
140
141    /*--------------------------------------------------------------------------*/
142    /* Comparison (friends) */
143    inline friend int operator != ( const badouble&, const badouble& );
144    inline friend int operator != ( double, const badouble& );
145    friend ADOLC_DLL_EXPORT int operator != ( const badouble&, double );
146    inline friend int operator == ( const badouble&, const badouble& );
147    inline friend int operator == ( double, const badouble& );
148    friend ADOLC_DLL_EXPORT int operator == ( const badouble&, double );
149    inline friend int operator <= ( const badouble&, const badouble& );
150    inline friend int operator <= ( double, const badouble& );
151    friend ADOLC_DLL_EXPORT int operator <= ( const badouble&, double );
152    inline friend int operator >= ( const badouble&, const badouble& );
153    inline friend int operator >= ( double, const badouble& );
154    friend ADOLC_DLL_EXPORT int operator >= ( const badouble&, double );
155    inline friend int operator >  ( const badouble&, const badouble& );
156    inline friend int operator >  ( double, const badouble& );
157    friend ADOLC_DLL_EXPORT int operator >  ( const badouble&, double );
158    inline friend int operator <  ( const badouble&, const badouble& );
159    inline friend int operator <  ( double, const badouble& );
160    friend ADOLC_DLL_EXPORT int operator <  ( const badouble&, double );
161
162
163    /*--------------------------------------------------------------------------*/
164    /* sign operators (friends) */
165    friend ADOLC_DLL_EXPORT adub operator + ( const badouble& x );
166    friend ADOLC_DLL_EXPORT adub operator - ( const badouble& x );
167
168    /*--------------------------------------------------------------------------*/
169    /* binary operators (friends) */
170    friend ADOLC_DLL_EXPORT adub operator + ( const badouble&, const badouble& );
171    friend ADOLC_DLL_EXPORT adub operator + ( double, const badouble& );
172    friend ADOLC_DLL_EXPORT adub operator + ( const badouble&, double );
173    friend ADOLC_DLL_EXPORT adub operator - ( const badouble&, const badouble& );
174    inline friend adub operator - ( const badouble&, double );
175    friend ADOLC_DLL_EXPORT adub operator - ( double, const badouble& );
176    friend ADOLC_DLL_EXPORT adub operator * ( const badouble&, const badouble& );
177    friend ADOLC_DLL_EXPORT adub operator * ( double, const badouble& );
178    inline friend adub operator * ( const badouble&, double );
179    inline friend adub operator / ( const badouble&, double );
180    friend ADOLC_DLL_EXPORT adub operator / ( const badouble&, const badouble& );
181    friend ADOLC_DLL_EXPORT adub operator / ( double, const badouble& );
182
183    /*--------------------------------------------------------------------------*/
184    /* unary operators (friends) */
185    friend ADOLC_DLL_EXPORT adub exp  ( const badouble& );
186    friend ADOLC_DLL_EXPORT adub log  ( const badouble& );
187    friend ADOLC_DLL_EXPORT adub sqrt ( const badouble& );
188    friend ADOLC_DLL_EXPORT adub sin  ( const badouble& );
189    friend ADOLC_DLL_EXPORT adub cos  ( const badouble& );
190    friend ADOLC_DLL_EXPORT adub tan  ( const badouble& );
191    friend ADOLC_DLL_EXPORT adub asin ( const badouble& );
192    friend ADOLC_DLL_EXPORT adub acos ( const badouble& );
193    friend ADOLC_DLL_EXPORT adub atan ( const badouble& );
194
195    /*--------------------------------------------------------------------------*/
196    /* special operators (friends) */
197    friend ADOLC_DLL_EXPORT adouble atan2 ( const badouble&, const badouble& );
198    /* no internal use of condassign: */
199    friend ADOLC_DLL_EXPORT adub    pow   ( const badouble&, double );
200    /* uses condassign internally */
201    friend ADOLC_DLL_EXPORT adouble pow   ( const badouble&, const badouble& );
202    friend ADOLC_DLL_EXPORT adouble pow   ( double, const badouble& );
203    friend ADOLC_DLL_EXPORT adub    log10 ( const badouble& );
204    /* User defined version of logarithm to test extend_quad macro */
205    friend ADOLC_DLL_EXPORT adouble myquad( const badouble& );
206
207    /*--------------------------------------------------------------------------*/
208    /* Additional ANSI C standard Math functions Added by DWJ on 8/6/90 */
209    friend ADOLC_DLL_EXPORT adub sinh  ( const badouble& );
210    friend ADOLC_DLL_EXPORT adub cosh  ( const badouble& );
211    friend ADOLC_DLL_EXPORT adub tanh  ( const badouble& );
212#if defined(ATRIG_ERF)
213    friend ADOLC_DLL_EXPORT adub asinh ( const badouble& );
214    friend ADOLC_DLL_EXPORT adub acosh ( const badouble& );
215    friend ADOLC_DLL_EXPORT adub atanh ( const badouble& );
216#endif
217
218    friend ADOLC_DLL_EXPORT adub fabs  ( const badouble& );
219    friend ADOLC_DLL_EXPORT adub ceil  ( const badouble& );
220    friend ADOLC_DLL_EXPORT adub floor ( const badouble& );
221
222    friend ADOLC_DLL_EXPORT adub fmax ( const badouble&, const badouble& );
223    friend ADOLC_DLL_EXPORT adub fmax ( double, const badouble& );
224    friend ADOLC_DLL_EXPORT adub fmax ( const badouble&, double );
225    friend ADOLC_DLL_EXPORT adub fmin ( const badouble&, const badouble& );
226    friend ADOLC_DLL_EXPORT adub fmin ( double, const badouble& );
227    friend ADOLC_DLL_EXPORT adub fmin ( const badouble&, double );
228
229    friend ADOLC_DLL_EXPORT adub ldexp ( const badouble&, int );
230    friend ADOLC_DLL_EXPORT adub frexp ( const badouble&, int* );
231    friend ADOLC_DLL_EXPORT adub erf   ( const badouble& );
232
233    /*--------------------------------------------------------------------------*/
234    /* Conditionals */
235    friend ADOLC_DLL_EXPORT void condassign( adouble &res, const adouble &cond,
236            const adouble &arg1, const adouble &arg2 );
237    friend ADOLC_DLL_EXPORT void condassign( adouble &res, const adouble &cond,
238            const adouble &arg );
239};
240
241
242
243/****************************************************************************/
244/*                                                               CLASS ADUB */
245
246/*
247   The class Adub
248   ---- Basically used as a temporary result.  The address for an
249        adub is usually generated within an operation.  That address
250        is "freed" when the adub goes out of scope (at destruction time).
251   ---- operates just like a badouble, but it has a destructor defined for it.
252*/
253
254class ADOLC_DLL_EXPORT adub:public badouble {
255    friend ADOLC_DLL_EXPORT class adouble;
256protected:
257    adub( locint lo ):badouble(lo) {};
258    adub( void ):badouble(0) {
259        fprintf(DIAG_OUT,"ADOL-C error: illegal default construction of adub"
260                " variable\n");
261        exit(-2);
262    };
263    adub( double ):badouble(0) {
264        fprintf(DIAG_OUT,"ADOL-C error: illegal  construction of adub variable"
265                " from double\n");
266        exit(-2);
267    };
268
269public:
270
271    /*--------------------------------------------------------------------------*/
272    /* sign operators (friends) */
273    friend ADOLC_DLL_EXPORT adub operator + ( const badouble& x );
274    friend ADOLC_DLL_EXPORT adub operator - ( const badouble& x );
275
276    /*--------------------------------------------------------------------------*/
277    /* binary operators (friends) */
278    friend ADOLC_DLL_EXPORT adub operator + ( const badouble&, const badouble& );
279    friend ADOLC_DLL_EXPORT adub operator + ( double, const badouble& );
280    friend ADOLC_DLL_EXPORT adub operator + ( const badouble&, double );
281    friend ADOLC_DLL_EXPORT adub operator - ( const badouble&, const badouble& );
282    inline friend adub operator - ( const badouble&, double );
283    friend ADOLC_DLL_EXPORT adub operator - ( double, const badouble& );
284    friend ADOLC_DLL_EXPORT adub operator * ( const badouble&, const badouble& );
285    friend ADOLC_DLL_EXPORT adub operator * ( double, const badouble& );
286    inline friend adub operator * ( const badouble&, double );
287    inline friend adub operator / ( const badouble&, double );
288    friend ADOLC_DLL_EXPORT adub operator / ( const badouble&, const badouble& );
289    friend ADOLC_DLL_EXPORT adub operator / ( double, const badouble& );
290
291    /*--------------------------------------------------------------------------*/
292    /* unary operators (friends) */
293    friend ADOLC_DLL_EXPORT adub exp  ( const badouble& );
294    friend ADOLC_DLL_EXPORT adub log  ( const badouble& );
295    friend ADOLC_DLL_EXPORT adub sqrt ( const badouble& );
296    friend ADOLC_DLL_EXPORT adub sin  ( const badouble& );
297    friend ADOLC_DLL_EXPORT adub cos  ( const badouble& );
298    friend ADOLC_DLL_EXPORT adub tan  ( const badouble& );
299    friend ADOLC_DLL_EXPORT adub asin ( const badouble& );
300    friend ADOLC_DLL_EXPORT adub acos ( const badouble& );
301    friend ADOLC_DLL_EXPORT adub atan ( const badouble& );
302
303    /*--------------------------------------------------------------------------*/
304    /* special operators (friends) */
305    /* no internal use of condassign: */
306    friend ADOLC_DLL_EXPORT adub    pow   ( const badouble&, double );
307    friend ADOLC_DLL_EXPORT adub    log10 ( const badouble& );
308
309    /*--------------------------------------------------------------------------*/
310    /* Additional ANSI C standard Math functions Added by DWJ on 8/6/90 */
311    friend ADOLC_DLL_EXPORT adub sinh  ( const badouble& );
312    friend ADOLC_DLL_EXPORT adub cosh  ( const badouble& );
313    friend ADOLC_DLL_EXPORT adub tanh  ( const badouble& );
314#if defined(ATRIG_ERF)
315    friend ADOLC_DLL_EXPORT adub asinh ( const badouble& );
316    friend ADOLC_DLL_EXPORT adub acosh ( const badouble& );
317    friend ADOLC_DLL_EXPORT adub atanh ( const badouble& );
318#endif
319
320    friend ADOLC_DLL_EXPORT adub fabs  ( const badouble& );
321    friend ADOLC_DLL_EXPORT adub ceil  ( const badouble& );
322    friend ADOLC_DLL_EXPORT adub floor ( const badouble& );
323
324    friend ADOLC_DLL_EXPORT adub fmax ( const badouble&, const badouble& );
325    friend ADOLC_DLL_EXPORT adub fmax ( double, const badouble& );
326    friend ADOLC_DLL_EXPORT adub fmax ( const badouble&, double );
327    friend ADOLC_DLL_EXPORT adub fmin ( const badouble&, const badouble& );
328    friend ADOLC_DLL_EXPORT adub fmin ( double, const badouble& );
329    friend ADOLC_DLL_EXPORT adub fmin ( const badouble&, double );
330
331    friend ADOLC_DLL_EXPORT adub ldexp ( const badouble&, int );
332    friend ADOLC_DLL_EXPORT adub frexp ( const badouble&, int* );
333    friend ADOLC_DLL_EXPORT adub erf   ( const badouble& );
334
335    ~adub();
336};
337
338
339/****************************************************************************/
340/*                                                            CLASS ADOUBLE */
341/*
342  The class adouble.
343  ---Derived from badouble.  Contains the standard constructors/destructors.
344  ---At construction, it is given a new address, and at destruction, that
345     address is freed.
346*/
347class ADOLC_DLL_EXPORT adouble:public badouble {
348public:
349    adouble( const adub& );
350    adouble( const adouble& );
351    adouble( void );
352    adouble( double );
353    /* adub prevents postfix operators to occur on the left
354       side of an assignment which would not work  */
355    adub operator++( int );
356    adub operator--( int );
357    badouble& operator++( void );
358    badouble& operator--( void );
359    /*   inline double value(); */
360    ~adouble();
361
362    adouble& operator = ( double );
363    adouble& operator = ( const badouble& );
364    /* adouble& operator = ( const adouble& );
365       !!! olvo 991210 was the same as badouble-assignment */
366    adouble& operator = ( const adub& );
367};
368
369
370/****************************************************************************/
371/*                                                       INLINE DEFINITIONS */
372
373/*--------------------------------------------------------------------------*/
374inline locint badouble::loc( void ) const {
375    return location;
376}
377
378/*--------------------------------------------------------------------------*/
379/* Comparison */
380inline int operator != ( const badouble& u, const badouble& v ) {
381    return (u-v != 0);
382}
383
384inline int operator != ( double coval, const badouble& v) {
385    if (coval)
386        return (-coval+v != 0);
387    else
388        return (v != 0);
389}
390
391inline int operator == ( const badouble& u, const badouble& v ) {
392    return (u-v == 0);
393}
394
395inline int operator == ( double coval, const badouble& v) {
396    if (coval)
397        return (-coval+v == 0);
398    else
399        return (v == 0);
400}
401
402inline int operator <= ( const badouble& u, const badouble& v ) {
403    return (u-v <= 0);
404}
405
406inline int operator <= ( double coval, const badouble& v ) {
407    if (coval)
408        return (-coval+v >= 0);
409    else
410        return (v >= 0);
411}
412
413inline int operator >= ( const badouble& u, const badouble& v ) {
414    return (u-v >= 0);
415}
416
417inline int operator >= ( double coval, const badouble& v ) {
418    if (coval)
419        return (-coval+v <= 0);
420    else
421        return (v <= 0);
422}
423
424inline int operator > ( const badouble& u, const badouble& v ) {
425    return (u-v > 0);
426}
427
428inline int operator > ( double coval, const badouble& v ) {
429    if (coval)
430        return (-coval+v < 0);
431    else
432        return (v < 0);
433}
434
435inline int operator < ( const badouble& u, const badouble& v ) {
436    return (u-v < 0);
437}
438
439inline int operator < ( double coval, const badouble& v ) {
440    if (coval)
441        return (-coval+v > 0);
442    else
443        return (v > 0);
444}
445
446/*--------------------------------------------------------------------------*/
447/* Subtract a floating point from an adouble  */
448inline adub operator - ( const badouble& x , double coval ) {
449    return (-coval) + x;
450}
451
452/*--------------------------------------------------------------------------*/
453/* Multiply an adouble by a floating point */
454inline adub operator * (const badouble& x, double coval) {
455    return coval * x;
456}
457
458/*--------------------------------------------------------------------------*/
459/* Divide an adouble by a floating point */
460inline adub operator / (const badouble& x, double coval) {
461    return (1.0/coval) * x;
462}
463
464/****************************************************************************/
465/* tapeless implementation                                                  */
466/****************************************************************************/
467#else
468
469namespace adtl {
470
471#if defined(NUMBER_DIRECTIONS)
472extern int ADOLC_numDir;
473#define ADOLC_TAPELESS_UNIQUE_INTERNALS int adtl::ADOLC_numDir = NUMBER_DIRECTIONS;
474#if !defined(DYNAMIC_DIRECTIONS)
475#  define ADVAL                adval[NUMBER_DIRECTIONS]
476#else
477#  define ADVAL                *adval;
478#endif
479#  define ADVAL_TYPE           const double *
480#  define FOR_I_EQ_0_LT_NUMDIR for (int _i=0; _i < ADOLC_numDir; ++_i)
481#  define ADVAL_I              adval[_i]
482#  define ADV_I                adv[_i]
483#  define V_I                  v[_i]
484#else
485#  define ADVAL                adval
486#  define ADVAL_TYPE           double
487#  define FOR_I_EQ_0_LT_NUMDIR
488#  define ADVAL_I              adval
489#  define ADV_I                adv
490#  define V_I                  v
491#endif
492
493inline double fmin( const double &x, const double &y ) {
494    if (x<y) return x;
495    else return y;
496}
497
498inline double fmax( const double &x, const double &y ) {
499    if (x>y) return x;
500    else return y;
501}
502
503inline double makeNaN() {
504#if defined(non_num)
505    double a,b;
506    a=non_num;
507    b=non_den;
508    return a/b;
509#else
510#  error Error: non_num undefined!
511#endif
512}
513
514class adouble {
515public:
516    // ctors
517    inline adouble();
518    inline adouble(const double v);
519    inline adouble(const double v, ADVAL_TYPE adv);
520    inline adouble(const adouble& a);
521#if defined(DYNAMIC_DIRECTIONS)
522    inline ~adouble();
523#endif
524
525    /*******************  temporary results  ******************************/
526    // sign
527    inline adouble operator - () const;
528    inline adouble operator + () const;
529
530    // addition
531    inline adouble operator + (const double v) const;
532    inline adouble operator + (const adouble& a) const;
533    inline friend
534    adouble operator + (const double v, const adouble& a);
535
536    // substraction
537    inline adouble operator - (const double v) const;
538    inline adouble operator - (const adouble& a) const;
539    inline friend
540    adouble operator - (const double v, const adouble& a);
541
542    // multiplication
543    inline adouble operator * (const double v) const;
544    inline adouble operator * (const adouble& a) const;
545    inline friend
546    adouble operator * (const double v, const adouble& a);
547
548    // division
549    inline adouble operator / (const double v) const;
550    inline adouble operator / (const adouble& a) const;
551    inline friend
552    adouble operator / (const double v, const adouble& a);
553
554    // inc/dec
555    inline adouble operator ++ ();
556    inline adouble operator ++ (int);
557    inline adouble operator -- ();
558    inline adouble operator -- (int);
559
560    // functions
561    inline friend adouble tan(const adouble &a);
562    inline friend adouble exp(const adouble &a);
563    inline friend adouble log(const adouble &a);
564    inline friend adouble sqrt(const adouble &a);
565    inline friend adouble sin(const adouble &a);
566    inline friend adouble cos(const adouble &a);
567    inline friend adouble asin(const adouble &a);
568    inline friend adouble acos(const adouble &a);
569    inline friend adouble atan(const adouble &a);
570
571    inline friend adouble atan2(const adouble &a, const adouble &b);
572    inline friend adouble pow(const adouble &a, double v);
573    inline friend adouble pow(const adouble &a, const adouble &b);
574    inline friend adouble pow(double v, const adouble &a);
575    inline friend adouble log10(const adouble &a);
576
577    inline friend adouble sinh (const adouble &a);
578    inline friend adouble cosh (const adouble &a);
579    inline friend adouble tanh (const adouble &a);
580#if defined(ATRIG_ERF)
581    inline friend adouble asinh (const adouble &a);
582    inline friend adouble acosh (const adouble &a);
583    inline friend adouble atanh (const adouble &a);
584#endif
585    inline friend adouble fabs (const adouble &a);
586    inline friend adouble ceil (const adouble &a);
587    inline friend adouble floor (const adouble &a);
588    inline friend adouble fmax (const adouble &a, const adouble &b);
589    inline friend adouble fmax (double v, const adouble &a);
590    inline friend adouble fmax (const adouble &a, double v);
591    inline friend adouble fmin (const adouble &a, const adouble &b);
592    inline friend adouble fmin (double v, const adouble &a);
593    inline friend adouble fmin (const adouble &a, double v);
594    inline friend adouble ldexp (const adouble &a, const adouble &b);
595    inline friend adouble ldexp (const adouble &a, const double v);
596    inline friend adouble ldexp (const double v, const adouble &a);
597    inline friend double frexp (const adouble &a, int* v);
598#if defined(ATRIG_ERF)
599    inline friend adouble erf (const adouble &a);
600#endif
601
602
603    /*******************  nontemporary results  ***************************/
604    // assignment
605    inline void operator = (const double v);
606    inline void operator = (const adouble& a);
607
608    // addition
609    inline void operator += (const double v);
610    inline void operator += (const adouble& a);
611
612    // substraction
613    inline void operator -= (const double v);
614    inline void operator -= (const adouble& a);
615
616    // multiplication
617    inline void operator *= (const double v);
618    inline void operator *= (const adouble& a);
619
620    // division
621    inline void operator /= (const double v);
622    inline void operator /= (const adouble& a);
623
624    // not
625    inline int operator ! () const;
626
627    // comparision
628    inline int operator != (const adouble&) const;
629    inline int operator != (const double) const;
630    inline friend int operator != (const double, const adouble&);
631
632    inline int operator == (const adouble&) const;
633    inline int operator == (const double) const;
634    inline friend int operator == (const double, const adouble&);
635
636    inline int operator <= (const adouble&) const;
637    inline int operator <= (const double) const;
638    inline friend int operator <= (const double, const adouble&);
639
640    inline int operator >= (const adouble&) const;
641    inline int operator >= (const double) const;
642    inline friend int operator >= (const double, const adouble&);
643
644    inline int operator >  (const adouble&) const;
645    inline int operator >  (const double) const;
646    inline friend int operator >  (const double, const adouble&);
647
648    inline int operator <  (const adouble&) const;
649    inline int operator <  (const double) const;
650    inline friend int operator <  (const double, const adouble&);
651
652    /*******************  getter / setter  ********************************/
653    inline double getValue() const;
654    inline void setValue(const double v);
655    inline ADVAL_TYPE getADValue() const;
656    inline void setADValue(ADVAL_TYPE v);
657#if defined(NUMBER_DIRECTIONS)
658    inline double getADValue(const unsigned int p) const;
659    inline void setADValue(const unsigned int p, const double v);
660#endif
661
662    /*******************  i/o operations  *********************************/
663    inline friend ostream& operator << ( ostream&, const adouble& );
664    inline friend istream& operator >> ( istream&, adouble& );
665
666
667private:
668    // internal variables
669    double val;
670    double ADVAL;
671};
672
673/*******************************  ctors  ************************************/
674adouble::adouble() {
675#if defined(DYNAMIC_DIRECTIONS)
676    adval = new double[ADOLC_numDir];
677#endif
678}
679
680adouble::adouble(const double v) : val(v) {
681#if defined(DYNAMIC_DIRECTIONS)
682    adval = new double[ADOLC_numDir];
683#endif
684    FOR_I_EQ_0_LT_NUMDIR
685    ADVAL_I = 0.0;
686}
687
688adouble::adouble(const double v, ADVAL_TYPE adv) : val(v) {
689#if defined(DYNAMIC_DIRECTIONS)
690    adval = new double[ADOLC_numDir];
691#endif
692    FOR_I_EQ_0_LT_NUMDIR
693    ADVAL_I=ADV_I;
694}
695
696adouble::adouble(const adouble& a) : val(a.val) {
697#if defined(DYNAMIC_DIRECTIONS)
698    adval = new double[ADOLC_numDir];
699#endif
700    FOR_I_EQ_0_LT_NUMDIR
701    ADVAL_I=a.ADVAL_I;
702}
703
704/*******************************  dtors  ************************************/
705#if defined(DYNAMIC_DIRECTIONS)
706adouble::~adouble() {
707    delete[] adval;
708}
709#endif
710
711/*************************  temporary results  ******************************/
712// sign
713adouble adouble::operator - () const {
714    adouble tmp;
715    tmp.val=-val;
716    FOR_I_EQ_0_LT_NUMDIR
717    tmp.ADVAL_I=-ADVAL_I;
718    return tmp;
719}
720
721adouble adouble::operator + () const {
722    return *this;
723}
724
725// addition
726adouble adouble::operator + (const double v) const {
727    return adouble(val+v, adval);
728}
729
730adouble adouble::operator + (const adouble& a) const {
731    adouble tmp;
732    tmp.val=val+a.val;
733    FOR_I_EQ_0_LT_NUMDIR
734    tmp.ADVAL_I=ADVAL_I+a.ADVAL_I;
735    return tmp;
736}
737
738adouble operator + (const double v, const adouble& a) {
739    return adouble(v+a.val, a.adval);
740}
741
742// subtraction
743adouble adouble::operator - (const double v) const {
744    return adouble(val-v, adval);
745}
746
747adouble adouble::operator - (const adouble& a) const {
748    adouble tmp;
749    tmp.val=val-a.val;
750    FOR_I_EQ_0_LT_NUMDIR
751    tmp.ADVAL_I=ADVAL_I-a.ADVAL_I;
752    return tmp;
753}
754
755adouble operator - (const double v, const adouble& a) {
756    adouble tmp;
757    tmp.val=v-a.val;
758    FOR_I_EQ_0_LT_NUMDIR
759    tmp.ADVAL_I=-a.ADVAL_I;
760    return tmp;
761}
762
763// multiplication
764adouble adouble::operator * (const double v) const {
765    adouble tmp;
766    tmp.val=val*v;
767    FOR_I_EQ_0_LT_NUMDIR
768    tmp.ADVAL_I=ADVAL_I*v;
769    return tmp;
770}
771
772adouble adouble::operator * (const adouble& a) const {
773    adouble tmp;
774    tmp.val=val*a.val;
775    FOR_I_EQ_0_LT_NUMDIR
776    tmp.ADVAL_I=ADVAL_I*a.val+val*a.ADVAL_I;
777    return tmp;
778}
779
780adouble operator * (const double v, const adouble& a) {
781    adouble tmp;
782    tmp.val=v*a.val;
783    FOR_I_EQ_0_LT_NUMDIR
784    tmp.ADVAL_I=v*a.ADVAL_I;
785    return tmp;
786}
787
788// division
789adouble adouble::operator / (const double v) const {
790    adouble tmp;
791    tmp.val=val/v;
792    FOR_I_EQ_0_LT_NUMDIR
793    tmp.ADVAL_I=ADVAL_I/v;
794    return tmp;
795}
796
797adouble adouble::operator / (const adouble& a) const {
798    adouble tmp;
799    tmp.val=val/a.val;
800    FOR_I_EQ_0_LT_NUMDIR
801    tmp.ADVAL_I=(ADVAL_I*a.val-val*a.ADVAL_I)/(a.val*a.val);
802    return tmp;
803}
804
805adouble operator / (const double v, const adouble& a) {
806    adouble tmp;
807    tmp.val=v/a.val;
808    FOR_I_EQ_0_LT_NUMDIR
809    tmp.ADVAL_I=(-v*a.ADVAL_I)/(a.val*a.val);
810    return tmp;
811}
812
813// inc/dec
814adouble adouble::operator ++ () {
815    ++val;
816    return *this;
817}
818
819adouble adouble::operator ++ (int) {
820    adouble tmp;
821    tmp.val=val++;
822    FOR_I_EQ_0_LT_NUMDIR
823    tmp.ADVAL_I=ADVAL_I;
824    return tmp;
825}
826
827adouble adouble::operator -- () {
828    --val;
829    return *this;
830}
831
832adouble adouble::operator -- (int) {
833    adouble tmp;
834    tmp.val=val--;
835    FOR_I_EQ_0_LT_NUMDIR
836    tmp.ADVAL_I=ADVAL_I;
837    return tmp;
838}
839
840// functions
841adouble tan(const adouble& a) {
842    adouble tmp;
843    double tmp2;
844    tmp.val=ADOLC_MATH_NSP::tan(a.val);
845    tmp2=ADOLC_MATH_NSP::cos(a.val);
846    tmp2*=tmp2;
847    FOR_I_EQ_0_LT_NUMDIR
848    tmp.ADVAL_I=a.ADVAL_I/tmp2;
849    return tmp;
850}
851
852adouble exp(const adouble &a) {
853    adouble tmp;
854    tmp.val=ADOLC_MATH_NSP::exp(a.val);
855    FOR_I_EQ_0_LT_NUMDIR
856    tmp.ADVAL_I=tmp.val*a.ADVAL_I;
857    return tmp;
858}
859
860adouble log(const adouble &a) {
861    adouble tmp;
862    tmp.val=ADOLC_MATH_NSP::log(a.val);
863    FOR_I_EQ_0_LT_NUMDIR
864    if (a.val>0 || a.val==0 && a.ADVAL_I>=0) tmp.ADVAL_I=a.ADVAL_I/a.val;
865    else tmp.ADVAL_I=makeNaN();
866    return tmp;
867}
868
869adouble sqrt(const adouble &a) {
870    adouble tmp;
871    tmp.val=ADOLC_MATH_NSP::sqrt(a.val);
872    FOR_I_EQ_0_LT_NUMDIR
873    if (a.val>0 || a.val==0 && a.ADVAL_I>=0) tmp.ADVAL_I=a.ADVAL_I/tmp.val/2;
874    else tmp.ADVAL_I=makeNaN();
875    return tmp;
876}
877
878adouble sin(const adouble &a) {
879    adouble tmp;
880    double tmp2;
881    tmp.val=ADOLC_MATH_NSP::sin(a.val);
882    tmp2=ADOLC_MATH_NSP::cos(a.val);
883    FOR_I_EQ_0_LT_NUMDIR
884    tmp.ADVAL_I=tmp2*a.ADVAL_I;
885    return tmp;
886}
887
888adouble cos(const adouble &a) {
889    adouble tmp;
890    double tmp2;
891    tmp.val=ADOLC_MATH_NSP::cos(a.val);
892    tmp2=-ADOLC_MATH_NSP::sin(a.val);
893    FOR_I_EQ_0_LT_NUMDIR
894    tmp.ADVAL_I=tmp2*a.ADVAL_I;
895    return tmp;
896}
897
898adouble asin(const adouble &a) {
899    adouble tmp;
900    tmp.val=ADOLC_MATH_NSP::asin(a.val);
901    double tmp2=ADOLC_MATH_NSP::sqrt(1-a.val*a.val);
902    FOR_I_EQ_0_LT_NUMDIR
903    tmp.ADVAL_I=a.ADVAL_I/tmp2;
904    return tmp;
905}
906
907adouble acos(const adouble &a) {
908    adouble tmp;
909    tmp.val=ADOLC_MATH_NSP::acos(a.val);
910    double tmp2=-ADOLC_MATH_NSP::sqrt(1-a.val*a.val);
911    FOR_I_EQ_0_LT_NUMDIR
912    tmp.ADVAL_I=a.ADVAL_I/tmp2;
913    return tmp;
914}
915
916adouble atan(const adouble &a) {
917    adouble tmp;
918    tmp.val=ADOLC_MATH_NSP::atan(a.val);
919    double tmp2=1+a.val*a.val;
920    tmp2=1/tmp2;
921    if (tmp2!=0)
922        FOR_I_EQ_0_LT_NUMDIR
923        tmp.ADVAL_I=a.ADVAL_I*tmp2;
924    else
925        FOR_I_EQ_0_LT_NUMDIR
926        tmp.ADVAL_I=0.0;
927    return tmp;
928}
929
930adouble atan2(const adouble &a, const adouble &b) {
931    adouble tmp;
932    tmp.val=ADOLC_MATH_NSP::atan2(a.val, b.val);
933    double tmp2=a.val*a.val;
934    double tmp3=b.val*b.val;
935    double tmp4=tmp3/(tmp2+tmp3);
936    if (tmp4!=0)
937        FOR_I_EQ_0_LT_NUMDIR
938        tmp.ADVAL_I=(a.ADVAL_I*b.val-a.val*b.ADVAL_I)/tmp3*tmp4;
939    else
940        FOR_I_EQ_0_LT_NUMDIR
941        tmp.ADVAL_I=0.0;
942    return tmp;
943}
944
945adouble pow(const adouble &a, double v) {
946    adouble tmp;
947    tmp.val=ADOLC_MATH_NSP::pow(a.val, v);
948    double tmp2=v*ADOLC_MATH_NSP::pow(a.val, v-1);
949    FOR_I_EQ_0_LT_NUMDIR
950    tmp.ADVAL_I=tmp2*a.ADVAL_I;
951    return tmp;
952}
953
954adouble pow(const adouble &a, const adouble &b) {
955    adouble tmp;
956    tmp.val=ADOLC_MATH_NSP::pow(a.val, b.val);
957    double tmp2=b.val*ADOLC_MATH_NSP::pow(a.val, b.val-1);
958    double tmp3=ADOLC_MATH_NSP::log(a.val)*tmp.val;
959    FOR_I_EQ_0_LT_NUMDIR
960    tmp.ADVAL_I=tmp2*a.ADVAL_I+tmp3*b.ADVAL_I;
961    return tmp;
962}
963
964adouble pow(double v, const adouble &a) {
965    adouble tmp;
966    tmp.val=ADOLC_MATH_NSP::pow(v, a.val);
967    double tmp2=tmp.val*ADOLC_MATH_NSP::log(v);
968    FOR_I_EQ_0_LT_NUMDIR
969    tmp.ADVAL_I=tmp2*a.ADVAL_I;
970    return tmp;
971}
972
973adouble log10(const adouble &a) {
974    adouble tmp;
975    tmp.val=ADOLC_MATH_NSP::log10(a.val);
976    double tmp2=ADOLC_MATH_NSP::log((double)10)*a.val;
977    FOR_I_EQ_0_LT_NUMDIR
978    tmp.ADVAL_I=a.ADVAL_I/tmp2;
979    return tmp;
980}
981
982adouble sinh (const adouble &a) {
983    adouble tmp;
984    tmp.val=ADOLC_MATH_NSP::sinh(a.val);
985    double tmp2=ADOLC_MATH_NSP::cosh(a.val);
986    FOR_I_EQ_0_LT_NUMDIR
987    tmp.ADVAL_I=a.ADVAL_I*tmp2;
988    return tmp;
989}
990
991adouble cosh (const adouble &a) {
992    adouble tmp;
993    tmp.val=ADOLC_MATH_NSP::cosh(a.val);
994    double tmp2=ADOLC_MATH_NSP::sinh(a.val);
995    FOR_I_EQ_0_LT_NUMDIR
996    tmp.ADVAL_I=a.ADVAL_I*tmp2;
997    return tmp;
998}
999
1000adouble tanh (const adouble &a) {
1001    adouble tmp;
1002    tmp.val=ADOLC_MATH_NSP::tanh(a.val);
1003    double tmp2=ADOLC_MATH_NSP::cosh(a.val);
1004    tmp2*=tmp2;
1005    FOR_I_EQ_0_LT_NUMDIR
1006    tmp.ADVAL_I=a.ADVAL_I/tmp2;
1007    return tmp;
1008}
1009
1010#if defined(ATRIG_ERF)
1011adouble asinh (const adouble &a) {
1012    adouble tmp;
1013    tmp.val=ADOLC_MATH_NSP::asinh(a.val);
1014    double tmp2=ADOLC_MATH_NSP::sqrt(a.val*a.val+1);
1015    FOR_I_EQ_0_LT_NUMDIR
1016    tmp.ADVAL_I=a.ADVAL_I/tmp2;
1017    return tmp;
1018}
1019
1020adouble acosh (const adouble &a) {
1021    adouble tmp;
1022    tmp.val=ADOLC_MATH_NSP::acosh(a.val);
1023    double tmp2=ADOLC_MATH_NSP::sqrt(a.val*a.val-1);
1024    FOR_I_EQ_0_LT_NUMDIR
1025    tmp.ADVAL_I=a.ADVAL_I/tmp2;
1026    return tmp;
1027}
1028
1029adouble atanh (const adouble &a) {
1030    adouble tmp;
1031    tmp.val=ADOLC_MATH_NSP::atanh(a.val);
1032    double tmp2=1-a.val*a.val;
1033    FOR_I_EQ_0_LT_NUMDIR
1034    tmp.ADVAL_I=a.ADVAL_I/tmp2;
1035    return tmp;
1036}
1037#endif
1038
1039adouble fabs (const adouble &a) {
1040    adouble tmp;
1041    tmp.val=ADOLC_MATH_NSP::fabs(a.val);
1042    int as=0;
1043    if (a.val>0) as=1;
1044    if (a.val<0) as=-1;
1045    if (as!=0)
1046        FOR_I_EQ_0_LT_NUMDIR
1047        tmp.ADVAL_I=a.ADVAL_I*as;
1048    else
1049        FOR_I_EQ_0_LT_NUMDIR {
1050            as=0;
1051            if (a.ADVAL_I>0) as=1;
1052            if (a.ADVAL_I<0) as=-1;
1053                tmp.ADVAL_I=a.ADVAL_I*as;
1054            }
1055            return tmp;
1056}
1057
1058adouble ceil (const adouble &a) {
1059    adouble tmp;
1060    tmp.val=ADOLC_MATH_NSP::ceil(a.val);
1061    FOR_I_EQ_0_LT_NUMDIR
1062    tmp.ADVAL_I=0.0;
1063    return tmp;
1064}
1065
1066adouble floor (const adouble &a) {
1067    adouble tmp;
1068    tmp.val=ADOLC_MATH_NSP::floor(a.val);
1069    FOR_I_EQ_0_LT_NUMDIR
1070    tmp.ADVAL_I=0.0;
1071    return tmp;
1072}
1073
1074adouble fmax (const adouble &a, const adouble &b) {
1075    adouble tmp;
1076    double tmp2=a.val-b.val;
1077    if (tmp2<0) {
1078        tmp.val=b.val;
1079        FOR_I_EQ_0_LT_NUMDIR
1080        tmp.ADVAL_I=b.ADVAL_I;
1081    } else {
1082        tmp.val=a.val;
1083        if (tmp2>0) {
1084            FOR_I_EQ_0_LT_NUMDIR
1085            tmp.ADVAL_I=a.ADVAL_I;
1086        } else {
1087            FOR_I_EQ_0_LT_NUMDIR
1088            {
1089                if (a.ADVAL_I<b.ADVAL_I) tmp.ADVAL_I=b.ADVAL_I;
1090                else tmp.ADVAL_I=a.ADVAL_I;
1091                }
1092            }
1093}
1094return tmp;
1095}
1096
1097adouble fmax (double v, const adouble &a) {
1098    adouble tmp;
1099    double tmp2=v-a.val;
1100    if (tmp2<0) {
1101        tmp.val=a.val;
1102        FOR_I_EQ_0_LT_NUMDIR
1103        tmp.ADVAL_I=a.ADVAL_I;
1104    } else {
1105        tmp.val=v;
1106        if (tmp2>0) {
1107            FOR_I_EQ_0_LT_NUMDIR
1108            tmp.ADVAL_I=0.0;
1109        } else {
1110            FOR_I_EQ_0_LT_NUMDIR
1111            {
1112                if (a.ADVAL_I>0) tmp.ADVAL_I=a.ADVAL_I;
1113                else tmp.ADVAL_I=0.0;
1114                }
1115            }
1116}
1117return tmp;
1118}
1119
1120adouble fmax (const adouble &a, double v) {
1121    adouble tmp;
1122    double tmp2=a.val-v;
1123    if (tmp2<0) {
1124        tmp.val=v;
1125        FOR_I_EQ_0_LT_NUMDIR
1126        tmp.ADVAL_I=0.0;
1127    } else {
1128        tmp.val=a.val;
1129        if (tmp2>0) {
1130            FOR_I_EQ_0_LT_NUMDIR
1131            tmp.ADVAL_I=a.ADVAL_I;
1132        } else {
1133            FOR_I_EQ_0_LT_NUMDIR
1134            {
1135                if (a.ADVAL_I>0) tmp.ADVAL_I=a.ADVAL_I;
1136                else tmp.ADVAL_I=0.0;
1137                }
1138            }
1139}
1140return tmp;
1141}
1142
1143adouble fmin (const adouble &a, const adouble &b) {
1144    adouble tmp;
1145    double tmp2=a.val-b.val;
1146    if (tmp2<0) {
1147        tmp.val=a.val;
1148        FOR_I_EQ_0_LT_NUMDIR
1149        tmp.ADVAL_I=a.ADVAL_I;
1150    } else {
1151        tmp.val=b.val;
1152        if (tmp2>0) {
1153            FOR_I_EQ_0_LT_NUMDIR
1154            tmp.ADVAL_I=b.ADVAL_I;
1155        } else {
1156            FOR_I_EQ_0_LT_NUMDIR
1157            {
1158                if (a.ADVAL_I<b.ADVAL_I) tmp.ADVAL_I=a.ADVAL_I;
1159                else tmp.ADVAL_I=b.ADVAL_I;
1160                }
1161            }
1162}
1163return tmp;
1164}
1165
1166adouble fmin (double v, const adouble &a) {
1167    adouble tmp;
1168    double tmp2=v-a.val;
1169    if (tmp2<0) {
1170        tmp.val=v;
1171        FOR_I_EQ_0_LT_NUMDIR
1172        tmp.ADVAL_I=0.0;
1173    } else {
1174        tmp.val=a.val;
1175        if (tmp2>0) {
1176            FOR_I_EQ_0_LT_NUMDIR
1177            tmp.ADVAL_I=a.ADVAL_I;
1178        } else {
1179            FOR_I_EQ_0_LT_NUMDIR
1180            {
1181                if (a.ADVAL_I<0) tmp.ADVAL_I=a.ADVAL_I;
1182                else tmp.ADVAL_I=0.0;
1183                }
1184            }
1185}
1186return tmp;
1187}
1188
1189adouble fmin (const adouble &a, double v) {
1190    adouble tmp;
1191    double tmp2=a.val-v;
1192    if (tmp2<0) {
1193        tmp.val=a.val;
1194        FOR_I_EQ_0_LT_NUMDIR
1195        tmp.ADVAL_I=a.ADVAL_I;
1196    } else {
1197        tmp.val=v;
1198        if (tmp2>0) {
1199            FOR_I_EQ_0_LT_NUMDIR
1200            tmp.ADVAL_I=0.0;
1201        } else {
1202            FOR_I_EQ_0_LT_NUMDIR
1203            {
1204                if (a.ADVAL_I<0) tmp.ADVAL_I=a.ADVAL_I;
1205                else tmp.ADVAL_I=0.0;
1206                }
1207            }
1208}
1209return tmp;
1210}
1211
1212adouble ldexp (const adouble &a, const adouble &b) {
1213    return a*pow(2.,b);
1214}
1215
1216adouble ldexp (const adouble &a, const double v) {
1217    return a*ADOLC_MATH_NSP::pow(2.,v);
1218}
1219
1220adouble ldexp (const double v, const adouble &a) {
1221    return v*pow(2.,a);
1222}
1223
1224double frexp (const adouble &a, int* v) {
1225    return ADOLC_MATH_NSP::frexp(a.val, v);
1226}
1227
1228#if defined(ATRIG_ERF)
1229adouble erf (const adouble &a) {
1230    adouble tmp;
1231    tmp.val=ADOLC_MATH_NSP::erf(a.val);
1232    double tmp2 = 2.0 /
1233        ADOLC_MATH_NSP::sqrt(ADOLC_MATH_NSP::acos(-1.0)) *
1234        ADOLC_MATH_NSP::exp(-a.val*a.val);
1235    FOR_I_EQ_0_LT_NUMDIR
1236    tmp.ADVAL_I=tmp2*a.ADVAL_I;
1237    return tmp;
1238}
1239#endif
1240
1241
1242/*******************  nontemporary results  *********************************/
1243void adouble::operator = (const double v) {
1244    val=v;
1245    FOR_I_EQ_0_LT_NUMDIR
1246    ADVAL_I=0.0;
1247}
1248
1249void adouble::operator = (const adouble& a) {
1250    val=a.val;
1251    FOR_I_EQ_0_LT_NUMDIR
1252    ADVAL_I=a.ADVAL_I;
1253}
1254
1255void adouble::operator += (const double v) {
1256    val+=v;
1257}
1258
1259void adouble::operator += (const adouble& a) {
1260    val=val+a.val;
1261    FOR_I_EQ_0_LT_NUMDIR
1262    ADVAL_I+=a.ADVAL_I;
1263}
1264
1265void adouble::operator -= (const double v) {
1266    val-=v;
1267}
1268
1269void adouble::operator -= (const adouble& a) {
1270    val=val-a.val;
1271    FOR_I_EQ_0_LT_NUMDIR
1272    ADVAL_I-=a.ADVAL_I;
1273}
1274
1275void adouble::operator *= (const double v) {
1276    val=val*v;
1277    FOR_I_EQ_0_LT_NUMDIR
1278    ADVAL_I*=v;
1279}
1280
1281void adouble::operator *= (const adouble& a) {
1282    FOR_I_EQ_0_LT_NUMDIR
1283    ADVAL_I=ADVAL_I*a.val+val*a.ADVAL_I;
1284    val*=a.val;
1285}
1286
1287void adouble::operator /= (const double v) {
1288    val/=v;
1289    FOR_I_EQ_0_LT_NUMDIR
1290    ADVAL_I/=v;
1291}
1292
1293void adouble::operator /= (const adouble& a) {
1294    FOR_I_EQ_0_LT_NUMDIR
1295    ADVAL_I=(ADVAL_I*a.val-val*a.ADVAL_I)/(a.val*a.val);
1296    val=val/a.val;
1297}
1298
1299// not
1300int adouble::operator ! () const {
1301    return val==0.0;
1302}
1303
1304// comparision
1305int adouble::operator != (const adouble &a) const {
1306    return val!=a.val;
1307}
1308
1309int adouble::operator != (const double v) const {
1310    return val!=v;
1311}
1312
1313int operator != (const double v, const adouble &a) {
1314    return v!=a.val;
1315}
1316
1317int adouble::operator == (const adouble &a) const {
1318    return val==a.val;
1319}
1320
1321int adouble::operator == (const double v) const {
1322    return val==v;
1323}
1324
1325int operator == (const double v, const adouble &a) {
1326    return v==a.val;
1327}
1328
1329int adouble::operator <= (const adouble &a) const {
1330    return val<=a.val;
1331}
1332
1333int adouble::operator <= (const double v) const {
1334    return val<=v;
1335}
1336
1337int operator <= (const double v, const adouble &a) {
1338    return v<=a.val;
1339}
1340
1341int adouble::operator >= (const adouble &a) const {
1342    return val>=a.val;
1343}
1344
1345int adouble::operator >= (const double v) const {
1346    return val>=v;
1347}
1348
1349int operator >= (const double v, const adouble &a) {
1350    return v>=a.val;
1351}
1352
1353int adouble::operator >  (const adouble &a) const {
1354    return val>a.val;
1355}
1356
1357int adouble::operator >  (const double v) const {
1358    return val>v;
1359}
1360
1361int operator >  (const double v, const adouble &a) {
1362    return v>a.val;
1363}
1364
1365int adouble::operator <  (const adouble &a) const {
1366    return val<a.val;
1367}
1368
1369int adouble::operator <  (const double v) const {
1370    return val<v;
1371}
1372
1373int operator <  (const double v, const adouble &a) {
1374    return v<a.val;
1375}
1376
1377/*******************  getter / setter  **************************************/
1378double adouble::getValue() const {
1379    return val;
1380}
1381
1382void adouble::setValue(const double v) {
1383    val=v;
1384}
1385
1386ADVAL_TYPE adouble::getADValue() const {
1387    return adval;
1388}
1389
1390void adouble::setADValue(ADVAL_TYPE v) {
1391    FOR_I_EQ_0_LT_NUMDIR
1392    ADVAL_I=V_I;
1393}
1394
1395#  if defined(NUMBER_DIRECTIONS)
1396double adouble::getADValue(const unsigned int p) const {
1397    if (p>=NUMBER_DIRECTIONS) {
1398        fprintf(DIAG_OUT, "Derivative array accessed out of bounds"\
1399                " while \"getADValue(...)\"!!!\n");
1400        exit(-1);
1401    }
1402    return adval[p];
1403}
1404
1405void adouble::setADValue(const unsigned int p, const double v) {
1406    if (p>=NUMBER_DIRECTIONS) {
1407        fprintf(DIAG_OUT, "Derivative array accessed out of bounds"\
1408                " while \"setADValue(...)\"!!!\n");
1409        exit(-1);
1410    }
1411    adval[p]=v;
1412}
1413#  endif
1414
1415#if defined(NUMBER_DIRECTIONS)
1416static void setNumDir(const unsigned int p) {
1417#if !defined(DYNAMIC_DIRECTIONS)
1418    if (p>NUMBER_DIRECTIONS) ADOLC_numDir=NUMBER_DIRECTIONS;
1419    else ADOLC_numDir=p;
1420#else
1421    ADOLC_numDir = p;
1422#endif
1423}
1424#endif
1425
1426/*******************  i/o operations  ***************************************/
1427ostream& operator << ( ostream& out, const adouble& a) {
1428    out << "Value: " << a.val;
1429#if !defined(NUMBER_DIRECTIONS)
1430    out << " ADValue: ";
1431#else
1432    out << " ADValues (" << ADOLC_numDir << "): ";
1433#endif
1434    FOR_I_EQ_0_LT_NUMDIR
1435    out << a.ADVAL_I << " ";
1436    out << "(a)";
1437    return out;
1438}
1439
1440istream& operator >> ( istream& in, adouble& a) {
1441    char c;
1442    do {
1443        in >> c;
1444    } while (c!=':' && !in.eof());
1445    in >> a.val;
1446#if !defined(NUMBER_DIRECTIONS)
1447    do in >> c;
1448    while (c!=':' && !in.eof());
1449#else
1450unsigned int num;
1451do in >> c;
1452while (c!='(' && !in.eof());
1453in >> num;
1454if (num>NUMBER_DIRECTIONS) {
1455    cout << "ADOL-C error: to many directions in input\n";
1456    exit(-1);
1457}
1458do in >> c;
1459while (c!=')' && !in.eof());
1460#endif
1461    FOR_I_EQ_0_LT_NUMDIR
1462    in >> a.ADVAL_I;
1463    do in >> c;
1464    while (c!=')' && !in.eof());
1465    return in;
1466}
1467}
1468
1469/****************************************************************************/
1470#endif /* ADOLC_TAPELESS */
1471
1472/****************************************************************************/
1473/*                                                                THAT'S ALL*/
1474#endif /* __cplusplus */
1475#endif /* ADOLC_ADOUBLE_H */
Note: See TracBrowser for help on using the repository browser.