source: trunk/ADOL-C/src/uni5_for.c @ 372

Last change on this file since 372 was 372, checked in by kulshres, 8 years ago

Reintroduce the old method of hessian pattern propagation

This reintroduces the hessian pattern propagation code that was in
previous releases under new interfaces

nonl_ind_old_forward_t() and nonl_ind_old_forward_s()

The standard interfaces nonl_ind_forward_t() and nonl_ind_forward_s()
rely on a newer approach introduced by
Andrea Walther <andrea.wather@…> since May 2012

Signed-off-by: Kshitij Kulshreshtha <kshitij@…>

  • Property svn:keywords set to Author Date Id Revision
File size: 162.4 KB
Line 
1
2/*----------------------------------------------------------------------------
3 ADOL-C -- Automatic Differentiation by Overloading in C++
4 File:     uni5_for.c
5
6
7 Revision: $Id: uni5_for.c 372 2012-12-07 12:24:46Z kulshres $
8
9 Contents: Contains the routines :
10           zos_forward (zero-order-scalar forward mode):      define _ZOS_
11           fos_forward (first-order-scalar forward mode):     define _FOS_
12           hos_forward (higher-order-scalar forward mode):    define _HOS_
13           fov_forward (first-order-vector forward mode):     define _FOV_
14           hov_forward (higher-order-vector forward mode):    define _HOV_
15           hov_wk_forward (higher-order-vector forward mode): define _HOV_WK_
16           int_forward_safe:                                  define _INT_FOR_ and _NTIGHT__
17
18           Uses the preprocessor to compile the 7 different object files
19           with/without "keep" parameter:                     define _KEEP_
20
21 Copyright (c) Andrea Walther, Andreas Griewank, Andreas Kowarz,
22               Hristo Mitev, Sebastian Schlenkrich, Jean Utke, Olaf Vogel,
23               Benjamin Letschert, Kshitij Kulshreshtha
24               
25
26 This file is part of ADOL-C. This software is provided as open source.
27 Any use, reproduction, or distribution of the software constitutes
28 recipient's acceptance of the terms of the accompanying license file.
29
30----------------------------------------------------------------------------*/
31
32#include <adolc/interfaces.h>
33#include <adolc/adalloc.h>
34#include <adolc/taping.h>
35#include "taping_p.h"
36#include "oplate.h"
37#include <adolc/externfcts.h>
38#include "externfcts_p.h"
39
40#include <math.h>
41
42#if defined(ADOLC_DEBUG)
43#include <string.h>
44#endif /* ADOLC_DEBUG */
45
46
47/****************************************************************************/
48/*                                                                   MACROS */
49#undef _ADOLC_VECTOR_
50#undef _HIGHER_ORDER_
51
52
53
54/*--------------------------------------------------------------------------*/
55#if defined(_ZOS_)
56#  define GENERATED_FILENAME "zos_forward"
57
58/*--------------------------------------------------------------------------*/
59#else
60#if defined(_FOS_)
61#define GENERATED_FILENAME "fos_forward"
62
63#define ARGUMENT(indexi,l,i) argument[indexi]
64#define TAYLORS(indexd,l,i)   taylors[indexd]
65
66/*--------------------------------------------------------------------------*/
67#else
68#if defined(_FOV_)
69#define GENERATED_FILENAME "fov_forward"
70
71#define _ADOLC_VECTOR_
72
73#if defined(_CHUNKED_)
74#define ARGUMENT(indexi,l,i) argument[indexi][l+offset]
75#define TAYLORS(indexd,l,i)   taylors[indexd][l+offset]
76#else
77#define ARGUMENT(indexi,l,i) argument[indexi][l]
78#define TAYLORS(indexd,l,i)   taylors[indexd][l]
79#endif
80
81/*--------------------------------------------------------------------------*/
82#else
83#if defined(_HOS_)
84#define GENERATED_FILENAME "hos_forward"
85
86#define _HIGHER_ORDER_
87
88#define ARGUMENT(indexi,l,i) argument[indexi][i]
89#define TAYLORS(indexd,l,i)   taylors[indexd][i]
90
91/*--------------------------------------------------------------------------*/
92#else
93#if defined(_HOV_)
94#define GENERATED_FILENAME "hov_forward"
95
96#define _ADOLC_VECTOR_
97#define _HIGHER_ORDER_
98
99#define ARGUMENT(indexi,l,i) argument[indexi][l][i]
100#define TAYLORS(indexd,l,i)   taylors[indexd][l][i]
101
102/*--------------------------------------------------------------------------*/
103#else
104#if defined(_HOV_WK_)
105#define GENERATED_FILENAME "hov_wk_forward"
106
107#define _ADOLC_VECTOR_
108#define _HIGHER_ORDER_
109
110#define ARGUMENT(indexi,l,i) argument[indexi][l][i]
111#define TAYLORS(indexd,l,i)   taylors[indexd][l][i]
112
113/*--------------------------------------------------------------------------*/
114#else
115#if defined(_INT_FOR_)
116#if defined(_TIGHT_)
117#define GENERATED_FILENAME "int_forward_t"
118#endif
119#if defined(_NTIGHT_)
120#define GENERATED_FILENAME "int_forward_s"
121#endif
122#define ARGUMENT(indexi,l,i) argument[indexi][l]
123#define TAYLORS(indexd,l,i)   taylors[indexd][l]
124/*--------------------------------------------------------------------------*/
125#else
126#if defined(_INDO_)
127void copy_index_domain(int res, int arg, locint **ind_dom);
128void merge_2_index_domains(int res, int arg, locint **ind_dom);
129void combine_2_index_domains(int res, int arg1, int arg2, locint **ind_dom);
130void merge_3_index_domains(int res, int arg1, int arg2, locint **ind_dom);
131
132#define NUMNNZ 20
133#define FMIN_ADOLC(x,y)  ((y<x)?y:x)
134
135#if defined(_INDOPRO_) && !defined(_NONLIND_OLD_)
136#if defined(_TIGHT_)
137#define GENERATED_FILENAME "indopro_forward_t"
138#endif
139#if defined(_NTIGHT_)
140#define GENERATED_FILENAME "indopro_forward_s"
141#endif
142#endif
143#if defined(_NONLIND_)
144
145typedef struct IndexElement {
146    locint  entry;
147    struct IndexElement* left;
148    struct IndexElement* right;
149}
150IndexElement;
151
152typedef struct IndexElement_sod {
153    locint  entry;
154    struct IndexElement_sod* left;
155}
156IndexElement_sod;
157
158void free_tree(IndexElement* nlf, int num);
159void traverse_crs(IndexElement* fod2, IndexElement_sod* sod, int num);
160void traverse_unary(IndexElement* fod, IndexElement* nonl_dom,  IndexElement* fodi, int num, int maxopind);
161
162#if defined(_TIGHT_)
163#define GENERATED_FILENAME "nonl_ind_forward_t"
164#elif defined(_NTIGHT_)
165#define GENERATED_FILENAME "nonl_ind_forward_s"
166#endif
167#endif
168#if defined(_NONLIND_OLD_)
169
170/*
171 * This is the type used for the list elements. The entry is either a counter
172 * (first element of the NID list) or the index of an independent variable.
173 */
174
175void extend_nonlinearity_domain_binary_step
176(int arg1, int arg2, locint **ind_dom, locint **nonl_dom);
177void extend_nonlinearity_domain_unary
178(int arg, locint **ind_dom, locint **nonl_dom);
179void extend_nonlinearity_domain_binary
180(int arg1, int arg2, locint **ind_dom, locint **nonl_dom);
181
182
183#if defined(_TIGHT_)
184#define GENERATED_FILENAME "nonl_ind_old_forward_t"
185#endif
186#if defined(_NTIGHT_)
187#define GENERATED_FILENAME "nonl_ind_old_forward_s"
188#endif
189#endif
190
191
192/*--------------------------------------------------------------------------*/
193#else
194#error Error ! Define [_ZOS_ | _FOS_ |\
195   _HOS_ | _FOV_ | _HOV_ | _HOV_WK_  | _INT_FOR_SAFE_ | _INT_FOR_TIGHT_ | _INDOPRO_ | _NONLIND_ ] [{_KEEP_}]
196#endif
197#endif
198#endif
199#endif
200#endif
201#endif
202#endif
203#endif
204
205/*--------------------------------------------------------------------------*/
206/*                                                               KEEP stuff */
207#if defined(_KEEP_)
208
209#if defined(_HOV_WK_) /* keep in this vector mode */
210#define IF_KEEP_TAYLOR_CLOSE \
211if (keep){\
212  fprintf(DIAG_OUT,"Succeeding reverse sweep will fail!\n");\
213  taylor_close(0);\
214}
215#define IF_KEEP_WRITE_TAYLOR(res,keep,k,p) \
216    { \
217        UPDATE_TAYLORWRITTEN(keep * k * p) \
218        if (keep) \
219        { \
220            ADOLC_WRITE_SCAYLOR(dp_T0[res]); \
221            if (keep > 1) \
222            write_taylors(res,(keep-1),k,p); \
223        } \
224    }
225#else
226#if defined(_ADOLC_VECTOR_) /* otherwise no keep */
227#define IF_KEEP_TAYLOR_CLOSE
228#define IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
229#else /* _ZOS_, _FOS_, _HOS_ */
230#define IF_KEEP_TAYLOR_CLOSE \
231if (keep){\
232  fprintf(DIAG_OUT,"Otherwise succeeding reverse sweep will fail!\n");\
233  taylor_close(0);\
234}
235#if defined(_ZOS_)
236#define IF_KEEP_WRITE_TAYLOR(res,keep,k,p) \
237    { \
238        UPDATE_TAYLORWRITTEN(keep) \
239        if (keep) \
240            ADOLC_WRITE_SCAYLOR(dp_T0[res]); \
241    }
242#else
243#if defined(_FOS_)
244#define IF_KEEP_WRITE_TAYLOR(res,keep,k,p) \
245    { \
246        UPDATE_TAYLORWRITTEN(keep) \
247        if (keep) \
248        { \
249            ADOLC_WRITE_SCAYLOR(dp_T0[res]); \
250            if (keep > 1) \
251                ADOLC_WRITE_SCAYLOR(dp_T[res]); \
252        } \
253    }
254#else
255#if defined(_HOS_)
256#define IF_KEEP_WRITE_TAYLOR(res,keep,k,p) \
257    { \
258        UPDATE_TAYLORWRITTEN(keep) \
259        if (keep) \
260        { \
261            ADOLC_WRITE_SCAYLOR(dp_T0[res]); \
262            if (keep > 1) \
263                write_taylor(res,keep-1); \
264        } \
265    }
266#endif
267#endif
268#endif
269#endif
270#endif
271
272#else  /* no _KEEP_ */
273#define IF_KEEP_TAYLOR_CLOSE
274#define IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
275#endif
276
277/*--------------------------------------------------------------------------*/
278/*                                                      access to variables */
279#if !defined(_ZOS_)
280#if defined(_FOS_)
281#define TRES         *Tres
282#define TARG         *Targ
283#define TARG1        *Targ1
284#define TARG2        *Targ2
285#define TQO          *Tqo
286
287#define TRES_INC     *Tres
288#define TARG_INC     *Targ
289#define TARG1_INC    *Targ1
290#define TARG2_INC    *Targ2
291#define TQO_INC      *Tqo
292
293#define TRES_DEC     *Tres
294#define TARG_DEC     *Targ
295#define TARG1_DEC    *Targ1
296#define TARG2_DEC    *Targ2
297#define TQO_DEC      *Tqo
298
299#define TRES_FOINC   *Tres
300#define TARG_FOINC   *Targ
301#define TARG1_FOINC  *Targ1
302#define TARG2_FOINC  *Targ2
303#define TQO_FOINC    *Tqo
304
305#define TRES_FODEC   *Tres
306#define DEC_TRES_FO
307#define TARG_FODEC   *Targ
308#define TARG1_FODEC  *Targ1
309#define TARG2_FODEC  *Targ2
310#define TQO_FODEC    *Tqo
311
312#define ASSIGN_T(a,b)  a = &b;
313
314#else
315#if defined(_INT_FOR_)
316#define TRES         *Tres
317#define TARG         *Targ
318#define TARG1        *Targ1
319#define TARG2        *Targ2
320#define TQO          *Tqo
321
322#define TRES_INC     *Tres++
323#define TARG_INC     *Targ++
324#define TARG1_INC    *Targ1++
325#define TARG2_INC    *Targ2++
326#define TQO_INC      *Tqo++
327
328#define TRES_DEC     *Tres--
329#define TARG_DEC     *Targ--
330#define TARG1_DEC    *Targ1--
331#define TARG2_DEC    *Targ2--
332#define TQO_DEC      *Tqo--
333
334#define TRES_FOINC   *Tres++
335#define TARG_FOINC   *Targ++
336#define TARG1_FOINC  *Targ1++
337#define TARG2_FOINC  *Targ2++
338#define TQO_FOINC    *Tqo++
339
340#define TRES_FODEC   *Tres--
341#define TARG_FODEC   *Targ--
342#define TARG1_FODEC  *Targ1--
343#define TARG2_FODEC  *Targ2--
344#define TQO_FODEC    *Tqo--
345
346
347#define ASSIGN_T(a,b)  a = b;
348
349#else  /* _HOS_, _FOV_, _HOV_, _HOV_WK */
350#define TRES         *Tres
351#define TARG         *Targ
352#define TARG1        *Targ1
353#define TARG2        *Targ2
354#define TQO          *Tqo
355
356#define TRES_INC     *Tres++
357#define TARG_INC     *Targ++
358#define TARG1_INC    *Targ1++
359#define TARG2_INC    *Targ2++
360#define TQO_INC      *Tqo++
361
362#define TRES_DEC     *Tres--
363#define TARG_DEC     *Targ--
364#define TARG1_DEC    *Targ1--
365#define TARG2_DEC    *Targ2--
366#define TQO_DEC      *Tqo--
367
368#if defined(_FOV_)
369#define TRES_FOINC   *Tres++
370#define TARG_FOINC   *Targ++
371#define TARG1_FOINC  *Targ1++
372#define TARG2_FOINC  *Targ2++
373#define TQO_FOINC    *Tqo++
374
375#define TRES_FODEC   *Tres
376#define DEC_TRES_FO  Tres--;
377#define TARG_FODEC   *Targ--
378#define TARG1_FODEC  *Targ1--
379#define TARG2_FODEC  *Targ2--
380#define TQO_FODEC    *Tqo--
381#else /* _HOS_, _HOV_, _HOV_WK */
382#define TRES_FOINC   *Tres
383#define TARG_FOINC   *Targ
384#define TARG1_FOINC  *Targ1
385#define TARG2_FOINC  *Targ2
386#define TQO_FOINC    *Tqo
387
388#define TRES_FODEC   *Tres
389#define DEC_TRES_FO
390#define TARG_FODEC   *Targ
391#define TARG1_FODEC  *Targ1
392#define TARG2_FODEC  *Targ2
393#define TQO_FODEC    *Tqo
394#endif
395#endif
396
397#define ASSIGN_T(a,b)  a = b;
398#endif
399#endif
400
401
402/*--------------------------------------------------------------------------*/
403/*                                                               loop stuff */
404#if defined(_ADOLC_VECTOR_)
405#define FOR_0_LE_l_LT_p for (l=0; l<p; l++)
406#define FOR_p_GT_l_GE_0 for (l=p-1; l>=0; l--)
407#else
408#if defined(_INT_FOR_)
409#define FOR_0_LE_l_LT_p for (l=0; l<p; l++)
410#define FOR_p_GT_l_GE_0 for (l=p-1; l>=0; l--)
411#else
412#define FOR_0_LE_l_LT_p
413#define FOR_p_GT_l_GE_0
414#endif
415#endif
416
417#if defined(_HIGHER_ORDER_)
418#define FOR_0_LE_i_LT_k for (i=0; i<k; i++)
419#define FOR_k_GT_i_GE_0 for (i=k-1; i>=0; i--)
420#else
421#define FOR_0_LE_i_LT_k
422#define FOR_k_GT_i_GE_0
423#endif
424
425#if defined(_HOV_)
426#define FOR_0_LE_l_LT_pk for (l=0; l<pk; l++)
427#define INC_pk_1(T)      T += pk-1;
428#define VEC_INC(T,inc)   T += inc;
429#define HOV_INC(T,inc)   T += inc;
430#else
431#if defined(_HOV_WK_)
432#define FOR_0_LE_l_LT_pk for (l=0; l<pk; l++)
433#define INC_pk_1(T)      T += pk-1;
434#define VEC_INC(T,inc)   T += inc;
435#define HOV_INC(T,inc)   T += inc;
436#else
437#if defined(_FOV_)
438#define FOR_0_LE_l_LT_pk for (l=0; l<p; l++)
439#define INC_pk_1(T)      T += p-1;
440#define VEC_INC(T,inc)   T++;
441#define HOV_INC(T,inc)
442#else
443#if defined(_HOS_)
444#define FOR_0_LE_l_LT_pk for (l=0; l<k; l++)
445#define INC_pk_1(T)      T += k-1;
446#define VEC_INC(T,inc)
447#define HOV_INC(T,inc)
448#else
449#if defined(_INT_FOR_)
450#define FOR_0_LE_l_LT_pk for (l=0; l<p; l++)
451#define INC_pk_1(T)      T += p-1;
452#define VEC_INC(T,inc)   T++;
453#else
454#define FOR_0_LE_l_LT_pk
455#define INC_pk_1(T)
456#define VEC_INC(T,inc)
457#define HOV_INC(T,inc)
458#endif
459#endif
460#endif
461#endif
462#endif
463
464/*--------------------------------------------------------------------------*/
465/*                                                        higher order case */
466#if defined(_HIGHER_ORDER_)
467#define BREAK_FOR_I break;
468#else
469#define BREAK_FOR_I ;
470#endif
471
472/* END Macros */
473
474BEGIN_C_DECLS
475
476#if defined(_ZOS_)
477/****************************************************************************/
478/* Zero Order Scalar version of the forward mode.                           */
479/****************************************************************************/
480#if defined(_KEEP_)
481int  zos_forward(
482#else
483int  zos_forward_nk(
484#endif
485    short  tnum,              /* tape id */
486    int    depcheck,          /* consistency chk on # of deps */
487    int    indcheck,          /* consistency chk on # of indeps */
488#if defined(_KEEP_)
489    int    keep,              /* flag for reverse sweep */
490#endif
491    const double *basepoint,  /* independant variable values */
492    double       *valuepoint) /* dependent variable values */
493
494#else
495#if defined(_FOS_)
496/****************************************************************************/
497/* First Order Scalar version of the forward mode.                          */
498/****************************************************************************/
499#if defined(_KEEP_)
500int  fos_forward(
501#else
502int  fos_forward_nk(
503#endif
504    short  tnum,        /* tape id */
505    int    depcheck,    /* consistency chk on # of deps */
506    int    indcheck,    /* consistency chk on # of indeps */
507#if defined(_KEEP_)
508    int    keep,        /* flag for reverse sweep */
509#endif
510    const double *basepoint,  /* independent variable values */
511    double *argument,   /* Taylor coefficients (input) */
512    double *valuepoint, /* Taylor coefficients (output) */
513    double *taylors)    /* matrix of coefficient vectors */
514/* the order of the indices in argument and taylors is [var][taylor] */
515
516#else
517#if defined(_INT_FOR_)
518#if defined(_TIGHT_)
519/****************************************************************************/
520/* First Order Vector version of the forward mode for bit patterns, tight   */
521/****************************************************************************/
522int int_forward_tight(
523    short               tnum,     /* tape id                              */
524    int                 depcheck, /* consistency chk on # of dependents   */
525    int                 indcheck, /* consistency chk on # of independents */
526    int                 p,        /* # of taylor series, bit pattern      */
527    const double       *basepoint,  /* independent variable values   (in)*/
528    unsigned long int **argument,  /* Taylor coeff.                 (in)*/
529    double             *valuepoint, /* dependent variable values    (out)*/
530    unsigned long int **taylors)   /* matrix of coefficient vectors(out)*/
531
532/* int_forward_tight( tag, m, n, p, x[n], X[n][p], y[m], Y[m][p]),
533
534     nBV = number of Boolean Vectors to be packed
535                      (see Chapter Dependence Analysis, ADOL-C Documentation)
536     bits_per_long = 8*sizeof(unsigned long int)
537     p = nBV / bits_per_long + ( (nBV % bits_per_long) != 0 )
538
539     The order of the indices in argument and taylors is [var][taylor]
540
541     For the full Jacobian matrix set
542     p = indep / bits_per_long + ((indep % bits_per_long) != 0)
543     and pass a bit pattern version of the identity matrix as an argument   */
544
545
546#endif
547#if defined (_NTIGHT_)
548/****************************************************************************/
549/* First Order Vector version of the forward mode, bit pattern, safe        */
550/****************************************************************************/
551int int_forward_safe(
552    short             tnum,     /* tape id                              */
553    int               depcheck, /* consistency chk on # of dependents   */
554    int               indcheck, /* consistency chk on # of independents */
555    int               p,        /* # of taylor series, bit pattern      */
556    unsigned long int **argument, /* Taylor coeff.                  (in)*/
557    unsigned long int **taylors)  /* matrix of coefficient vectors (out)*/
558
559/* int_forward_safe( tag, m, n, p, X[n][p], Y[m][p]),
560
561nBV = number of Boolean Vectors to be packed
562(see Chapter Dependence Analysis, ADOL-C Documentation)
563bits_per_long = 8*sizeof(unsigned long int)
564p = nBV / bits_per_long + ( (nBV % bits_per_long) != 0 )
565
566The order of the indices in argument and taylors is [var][taylor]
567
568For the full Jacobian matrix set
569p = indep / bits_per_long + ((indep % bits_per_long) != 0)
570and pass a bit pattern version of the identity matrix as an argument    */
571#endif
572#else
573#if defined(_INDOPRO_) && !defined(_NONLIND_OLD_)
574#if defined(_TIGHT_)
575/****************************************************************************/
576/* First Order Vector version of the forward mode for bit patterns, tight   */
577/****************************************************************************/
578int indopro_forward_tight(
579    short             tnum,        /* tape id                              */
580    int               depcheck,    /* consistency chk on # of dependents   */
581    int               indcheck,    /* consistency chk on # of independents */
582    const double     *basepoint,  /* independent variable values   (in)   */
583    unsigned int    **crs)        /* returned row index storage (out)     */
584
585/* indopro_forward_tight( tag, m, n, x[n], *crs[m]),
586
587  */
588
589
590#endif
591#if defined (_NTIGHT_)
592/****************************************************************************/
593/* First Order Vector version of the forward mode, bit pattern, safe        */
594/****************************************************************************/
595int indopro_forward_safe(
596    short             tnum,        /* tape id                              */
597    int               depcheck,    /* consistency chk on # of dependents   */
598    int               indcheck,    /* consistency chk on # of independents */
599    const double     *basepoint,   /* independent variable values   (in)   */
600    unsigned int    **crs)         /* returned row index storage (out)     */
601
602/* indopro_forward_safe( tag, m, n, x[n], *crs[m]),
603
604  */
605#endif
606#else
607#if defined(_NONLIND_)
608#if defined(_TIGHT_)
609/****************************************************************************/
610/* First Order Vector version of the forward mode for bit patterns, tight   */
611/****************************************************************************/
612int nonl_ind_forward_tight(
613    short             tnum,        /* tape id                              */
614    int               depcheck,    /* consistency chk on # of dependents   */
615    int               indcheck,    /* consistency chk on # of independents */
616    const double     *basepoint,  /* independent variable values   (in)   */
617    unsigned int     **crs)        /* returned row index storage (out)     */
618
619#endif
620#if defined (_NTIGHT_)
621/****************************************************************************/
622/* First Order Vector version of the forward mode, bit pattern, safe        */
623/****************************************************************************/
624int nonl_ind_forward_safe(
625    short             tnum,        /* tape id                              */
626    int               depcheck,    /* consistency chk on # of dependents   */
627    int               indcheck,    /* consistency chk on # of independents */
628    const double      *basepoint,  /* independent variable values   (in)   */
629    unsigned int    **crs)        /* returned row index storage (out)     */
630
631/* indopro_forward_safe( tag, m, n, x[n], *crs[m]),
632
633  */
634#endif
635#else
636#if defined(_NONLIND_OLD_)
637#if defined(_TIGHT_)
638/****************************************************************************/
639/* First Order Vector version of the forward mode for bit patterns, tight   */
640/****************************************************************************/
641int nonl_ind_old_forward_tight(
642    short             tnum,        /* tape id                              */
643    int               depcheck,    /* consistency chk on # of dependents   */
644    int               indcheck,    /* consistency chk on # of independents */
645    const double     *basepoint,  /* independent variable values   (in)   */
646    unsigned int     **crs)        /* returned row index storage (out)     */
647
648#endif
649#if defined (_NTIGHT_)
650/****************************************************************************/
651/* First Order Vector version of the forward mode, bit pattern, safe        */
652/****************************************************************************/
653int nonl_ind_old_forward_safe(
654    short             tnum,        /* tape id                              */
655    int               depcheck,    /* consistency chk on # of dependents   */
656    int               indcheck,    /* consistency chk on # of independents */
657    const double      *basepoint,  /* independent variable values   (in)   */
658    unsigned int    **crs)        /* returned row index storage (out)     */
659
660/* indopro_forward_safe( tag, m, n, x[n], *crs[m]),
661
662  */
663#endif
664#else
665#if defined(_FOV_)
666#if defined(_CHUNKED_)
667/****************************************************************************/
668/* First Order Vector version of the forward mode with p-offset in          */
669/* **argument and **taylors                                                 */
670/****************************************************************************/
671int  fov_offset_forward(
672    short  tnum,        /* tape id */
673    int    depcheck,    /* consistency chk on # of deps */
674    int    indcheck,    /* consistency chk on # of indeps */
675    int    p,           /* # of taylor series */
676    int    offset,      /* offset for assignments */
677    double *basepoint,  /* independent variable values */
678    double **argument,  /* Taylor coefficients (input) */
679    double *valuepoint, /* Taylor coefficients (output) */
680    double **taylors)   /* matrix of coifficient vectors */
681/* the order of the indices in argument and taylors is [var][taylor] */
682#else
683/****************************************************************************/
684/* First Order Vector version of the forward mode.                          */
685/****************************************************************************/
686int  fov_forward(
687    short         tnum,        /* tape id */
688    int           depcheck,    /* consistency chk on # of deps */
689    int           indcheck,    /* consistency chk on # of indeps */
690    int           p,           /* # of taylor series */
691    const double *basepoint,   /* independent variable values */
692    double      **argument,    /* Taylor coefficients (input) */
693    double       *valuepoint,  /* Taylor coefficients (output) */
694    double      **taylors)     /* matrix of coifficient vectors */
695/* the order of the indices in argument and taylors is [var][taylor] */
696#endif
697
698#else
699#if defined(_HOS_)
700/****************************************************************************/
701/* Higher Order Scalar version of the forward mode.                         */
702/****************************************************************************/
703#if defined(_KEEP_)
704int  hos_forward(
705#else
706int  hos_forward_nk(
707#endif
708    short  tnum,        /* tape id */
709    int    depcheck,    /* consistency chk on # of dependents */
710    int    indcheck,    /* consistency chk on # of independents */
711    int    gdegree,     /* highest derivative degree */
712#if defined(_KEEP_)
713    int    keep,        /* flag for reverse sweep */
714#endif
715    double *basepoint,  /* independent variable values */
716    double **argument,  /* independant variable values */
717    double *valuepoint, /* Taylor coefficients (output) */
718    double **taylors)   /* matrix of coifficient vectors */
719
720
721#else
722/****************************************************************************/
723/* Higher Order Vector version of the forward mode.                         */
724/****************************************************************************/
725#if defined(_KEEP_)
726int  hov_wk_forward(
727#else
728int  hov_forward(
729#endif
730    short  tnum,        /* tape id */
731    int    depcheck,    /* consistency chk on # of deps */
732    int    indcheck,    /* consistency chk on # of indeps */
733    int    gdegree,     /* highest derivative degree */
734#if defined(_KEEP_)
735    int    keep,        /* flag for reverse sweep */
736#endif
737    int    p,           /* # of taylor series */
738    double *basepoint,  /* independent variable values */
739    double ***argument, /* Taylor coefficients (input) */
740    double *valuepoint, /* Taylor coefficients (output) */
741    double ***taylors)  /* matrix of coifficient vectors */
742/* the order of the indices in argument and taylors is [var][taylor][deriv] */
743
744#endif
745#endif
746#endif
747#endif
748#endif
749#endif
750#endif
751#endif
752{
753    /****************************************************************************/
754    /*                                                            ALL VARIABLES */
755
756    unsigned char operation;   /* operation code */
757    int ret_c =3;              /* return value */
758
759    locint size = 0;
760    locint res  = 0;
761    locint arg  = 0;
762    locint arg1 = 0;
763    locint arg2 = 0;
764
765    double coval = 0, *d = 0;
766
767    int indexi = 0,  indexd = 0;
768
769    /* loop indices */
770#if !defined (_ZOS_)
771#if !defined (_INT_FOR_)
772    int i;
773#if !defined (_INDO_)
774    int ii;
775#endif
776#endif
777#endif
778#if defined (_HIGHER_ORDER_)
779    int j, l=0;
780#endif
781    int ls;
782#if defined(_ADOLC_VECTOR_)
783#if !defined (_HIGHER_ORDER_)
784    int l=0;
785#endif
786#endif
787#if defined (_INT_FOR_)
788    int l=0;
789#endif
790#if defined (_INDO_)
791    int l=0;
792    int max_ind_dom;
793#if defined(_INDOPRO_)
794    /* index domains */
795    locint** ind_dom;
796#endif
797#if defined(_NONLIND_)
798    /* nonlinear interaction domains */
799    locint* arg_index;
800    IndexElement* fod;
801    IndexElement* nonl_dom;
802    IndexElement_sod* sod;
803    IndexElement_sod* temp;
804    IndexElement_sod* temp1;
805    int maxopind;
806    int opind;
807    int ii;
808#endif
809#if defined(_NONLIND_OLD_)
810    /* nonlinear interaction domains */
811    locint** nonl_dom;
812    locint*  temp;
813    locint*  temp1;
814#endif
815#endif
816
817    /* other necessary variables */
818#if !defined (_ZOS_)
819#if !defined (_INDO_)
820#if !defined (_INT_FOR_)
821    double r0=0.0, x, y, divs;
822    int even;
823#endif
824#endif
825#endif
826
827#if defined(_INT_FOR_)
828#ifdef _TIGHT_
829    double  *dp_T0;
830    double y, divs;
831#endif /* _TIGHT_ */
832
833    /* Taylor stuff */
834    unsigned long int  **up_T;
835
836    unsigned long int         *Tres, *Targ, *Targ1, *Targ2;
837#ifdef _TIGHT_
838    unsigned long int         *Tqo;
839    unsigned long int         *Targ1OP, *Targ2OP;
840#endif
841
842#define T0res  T0temp
843#else
844#if defined(_INDO_)
845#ifdef _TIGHT_
846    double  *dp_T0;
847    double  T0temp;
848    double divs;
849#endif /* _TIGHT_ */
850#define T0res  T0temp
851#define T0arg  T0temp
852
853#else
854    double *dp_T0;
855#if !defined(_ZOS_)
856#if  defined(_FOS_)
857    double  *dp_T;
858# define T_TEMP Ttemp;
859# else
860    double *dp_Ttemp, **dpp_T;
861#endif
862    double         *Tres, *Targ, *Targ1, *Targ2, *Tqo;
863
864#if defined (_HIGHER_ORDER_)
865    double         *TresOP, *TresOP2, *zOP;
866    double *dp_z;
867#endif
868   double         *TargOP, *Targ1OP, *Targ2OP;
869   double         T0temp;
870#endif
871#define T0res  T0temp
872#define T0arg  T0temp
873#endif
874#endif
875
876#if defined(_HIGHER_ORDER_)
877    int k = gdegree;
878#endif
879
880#if defined(_KEEP_)
881    int taylbuf=0;
882#endif
883
884#if defined(_HOV_)
885    int pk = k*p;
886#else
887#if defined(_HOV_WK_)
888    int pk = k*p;
889#endif
890#endif
891
892    /* extern diff. function variables */
893#if defined(_EXTERN_)
894#  undef (_EXTERN_)
895#endif
896    /* ZOS_FORWARD */
897#if defined(_ZOS_)
898#   define _EXTERN_ 1
899#   define ADOLC_EXT_FCT_POINTER zos_forward
900#   define ADOLC_EXT_FCT_COMPLETE \
901    zos_forward(n, edfct->dp_x, m, edfct->dp_y)
902#   define ADOLC_EXT_LOOP
903#   define ADOLC_EXT_SUBSCRIPT
904#endif
905    /* FOS_FORWARD */
906#if defined(_FOS_)
907#   define _EXTERN_ 1
908#   define ADOLC_EXT_FCT_POINTER fos_forward
909#   define ADOLC_EXT_FCT_COMPLETE \
910    fos_forward(n, edfct->dp_x, edfct->dp_X, m, edfct->dp_y, edfct->dp_Y)
911#   define ADOLC_EXT_POINTER_X edfct->dp_X
912#   define ADOLC_EXT_POINTER_Y edfct->dp_Y
913#   define ADOLC_EXT_LOOP
914#   define ADOLC_EXT_SUBSCRIPT
915#endif
916    /* FOV_FORWARD */
917#if defined(_FOV_)
918#   define _EXTERN_ 1
919#   define ADOLC_EXT_FCT_POINTER fov_forward
920#   define ADOLC_EXT_FCT_COMPLETE \
921    fov_forward(n, edfct->dp_x,p, edfct->dpp_X, m, edfct->dp_y, edfct->dpp_Y)
922#   define ADOLC_EXT_POINTER_X edfct->dpp_X
923#   define ADOLC_EXT_POINTER_Y edfct->dpp_Y
924#   define ADOLC_EXT_LOOP for (loop2 = 0; loop2 < p; ++loop2)
925#   define ADOLC_EXT_SUBSCRIPT [loop2]
926#endif
927
928#if defined(_EXTERN_)
929    locint n, m;
930    ext_diff_fct *edfct;
931    int loop;
932#   if defined(_FOV_)
933        int loop2;
934#   endif
935    int ext_retc;
936#endif
937
938    ADOLC_OPENMP_THREAD_NUMBER;
939
940#if defined(ADOLC_DEBUG)
941    /****************************************************************************/
942    /*                                                           DEBUG MESSAGES */
943    fprintf(DIAG_OUT,"Call of %s(..) with tag: %d, n: %d, m %d,\n",
944            GENERATED_FILENAME, tnum, indcheck, depcheck);
945#if defined(_KEEP_)
946    fprintf(DIAG_OUT,"                    keep: %d\n", keep);
947#endif
948#if defined(_HIGHER_ORDER_)
949    fprintf(DIAG_OUT,"                    degree: %d\n",gdegree);
950#endif
951#if defined(_ADOLC_VECTOR_)
952    fprintf(DIAG_OUT,"                    p: %d\n\n",p);
953#endif
954
955#endif
956
957    /****************************************************************************/
958    /*                                                                    INITs */
959
960   /* Set up stuff for the tape */
961    ADOLC_OPENMP_GET_THREAD_NUMBER;
962
963    /* Initialize the Forward Sweep */
964
965    init_for_sweep(tnum);
966
967      if ( (depcheck != ADOLC_CURRENT_TAPE_INFOS.stats[NUM_DEPENDENTS]) ||
968            (indcheck != ADOLC_CURRENT_TAPE_INFOS.stats[NUM_INDEPENDENTS]) ) {
969        fprintf(DIAG_OUT,"ADOL-C error: forward sweep on tape %d  aborted!\n"
970                "Number of dependent(%u) and/or independent(%u) variables passed"
971                " to forward is\ninconsistent with number "
972                "recorded on tape (%d, %d) \n", tnum,
973                depcheck, indcheck,
974                ADOLC_CURRENT_TAPE_INFOS.stats[NUM_DEPENDENTS],
975                ADOLC_CURRENT_TAPE_INFOS.stats[NUM_INDEPENDENTS]);
976        exit (-1);
977    }
978
979    /****************************************************************************/
980    /*                                                        MEMORY ALLOCATION */
981    /* olvo 980626 has to be revised for common blocks */
982
983    /*--------------------------------------------------------------------------*/
984#if !defined(_NTIGHT_)
985    dp_T0 = myalloc1(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES]);
986    ADOLC_CURRENT_TAPE_INFOS.dp_T0 = dp_T0;
987#endif /* !_NTIGHT_ */
988#if defined(_ZOS_)                                                   /* ZOS */
989
990#if defined(_KEEP_)
991    if (keep>1) {
992        fprintf(DIAG_OUT,"\n ADOL-C error: zero order scalar forward cannot save"
993                " more\nthan zero order taylor coefficients!\n");
994        exit (-1);
995    }
996#endif
997#if defined(_KEEP_)
998    if (keep) {
999      taylbuf = ADOLC_CURRENT_TAPE_INFOS.stats[TAY_BUFFER_SIZE];
1000
1001        taylor_begin(taylbuf,&dp_T0,keep-1);
1002    }
1003#endif
1004
1005    /*--------------------------------------------------------------------------*/
1006#else                                                                /* FOS */
1007#if defined(_FOS_)
1008#if defined(_KEEP_)
1009    if (keep>2) {
1010        fprintf(DIAG_OUT,"\n ADOL-C error: first order scalar forward cannot save"
1011                " more  \nthan first order taylor coefficients!\n");
1012        exit (-1);
1013    }
1014#endif
1015    dp_T = myalloc1(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES]);
1016# define TAYLOR_BUFFER dp_T
1017#if defined(_KEEP_)
1018    if (keep) {
1019        taylbuf = ADOLC_CURRENT_TAPE_INFOS.stats[TAY_BUFFER_SIZE];
1020        taylor_begin(taylbuf,&dp_T,keep-1);
1021    }
1022#endif
1023
1024    /*--------------------------------------------------------------------------*/
1025#else                                                                /* INF_FOR */
1026#if defined(_INT_FOR_)
1027        up_T     = myalloc2_ulong(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES],p);
1028#define TAYLOR_BUFFER up_T
1029
1030    /*--------------------------------------------------------------------------*/
1031#else                                                                /* INDOPRO */
1032#if defined(_INDO_)
1033#if defined(_INDOPRO_)
1034    /* index domains */
1035    ind_dom = (locint **)  malloc(sizeof(locint*) * ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES]);
1036    max_ind_dom = ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES];
1037
1038    for(i=0;i<max_ind_dom;i++)
1039    {
1040        ind_dom[i] = (locint *)  malloc(sizeof(locint) * (NUMNNZ+2));
1041        ind_dom[i][0] = 0;
1042        ind_dom[i][1] = NUMNNZ;
1043    }
1044#endif
1045#if defined(_NONLIND_)
1046        maxopind=ADOLC_CURRENT_TAPE_INFOS.stats[NUM_OPERATIONS]+ADOLC_CURRENT_TAPE_INFOS.stats[NUM_EQ_PROD];
1047        fod = (struct IndexElement*) malloc(sizeof(struct IndexElement) * maxopind);
1048        nonl_dom = (struct IndexElement*) malloc(sizeof(struct IndexElement) * indcheck);
1049        arg_index = (locint *)  malloc(sizeof(locint) * (ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES]));
1050        opind = 0;
1051        sod = (struct IndexElement_sod*) malloc(sizeof(struct IndexElement_sod) * indcheck);
1052
1053        for(i=0;i<indcheck;i++)
1054          {
1055            sod[i].entry = 0;
1056            sod[i].left = NULL;
1057            nonl_dom[i].entry = indcheck+2;
1058            nonl_dom[i].left = NULL;
1059            nonl_dom[i].right = NULL;
1060          }
1061#endif
1062#if defined(_NONLIND_OLD_)
1063
1064    nonl_dom = (locint**) malloc(sizeof(locint*) * indcheck);
1065    for(i=0;i<indcheck;i++){
1066          nonl_dom[i] = (locint*) malloc(sizeof(locint)*(NUMNNZ+2));
1067          nonl_dom[i][0]=0;
1068          nonl_dom[i][1]=NUMNNZ;
1069       }
1070#endif
1071
1072    /*--------------------------------------------------------------------------*/
1073#else                                                                /* FOV */
1074#if defined(_FOV_)
1075    dpp_T = myalloc2(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES],p);
1076# define TAYLOR_BUFFER dpp_T
1077    dp_Ttemp = myalloc1(p);
1078# define T_TEMP dp_Ttemp;
1079
1080    /*--------------------------------------------------------------------------*/
1081#else                                                                /* HOS */
1082#if defined(_HOS_)
1083    dpp_T = myalloc2(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES],k);
1084# define TAYLOR_BUFFER dpp_T
1085    dp_z  = myalloc1(k);
1086    dp_Ttemp = myalloc1(k);
1087# define T_TEMP dp_Ttemp;
1088#if defined(_KEEP_)
1089    if (keep) {
1090        taylbuf = ADOLC_CURRENT_TAPE_INFOS.stats[TAY_BUFFER_SIZE];
1091        taylor_begin(taylbuf,dpp_T,keep-1);
1092    }
1093#endif
1094
1095    /*--------------------------------------------------------------------------*/
1096#else                                                     /* HOV and HOV_WK */
1097    dpp_T = myalloc2(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES],p*k);
1098# define TAYLOR_BUFFER dpp_T
1099    dp_z  = myalloc1(k);
1100    dp_Ttemp = myalloc1(p*k);
1101# define T_TEMP dp_Ttemp;
1102#if defined(_KEEP_)
1103    if (keep) {
1104        taylbuf = ADOLC_CURRENT_TAPE_INFOS.stats[TAY_BUFFER_SIZE];
1105        taylor_begin(taylbuf,dpp_T,keep-1);
1106    }
1107#endif
1108#endif
1109#endif
1110#endif
1111#endif
1112#endif
1113#endif
1114    /****************************************************************************/
1115    /*                                                            FORWARD SWEEP */
1116
1117#if defined(ADOLC_DEBUG)
1118/* #include <string.h> */
1119    int v = 0;
1120    unsigned int countPerOperation[256], taylorPerOperation[256];
1121    memset(countPerOperation, 0, 1024);
1122    memset(taylorPerOperation, 0, 1024);
1123#   define UPDATE_TAYLORWRITTEN(X) taylorPerOperation[operation] += X;
1124#else
1125#   define UPDATE_TAYLORWRITTEN(X)
1126#endif /* ADOLC_DEBUG */
1127
1128    operation=get_op_f();
1129#if defined(ADOLC_DEBUG)
1130    ++countPerOperation[operation];
1131#endif /* ADOLC_DEBUG */
1132
1133    while (operation !=end_of_tape) {
1134     
1135      switch (operation) {
1136
1137
1138                /****************************************************************************/
1139                /*                                                                  MARKERS */
1140
1141                /*--------------------------------------------------------------------------*/
1142            case end_of_op:                                          /* end_of_op */
1143                get_op_block_f();
1144                operation=get_op_f();
1145                /* Skip next operation, it's another end_of_op */
1146                break;
1147
1148                /*--------------------------------------------------------------------------*/
1149            case end_of_int:                                        /* end_of_int */
1150                get_loc_block_f();
1151                break;
1152
1153                /*--------------------------------------------------------------------------*/
1154            case end_of_val:                                        /* end_of_val */
1155               get_val_block_f();
1156                break;
1157                /*--------------------------------------------------------------------------*/
1158            case start_of_tape:                                  /* start_of_tape */
1159            case end_of_tape:                                      /* end_of_tape */
1160                break;
1161
1162
1163                /****************************************************************************/
1164                /*                                                               COMPARISON */
1165
1166                /*--------------------------------------------------------------------------*/
1167            case eq_zero:                                              /* eq_zero */
1168                arg = get_locint_f();
1169
1170#if !defined(_NTIGHT_)
1171                if (dp_T0[arg] != 0) {
1172                    if (ADOLC_GLOBAL_TAPE_VARS.branchSwitchWarning)
1173                        fprintf(DIAG_OUT,
1174                                "ADOL-C Warning: Branch switch detected in comparison "
1175                                "(operator eq_zero).\n"
1176                                "Forward sweep aborted! Retaping recommended!\n");
1177                    ret_c = -1;
1178                    operation = end_of_tape;
1179                    continue;
1180                }
1181                ret_c = 0;
1182#endif /* !_NTIGHT_ */
1183                break;
1184
1185                /*--------------------------------------------------------------------------*/
1186            case neq_zero:                                            /* neq_zero */
1187                arg = get_locint_f();
1188
1189#if !defined(_NTIGHT_)
1190                if (dp_T0[arg] == 0) {
1191                    if (ADOLC_GLOBAL_TAPE_VARS.branchSwitchWarning)
1192                        fprintf(DIAG_OUT,
1193                                "ADOL-C Warning: Branch switch detected in comparison "
1194                                "(operator neq_zero).\n"
1195                                "Forward sweep aborted! Retaping recommended!\n");
1196                    ret_c = -1;
1197                    operation = end_of_tape;
1198                    continue;
1199                }
1200#endif /* !_NTIGHT_ */
1201                break;
1202
1203                /*--------------------------------------------------------------------------*/
1204            case le_zero:                                              /* le_zero */
1205                arg = get_locint_f();
1206
1207#if !defined(_NTIGHT_)
1208                if (dp_T0[arg] > 0) {
1209                    if (ADOLC_GLOBAL_TAPE_VARS.branchSwitchWarning)
1210                        fprintf(DIAG_OUT,
1211                                "ADOL-C Warning: Branch switch detected in comparison "
1212                                "(operator le_zero).\n"
1213                                "Forward sweep aborted! Retaping recommended!\n");
1214                    ret_c = -1;
1215                    operation = end_of_tape;
1216                    continue;
1217                }
1218                if (dp_T0[arg] == 0)
1219                    ret_c = 0;
1220#endif /* !_NTIGHT_ */
1221                break;
1222
1223                /*--------------------------------------------------------------------------*/
1224            case gt_zero:                                              /* gt_zero */
1225                arg = get_locint_f();
1226
1227#if !defined(_NTIGHT_)
1228                if (dp_T0[arg] <= 0) {
1229                    if (ADOLC_GLOBAL_TAPE_VARS.branchSwitchWarning)
1230                        fprintf(DIAG_OUT,
1231                                "ADOL-C Warning: Branch switch detected in comparison "
1232                                "(operator gt_zero).\n"
1233                                "Forward sweep aborted! Retaping recommended!\n");
1234                    ret_c = -1;
1235                    operation = end_of_tape;
1236                    continue;
1237                }
1238#endif /* !_NTIGHT_ */
1239                break;
1240
1241                /*--------------------------------------------------------------------------*/
1242            case ge_zero:                                              /* ge_zero */
1243                arg = get_locint_f();
1244
1245#if !defined(_NTIGHT_)
1246                if (dp_T0[arg] < 0) {
1247                    if (ADOLC_GLOBAL_TAPE_VARS.branchSwitchWarning)
1248                        fprintf(DIAG_OUT,
1249                                "ADOL-C Warning: Branch switch detected in comparison "
1250                                "(operator ge_zero).\n"
1251                                "Forward sweep aborted! Retaping recommended!\n");
1252                    ret_c = -1;
1253                    operation = end_of_tape;
1254                    continue;
1255                }
1256                if (dp_T0[arg] == 0)
1257                    ret_c = 0;
1258#endif /* !_NTIGHT_ */
1259                break;
1260
1261                /*--------------------------------------------------------------------------*/
1262            case lt_zero:                                              /* lt_zero */
1263                arg = get_locint_f();
1264
1265#if !defined(_NTIGHT_)
1266                if (dp_T0[arg] >= 0) {
1267                    if (ADOLC_GLOBAL_TAPE_VARS.branchSwitchWarning)
1268                        fprintf(DIAG_OUT,
1269                                "ADOL-C Warning: Branch switch detected in comparison "
1270                                "(operator lt_zero).\n"
1271                                "Forward sweep aborted! Retaping recommended!\n");
1272                    ret_c = -1;
1273                    operation = end_of_tape;
1274                    continue;
1275                }
1276#endif /* !_NTIGHT_ */
1277                break;
1278
1279
1280                /****************************************************************************/
1281                /*                                                              ASSIGNMENTS */
1282
1283                /*--------------------------------------------------------------------------*/
1284            case assign_a:           /* assign an adouble variable an    assign_a */
1285                /* adouble value. (=) */
1286                arg = get_locint_f();
1287                res = get_locint_f();
1288
1289                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1290
1291
1292#if !defined(_NTIGHT_)
1293                dp_T0[res] = dp_T0[arg];
1294#endif /* !_NTIGHT_ */
1295
1296#if defined(_INDO_)
1297#if defined(_INDOPRO_)
1298                copy_index_domain(res, arg, ind_dom);
1299#endif
1300#if defined(_NONLIND_)
1301                arg_index[res] = arg_index[arg];               
1302#endif           
1303#else
1304#if !defined(_ZOS_) /* BREAK_ZOS */
1305                ASSIGN_T(Targ,TAYLOR_BUFFER[arg])
1306                ASSIGN_T(Tres,TAYLOR_BUFFER[res])
1307
1308                FOR_0_LE_l_LT_pk
1309                TRES_INC = TARG_INC;
1310#endif
1311#endif /* ALL_TOGETHER_AGAIN */
1312                break;
1313
1314                /*--------------------------------------------------------------------------*/
1315            case assign_d:            /* assign an adouble variable a    assign_d */
1316                /* double value. (=) */
1317                res   = get_locint_f();
1318                coval = get_val_f();
1319
1320                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1321
1322#if !defined(_NTIGHT_)
1323                dp_T0[res] = coval;
1324#endif /* !_NTIGHT_ */
1325
1326#if defined(_INDO_)
1327#if defined(_INDOPRO_)
1328                ind_dom[res][0]=0;
1329#endif
1330#if defined(_NONLIND_)
1331                fod[opind].entry = maxopind+2;
1332                fod[opind].left = NULL;
1333                fod[opind].right = NULL;
1334                arg_index[res] = opind++;               
1335#endif
1336#else
1337#if !defined(_ZOS_) /* BREAK_ZOS */
1338                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
1339
1340                FOR_0_LE_l_LT_pk
1341                TRES_INC = 0;
1342#endif
1343#endif /* ALL_TOGETHER_AGAIN */
1344                break;
1345
1346                /*--------------------------------------------------------------------------*/
1347            case assign_d_zero:  /* assign an adouble variable a    assign_d_zero */
1348                /* double value. (0) (=) */
1349                res   = get_locint_f();
1350
1351                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1352
1353#if !defined(_NTIGHT_)
1354                dp_T0[res] = 0.0;
1355#endif /* !_NTIGHT_ */
1356
1357#if defined(_INDO_)
1358#if defined(_INDOPRO_)
1359                ind_dom[res][0]=0;
1360#endif
1361#if defined(_NONLIND_)
1362                fod[opind].entry = maxopind+2;
1363                fod[opind].left = NULL;
1364                fod[opind].right = NULL;
1365                arg_index[res] = opind++;               
1366#endif
1367#else
1368#if !defined(_ZOS_) /* BREAK_ZOS */
1369                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
1370
1371                FOR_0_LE_l_LT_pk
1372                TRES_INC = 0;
1373#endif
1374#endif /* ALL_TOGETHER_AGAIN */
1375                break;
1376
1377                /*--------------------------------------------------------------------------*/
1378            case assign_d_one:    /* assign an adouble variable a    assign_d_one */
1379                /* double value. (1) (=) */
1380                res   = get_locint_f();
1381
1382                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1383
1384#if !defined(_NTIGHT_)
1385                dp_T0[res] = 1.0;
1386#endif /* !_NTIGHT_ */
1387
1388#if defined(_INDO_)
1389#if defined(_INDOPRO_)
1390                ind_dom[res][0]=0;
1391#endif
1392#if defined(_NONLIND_)
1393                fod[opind].entry = maxopind+2;
1394                fod[opind].left = NULL;
1395                fod[opind].right = NULL;
1396                arg_index[res] = opind++;               
1397#endif
1398#else
1399#if !defined(_ZOS_) /* BREAK_ZOS */
1400                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
1401
1402                FOR_0_LE_l_LT_pk
1403                TRES_INC = 0;
1404
1405#endif
1406#endif /* ALL_TOGETHER_AGAIN */
1407                break;
1408
1409                /*--------------------------------------------------------------------------*/
1410            case assign_ind:       /* assign an adouble variable an    assign_ind */
1411                /* independent double value (<<=) */
1412                res = get_locint_f();
1413
1414                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1415
1416#if !defined(_NTIGHT_)
1417                dp_T0[res] = basepoint[indexi];
1418#endif /* !_NTIGHT_ */
1419
1420#if defined(_INDO_)
1421#if defined(_INDOPRO_)
1422                ind_dom[res][0] = 1;
1423                ind_dom[res][2] = indexi;
1424#endif         
1425#if defined(_NONLIND_)
1426                fod[opind].entry = indexi;
1427                fod[opind].left = NULL;
1428                fod[opind].right = NULL;
1429                arg_index[res] = opind++;               
1430#endif         
1431#else
1432#if !defined(_ZOS_) /* BREAK_ZOS */
1433                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
1434
1435#ifdef _INT_FOR_
1436                FOR_0_LE_l_LT_p
1437                TRES_INC = ARGUMENT(indexi,l,i);
1438#else
1439                FOR_0_LE_l_LT_p
1440                FOR_0_LE_i_LT_k
1441                TRES_INC = ARGUMENT(indexi,l,i);
1442#endif
1443#endif
1444#endif /* ALL_TOGETHER_AGAIN */
1445                ++indexi;
1446                break;
1447
1448                /*--------------------------------------------------------------------------*/
1449            case assign_dep:           /* assign a float variable a    assign_dep */
1450                /* dependent adouble value. (>>=) */
1451                res = get_locint_f();
1452
1453#if !defined(_INDO_)
1454#if !defined(_NTIGHT_)
1455                if ( valuepoint != NULL )
1456                  valuepoint[indexd] = dp_T0[res];
1457#endif /* !_NTIGHT_ */
1458#endif
1459
1460#if defined(_INDO_)
1461#if defined(_INDOPRO_) && !defined(_NONLIND_OLD_)
1462          if (ind_dom[res][0] != 0) {
1463            crs[indexd] = (unsigned int*) malloc(sizeof(unsigned int) * (ind_dom[res][0]+1));
1464            crs[indexd][0] = ind_dom[res][0];
1465            for(l=1;l<=crs[indexd][0];l++) {
1466              crs[indexd][l] = ind_dom[res][l+1];
1467            }
1468          }
1469          else {
1470            crs[indexd] = (unsigned int*) malloc(sizeof(unsigned int));
1471            crs[indexd][0] =0;
1472          }
1473#endif
1474#else
1475#if !defined(_ZOS_) /* BREAK_ZOS */
1476                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
1477
1478#ifdef _INT_FOR_
1479                if (taylors != 0 )  /* ??? question: why here? */
1480                    FOR_0_LE_l_LT_p
1481                    TAYLORS(indexd,l,i) = TRES_INC;
1482#else
1483                if (taylors != 0 )  /* ??? question: why here? */
1484                    FOR_0_LE_l_LT_p
1485                    FOR_0_LE_i_LT_k
1486                    TAYLORS(indexd,l,i) = TRES_INC;
1487#endif
1488#endif
1489#endif /* ALL_TOGETHER_AGAIN */
1490                indexd++;
1491                break;
1492
1493
1494                /****************************************************************************/
1495                /*                                                   OPERATION + ASSIGNMENT */
1496
1497                /*--------------------------------------------------------------------------*/
1498            case eq_plus_d:            /* Add a floating point to an    eq_plus_d */
1499                /* adouble. (+=) */
1500                res   = get_locint_f();
1501                coval = get_val_f();
1502
1503                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1504
1505#if !defined(_NTIGHT_)
1506                dp_T0[res] += coval;
1507#endif /* !_NTIGHT_ */
1508                break;
1509
1510                /*--------------------------------------------------------------------------*/
1511            case eq_plus_a:             /* Add an adouble to another    eq_plus_a */
1512                /* adouble. (+=) */
1513                arg = get_locint_f();
1514                res = get_locint_f();
1515
1516                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1517
1518#if !defined(_NTIGHT_)
1519                dp_T0[res] += dp_T0[arg];
1520#endif /* !_NTIGHT_ */
1521
1522#if defined(_INDO_)
1523#if defined(_INDOPRO_)
1524                merge_2_index_domains(res, arg, ind_dom);
1525#endif
1526#if defined(_NONLIND_)
1527                fod[opind].entry = maxopind+2;
1528                fod[opind].left = &fod[arg_index[res]];
1529                fod[opind].right = &fod[arg_index[arg]];
1530                arg_index[res] = opind++;               
1531#endif
1532#else
1533#if !defined(_ZOS_) /* BREAK_ZOS */
1534                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
1535                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
1536
1537#ifdef _INT_FOR_
1538                FOR_0_LE_l_LT_pk
1539                TRES_INC |= TARG_INC;
1540#else
1541                FOR_0_LE_l_LT_pk
1542                TRES_INC += TARG_INC;
1543#endif
1544#endif
1545#endif /* ALL_TOGETHER_AGAIN */
1546                break;
1547
1548                /*--------------------------------------------------------------------------*/
1549            case eq_min_d:       /* Subtract a floating point from an    eq_min_d */
1550                /* adouble. (-=) */
1551                res = get_locint_f();
1552                coval = get_val_f();
1553
1554                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1555
1556#if !defined(_NTIGHT_)
1557                dp_T0[res] -= coval;
1558#endif /* !_NTIGHT_ */
1559                break;
1560
1561                /*--------------------------------------------------------------------------*/
1562            case eq_min_a:        /* Subtract an adouble from another    eq_min_a */
1563                /* adouble. (-=) */
1564                arg = get_locint_f();
1565                res = get_locint_f();
1566
1567                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1568
1569#if !defined(_NTIGHT_)
1570                dp_T0[res] -= dp_T0[arg];
1571#endif /* !_NTIGHT_ */
1572
1573#if defined(_INDO_)
1574#if defined(_INDOPRO_)
1575                merge_2_index_domains(res, arg, ind_dom);
1576#endif
1577#if defined(_NONLIND_)
1578                fod[opind].entry = maxopind+2;
1579                fod[opind].left = &fod[arg_index[res]];
1580                fod[opind].right = &fod[arg_index[arg]];
1581                arg_index[res] = opind++;               
1582#endif
1583#else
1584#if !defined(_ZOS_) /* BREAK_ZOS */
1585                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
1586                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
1587
1588#ifdef _INT_FOR_
1589                FOR_0_LE_l_LT_pk
1590                TRES_INC |= TARG_INC;
1591#else
1592                FOR_0_LE_l_LT_pk
1593                TRES_INC -= TARG_INC;
1594#endif
1595#endif
1596#endif /* ALL_TOGETHER_AGAIN */
1597                break;
1598
1599                /*--------------------------------------------------------------------------*/
1600            case eq_mult_d:              /* Multiply an adouble by a    eq_mult_d */
1601                /* flaoting point. (*=) */
1602                res   = get_locint_f();
1603                coval = get_val_f();
1604
1605                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1606
1607#if !defined(_NTIGHT_)
1608                dp_T0[res] *= coval;
1609#endif /* !_NTIGHT_ */
1610
1611#if !defined(_INDO_)
1612#if !defined(_ZOS_) /* BREAK_ZOS */
1613#if !defined( _INT_FOR_)
1614
1615                FOR_0_LE_l_LT_pk
1616                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
1617
1618                FOR_0_LE_l_LT_pk
1619                TRES_INC *= coval;
1620#endif
1621#endif
1622#endif /* ALL_TOGETHER_AGAIN */
1623                break;
1624
1625                /*--------------------------------------------------------------------------*/
1626            case eq_mult_a:       /* Multiply one adouble by another    eq_mult_a */
1627                /* (*=) */
1628                arg = get_locint_f();
1629                res = get_locint_f();
1630
1631                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1632
1633#if defined(_INDO_)
1634#if defined(_INDOPRO_)
1635                merge_2_index_domains(res, arg, ind_dom);
1636#endif
1637#if defined(_NONLIND_)
1638                fod[opind].entry = maxopind+2;
1639                fod[opind].left = &fod[arg_index[res]];
1640                fod[opind].right = &fod[arg_index[arg]];
1641                traverse_unary(&fod[arg_index[res]], nonl_dom, &fod[arg_index[arg]], indcheck+1,maxopind+2);
1642                traverse_unary(&fod[arg_index[arg]], nonl_dom, &fod[arg_index[res]], indcheck+1,maxopind+2);
1643                arg_index[res] = opind++;               
1644#endif
1645#if defined(_NONLIND_OLD_)
1646                extend_nonlinearity_domain_binary(res, arg, ind_dom, nonl_dom);
1647#endif
1648#else
1649#if !defined(_ZOS_) /* BREAK_ZOS */
1650                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
1651                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
1652
1653                INC_pk_1(Tres)
1654                INC_pk_1(Targ)
1655
1656#ifdef _INT_FOR_
1657                FOR_p_GT_l_GE_0
1658                TRES_FODEC |= TARG_DEC;
1659#else
1660                FOR_p_GT_l_GE_0
1661                FOR_k_GT_i_GE_0
1662                { TRES_FODEC = dp_T0[res]*TARG_DEC +
1663                               TRES*dp_T0[arg];
1664                  DEC_TRES_FO
1665#ifdef _HIGHER_ORDER_
1666                  TresOP = Tres-i;
1667                  TargOP = Targ;
1668
1669                  for (j=0;j<i;j++)
1670                  *Tres += (*TresOP++) * (*TargOP--);
1671                  Tres--;
1672#endif /* _HIGHER_ORDER_ */
1673                }
1674#endif
1675#endif
1676#endif /* ALL_TOGETHER_AGAIN */
1677#if !defined(_NTIGHT_)
1678               dp_T0[res] *= dp_T0[arg];
1679#endif /* !_NTIGHT_ */
1680                break;
1681
1682                /*--------------------------------------------------------------------------*/
1683            case incr_a:                        /* Increment an adouble    incr_a */
1684                res   = get_locint_f();
1685
1686                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1687
1688#if !defined(_NTIGHT_)
1689                dp_T0[res]++;
1690#endif /* !_NTIGHT_ */
1691                break;
1692
1693                /*--------------------------------------------------------------------------*/
1694            case decr_a:                        /* Increment an adouble    decr_a */
1695                res   = get_locint_f();
1696
1697                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1698
1699#if !defined(_NTIGHT_)
1700                dp_T0[res]--;
1701#endif /* !_NTIGHT_ */
1702                break;
1703
1704
1705                /****************************************************************************/
1706                /*                                                        BINARY OPERATIONS */
1707
1708                /*--------------------------------------------------------------------------*/
1709            case plus_a_a:                 /* : Add two adoubles. (+)    plus a_a */
1710                arg1 = get_locint_f();
1711                arg2 = get_locint_f();
1712                res  = get_locint_f();
1713
1714                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1715
1716#if !defined(_NTIGHT_)
1717                dp_T0[res] = dp_T0[arg1] + dp_T0[arg2];
1718#endif /* !_NTIGHT_ */
1719
1720#if defined(_INDO_)
1721#if defined(_INDOPRO_)
1722                combine_2_index_domains(res, arg1, arg2, ind_dom);
1723#endif     
1724#if defined(_NONLIND_)
1725                fod[opind].entry = maxopind+2;
1726                fod[opind].left = &fod[arg_index[arg1]];
1727                fod[opind].right = &fod[arg_index[arg2]];
1728                arg_index[res] = opind++;               
1729#endif
1730#else
1731#if !defined(_ZOS_) /* BREAK_ZOS */
1732                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
1733                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
1734                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
1735
1736#ifdef _INT_FOR_
1737                FOR_0_LE_l_LT_pk
1738                TRES_INC = TARG1_INC | TARG2_INC;
1739#else
1740                FOR_0_LE_l_LT_pk
1741                TRES_INC = TARG1_INC + TARG2_INC;
1742#endif
1743#endif
1744#endif /* ALL_TOGETHER_AGAIN */
1745                break;
1746
1747                /*--------------------------------------------------------------------------*/
1748            case plus_d_a:             /* Add an adouble and a double    plus_d_a */
1749                /* (+) */
1750                arg   = get_locint_f();
1751                res   = get_locint_f();
1752                coval = get_val_f();
1753
1754                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1755
1756#if !defined(_NTIGHT_)
1757                dp_T0[res] = dp_T0[arg] + coval;
1758#endif /* !_NTIGHT_ */
1759
1760#if defined(_INDO_)
1761#if defined(_INDOPRO_)
1762                copy_index_domain(res, arg, ind_dom);
1763#endif
1764#if defined(_NONLIND_)
1765                arg_index[res] = arg_index[arg];
1766#endif               
1767#else
1768#if !defined(_ZOS_) /* BREAK_ZOS */
1769                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
1770                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
1771
1772                FOR_0_LE_l_LT_pk
1773                TRES_INC = TARG_INC;
1774#endif
1775#endif /* ALL_TOGETHER_AGAIN */
1776                break;
1777
1778                /*--------------------------------------------------------------------------*/
1779            case min_a_a:              /* Subtraction of two adoubles     min_a_a */
1780                /* (-) */
1781                arg1 = get_locint_f();
1782                arg2 = get_locint_f();
1783                res  = get_locint_f();
1784
1785                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1786
1787#if !defined(_NTIGHT_)
1788                dp_T0[res] = dp_T0[arg1] -
1789                                               dp_T0[arg2];
1790#endif /* !_NTIGHT_ */
1791
1792
1793#if defined(_INDO_)   
1794#if defined(_INDOPRO_)
1795                combine_2_index_domains(res, arg1, arg2, ind_dom);
1796#endif
1797#if defined(_NONLIND_)
1798                fod[opind].entry = maxopind+2;
1799                fod[opind].left = &fod[arg_index[arg1]];
1800                fod[opind].right = &fod[arg_index[arg2]];
1801                arg_index[res] = opind++;               
1802#endif
1803#else
1804#if !defined(_ZOS_) /* BREAK_ZOS */
1805                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
1806                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
1807                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
1808
1809#ifdef _INT_FOR_
1810                FOR_0_LE_l_LT_pk
1811                TRES_INC = TARG1_INC | TARG2_INC;
1812#else
1813                 FOR_0_LE_l_LT_pk
1814                TRES_INC = TARG1_INC - TARG2_INC;
1815#endif
1816#endif
1817#endif /* ALL_TOGETHER_AGAIN */
1818                break;
1819
1820                /*--------------------------------------------------------------------------*/
1821            case min_d_a:                /* Subtract an adouble from a    min_d_a */
1822                /* double (-) */
1823                arg =get_locint_f();
1824                res = get_locint_f();
1825                coval = get_val_f();
1826
1827                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1828
1829#if !defined(_NTIGHT_)
1830                dp_T0[res] = coval - dp_T0[arg];
1831#endif /* !_NTIGHT_ */
1832
1833#if defined(_INDO_)
1834#if defined(_INDOPRO_)
1835                copy_index_domain(res, arg, ind_dom);
1836#endif
1837#if defined(_NONLIND_)
1838                arg_index[res] = arg_index[arg];               
1839#endif
1840#else
1841#if !defined(_ZOS_) /* BREAK_ZOS */
1842                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
1843                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
1844
1845#ifdef _INT_FOR_
1846                FOR_0_LE_l_LT_pk
1847                TRES_INC = TARG_INC;
1848#else
1849                FOR_0_LE_l_LT_pk
1850                TRES_INC = -TARG_INC;
1851#endif
1852#endif
1853#endif /* ALL_TOGETHER_AGAIN */
1854                break;
1855
1856                /*--------------------------------------------------------------------------*/
1857            case mult_a_a:               /* Multiply two adoubles (*)    mult_a_a */
1858                arg1 = get_locint_f();
1859                arg2 = get_locint_f();
1860                res  = get_locint_f();
1861
1862                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
1863
1864#if defined(_INDO_)
1865#if defined(_INDOPRO_)
1866                combine_2_index_domains(res, arg1, arg2, ind_dom);
1867#endif
1868#if defined(_NONLIND_)
1869                fod[opind].entry = maxopind+2;
1870                fod[opind].left = &fod[arg_index[arg1]];
1871                fod[opind].right = &fod[arg_index[arg2]];
1872                traverse_unary(&fod[arg_index[arg1]], nonl_dom, &fod[arg_index[arg2]], indcheck+1,maxopind+2);
1873                traverse_unary(&fod[arg_index[arg2]], nonl_dom, &fod[arg_index[arg1]], indcheck+1,maxopind+2);
1874                arg_index[res] = opind++;               
1875#endif
1876#if defined(_NONLIND_OLD_)
1877                extend_nonlinearity_domain_binary(arg1, arg2, ind_dom, nonl_dom);
1878#endif
1879#else
1880#if !defined(_ZOS_) /* BREAK_ZOS */
1881                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
1882                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
1883                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
1884
1885#ifdef _INT_FOR_
1886                FOR_0_LE_l_LT_p
1887                TRES_FOINC = TARG2_INC | TARG1_INC;
1888#else
1889                /* olvo 980915 now in reverse order to allow x = x*x etc. */
1890                INC_pk_1(Tres)
1891                INC_pk_1(Targ1)
1892                INC_pk_1(Targ2)
1893
1894                FOR_p_GT_l_GE_0
1895                FOR_k_GT_i_GE_0
1896                { TRES_FODEC = dp_T0[arg1]*TARG2_DEC +
1897                               TARG1_DEC*dp_T0[arg2];
1898                  DEC_TRES_FO
1899#if defined(_HIGHER_ORDER_)
1900                  Targ1OP = Targ1-i+1;
1901                  Targ2OP = Targ2;
1902
1903                  for (j=0;j<i;j++) {
1904                  *Tres += (*Targ1OP++) * (*Targ2OP--);
1905                  }
1906                  Tres--;
1907#endif /* _HIGHER_ORDER_ */
1908            }
1909#endif
1910#endif
1911#endif /* ALL_TOGETHER_AGAIN */
1912#if !defined(_NTIGHT_)
1913                dp_T0[res] = dp_T0[arg1] *
1914                                               dp_T0[arg2];
1915#endif /* !_NTIGHT_ */
1916                break;
1917
1918                /*--------------------------------------------------------------------------*/
1919                /* olvo 991122: new op_code with recomputation */
1920            case eq_plus_prod:   /* increment a product of           eq_plus_prod */
1921                /* two adoubles (*) */
1922                arg1 = get_locint_f();
1923                arg2 = get_locint_f();
1924                res  = get_locint_f();
1925
1926#if defined(_INDO_)
1927#if defined(_INDOPRO_)
1928                merge_3_index_domains(res, arg1, arg2, ind_dom);
1929#endif
1930#if defined(_NONLIND_)
1931                // operation: v = v+u*w
1932                // first step: z = u*w, index domains
1933                fod[opind].entry = maxopind+2;
1934                fod[opind].left = &fod[arg_index[arg1]];
1935                fod[opind].right = &fod[arg_index[arg2]];
1936                // first step: z = u*w,
1937                traverse_unary(&fod[arg_index[arg1]], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
1938                traverse_unary(&fod[arg_index[arg2]], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
1939                opind++;
1940                // second step: v = v+z, index domains
1941                fod[opind].entry = maxopind+2;
1942                fod[opind].left = &fod[arg_index[res]];
1943                fod[opind].right = &fod[opind-1];
1944                // second step: v = v+z,
1945                arg_index[res] = opind++;               
1946#endif
1947#if defined(_NONLIND_OLD_)
1948                extend_nonlinearity_domain_binary(arg1, arg2, ind_dom, nonl_dom);
1949#endif
1950#else
1951#if !defined(_ZOS_) /* BREAK_ZOS */
1952                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
1953                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
1954                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
1955
1956#ifdef _INT_FOR_
1957                FOR_0_LE_l_LT_p
1958                TRES_FOINC |= TARG2_INC | TARG1_INC;
1959#else
1960                /* olvo 980915 now in reverse order to allow x = x*x etc. */
1961                INC_pk_1(Tres)
1962                INC_pk_1(Targ1)
1963                INC_pk_1(Targ2)
1964
1965                FOR_p_GT_l_GE_0
1966                FOR_k_GT_i_GE_0
1967                { TRES_FODEC += dp_T0[arg1]*TARG2_DEC +
1968                                TARG1_DEC*dp_T0[arg2];
1969                  DEC_TRES_FO
1970#if defined(_HIGHER_ORDER_)
1971                  Targ1OP = Targ1-i+1;
1972                  Targ2OP = Targ2;
1973
1974                  for (j=0;j<i;j++)
1975                  *Tres += (*Targ1OP++) * (*Targ2OP--);
1976                  Tres--;
1977#endif /* _HIGHER_ORDER_ */
1978                }
1979#endif
1980#endif
1981#endif /* ALL_TOGETHER_AGAIN */
1982#if !defined(_NTIGHT_)
1983                dp_T0[res] += dp_T0[arg1] *
1984                                                    dp_T0[arg2];
1985#endif /* !_NTIGHT_ */
1986                break;
1987
1988                /*--------------------------------------------------------------------------*/
1989                /* olvo 991122: new op_code with recomputation */
1990            case eq_min_prod:    /* decrement a product of            eq_min_prod */
1991                /* two adoubles (*) */
1992                arg1 = get_locint_f();
1993                arg2 = get_locint_f();
1994                res  = get_locint_f();
1995
1996#if defined(_INDO_)
1997#if defined(_INDOPRO_)
1998                merge_3_index_domains(res, arg1, arg2, ind_dom);
1999#endif
2000#if defined(_NONLIND_)
2001                // operation: v = v-u*w
2002                // first step: z = u*w, index domains
2003                fod[opind].entry = maxopind+2;
2004                fod[opind].left = &fod[arg_index[arg1]];
2005                fod[opind].right = &fod[arg_index[arg2]];
2006                // first step: z = u*w,
2007                traverse_unary(&fod[arg_index[arg1]], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2008                traverse_unary(&fod[arg_index[arg2]], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2009                opind++;
2010                // second step: v = v-z, index domains
2011                fod[opind].entry = maxopind+2;
2012                fod[opind].left = &fod[arg_index[res]];
2013                fod[opind].right = &fod[opind-1];
2014                // second step: v = v-z,
2015                arg_index[res] = opind++;       
2016#endif
2017#if defined(_NONLIND_OLD_)
2018                extend_nonlinearity_domain_binary(arg1, arg2, ind_dom, nonl_dom);
2019#endif
2020#else
2021#if !defined(_ZOS_) /* BREAK_ZOS */
2022                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
2023                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
2024                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
2025
2026#ifdef _INT_FOR_
2027                FOR_0_LE_l_LT_p
2028                TRES_FOINC |= TARG2_INC | TARG1_INC;
2029#else
2030                /* olvo 980915 now in reverse order to allow x = x*x etc. */
2031                INC_pk_1(Tres)
2032                INC_pk_1(Targ1)
2033                INC_pk_1(Targ2)
2034
2035                FOR_p_GT_l_GE_0
2036                FOR_k_GT_i_GE_0
2037                { TRES_FODEC -= dp_T0[arg1]*TARG2_DEC +
2038                                TARG1_DEC*dp_T0[arg2];
2039                  DEC_TRES_FO
2040#if defined(_HIGHER_ORDER_)
2041                  Targ1OP = Targ1-i+1;
2042                  Targ2OP = Targ2;
2043
2044                  for (j=0;j<i;j++)
2045                  *Tres -= (*Targ1OP++) * (*Targ2OP--);
2046                  Tres--;
2047#endif /* _HIGHER_ORDER_ */
2048                }
2049#endif
2050#endif
2051#endif /* ALL_TOGETHER_AGAIN */
2052
2053#if !defined(_NTIGHT_)
2054                dp_T0[res] -= dp_T0[arg1] *
2055                                                    dp_T0[arg2];
2056#endif /* !_NTIGHT_ */
2057                break;
2058
2059                /*--------------------------------------------------------------------------*/
2060            case mult_d_a:         /* Multiply an adouble by a double    mult_d_a */
2061                /* (*) */
2062                arg   = get_locint_f();
2063                res   = get_locint_f();
2064                coval = get_val_f();
2065
2066                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2067
2068#if !defined(_NTIGHT_)
2069                dp_T0[res] = dp_T0[arg] * coval;
2070#endif /* !_NTIGHT_ */
2071
2072#if defined(_INDO_)
2073#if defined(_INDOPRO_)
2074                copy_index_domain(res, arg, ind_dom);
2075#endif
2076#if defined(_NONLIND_)
2077                arg_index[res] = arg_index[arg];               
2078#endif
2079#else
2080#if !defined(_ZOS_) /* BREAK_ZOS */
2081                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
2082                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
2083
2084#ifdef _INT_FOR_
2085                FOR_0_LE_l_LT_pk
2086                TRES_INC = TARG_INC;
2087#else
2088                FOR_0_LE_l_LT_pk
2089                TRES_INC = TARG_INC * coval;
2090#endif
2091#endif
2092#endif /* ALL_TOGETHER_AGAIN */
2093                break;
2094
2095                /*--------------------------------------------------------------------------*/
2096            case div_a_a:           /* Divide an adouble by an adouble    div_a_a */
2097                /* (/) */
2098                arg1 = get_locint_f();
2099                arg2 = get_locint_f();
2100                res  = get_locint_f();
2101
2102                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2103
2104#if !defined(_NTIGHT_)
2105#if !defined(_ZOS_) /* BREAK_ZOS */
2106                divs = 1.0 / dp_T0[arg2];
2107#endif /* ALL_TOGETHER_AGAIN */
2108
2109                dp_T0[res] = dp_T0[arg1] /
2110                                               dp_T0[arg2];
2111#endif /* !_NTIGHT_ */
2112
2113#if defined(_INDO_)
2114#if defined(_INDOPRO_)
2115                combine_2_index_domains(res, arg1, arg2, ind_dom);
2116#endif
2117#if defined(_NONLIND_)
2118                fod[opind].entry = maxopind+2;
2119                fod[opind].left = &fod[arg_index[arg1]];
2120                fod[opind].right = &fod[arg_index[arg2]];
2121                traverse_unary(&fod[arg_index[arg1]], nonl_dom, &fod[arg_index[arg2]], indcheck+1,maxopind+2);
2122                traverse_unary(&fod[arg_index[arg2]], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2123                arg_index[res] = opind++;               
2124#endif
2125#if defined(_NONLIND_OLD_)
2126                extend_nonlinearity_domain_binary(arg1, arg2, ind_dom, nonl_dom);
2127                extend_nonlinearity_domain_unary(arg2, ind_dom, nonl_dom);
2128#endif
2129#else
2130#if !defined(_ZOS_) /* BREAK_ZOS */
2131                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
2132                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
2133                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
2134
2135#ifdef _INT_FOR_
2136                FOR_0_LE_l_LT_p
2137                TRES_FOINC = TARG1_INC | TARG2_FOINC;
2138#else
2139                FOR_0_LE_l_LT_p
2140                FOR_0_LE_i_LT_k
2141                { /* olvo 980922 changed order to allow x = y/x */
2142#if defined(_HIGHER_ORDER_)
2143                    zOP      = dp_z+i;
2144                    (*zOP--) = -(*Targ2) * divs;
2145#endif /* _HIGHER_ORDER_ */
2146
2147                    TRES_FOINC = TARG1_INC * divs + dp_T0[res] *
2148                                 (-TARG2_INC * divs);
2149
2150#if defined(_HIGHER_ORDER_)
2151                    TresOP = Tres-i;
2152
2153                    for (j=0;j<i;j++)
2154                    *Tres += (*TresOP++) * (*zOP--);
2155                    Tres++;
2156#endif /* _HIGHER_ORDER_ */
2157                }
2158#endif
2159#endif
2160#endif /* ALL_TOGETHER_AGAIN */
2161                break;
2162
2163            /*--------------------------------------------------------------------------*/
2164        case div_d_a:             /* Division double - adouble (/)    div_d_a */
2165            arg   = get_locint_f();
2166                res   = get_locint_f();
2167                coval = get_val_f();
2168
2169                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2170
2171                /* olvo 980922 necessary for reverse */
2172                if (arg == res) {
2173                    IF_KEEP_WRITE_TAYLOR(arg,keep,k,p)
2174                }
2175
2176#if !defined(_NTIGHT_)
2177#if !defined(_ZOS_) /* BREAK_ZOS */
2178                divs = 1.0 / dp_T0[arg];
2179#endif /* ALL_TOGETHER_AGAIN */
2180
2181                dp_T0[res] = coval / dp_T0[arg];
2182#endif /* !_NTIGHT_ */
2183
2184#if defined(_INDO_)
2185#if defined(_INDOPRO_)
2186                copy_index_domain(res, arg, ind_dom);
2187#endif
2188#if defined(_NONLIND_)
2189                fod[opind].entry = maxopind+2;
2190                fod[opind].left = &fod[arg_index[arg]];
2191                fod[opind].right = NULL;
2192                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2193                arg_index[res] = opind++;               
2194#endif
2195#if defined(_NONLIND_OLD_)
2196                extend_nonlinearity_domain_unary(arg, ind_dom, nonl_dom);
2197#endif
2198#else
2199#if !defined(_ZOS_) /* BREAK_ZOS */
2200                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
2201                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
2202
2203#ifdef _INT_FOR_
2204                FOR_0_LE_l_LT_p
2205                TRES_FOINC = TARG_FOINC;
2206#else
2207                FOR_0_LE_l_LT_p
2208                FOR_0_LE_i_LT_k
2209                { /* olvo 980922 changed order to allow x = d/x */
2210#if defined(_HIGHER_ORDER_)
2211                    zOP      = dp_z+i;
2212                    (*zOP--) = -(*Targ) * divs;
2213#endif /* _HIGHER_ORDER_ */
2214
2215                    TRES_FOINC = dp_T0[res] * (-TARG_INC * divs);
2216
2217#if defined(_HIGHER_ORDER_)
2218                    TresOP = Tres-i;
2219
2220                    for (j=0;j<i;j++)
2221                    *Tres += (*TresOP++) * (*zOP--);
2222                    Tres++;
2223#endif /* _HIGHER_ORDER_ */
2224                }
2225#endif
2226#endif
2227#endif /* ALL_TOGETHER_AGAIN */
2228                break;
2229
2230
2231            /****************************************************************************/
2232            /*                                                         SIGN  OPERATIONS */
2233
2234            /*--------------------------------------------------------------------------*/
2235        case pos_sign_a:                                        /* pos_sign_a */
2236            arg   = get_locint_f();
2237                res   = get_locint_f();
2238
2239                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2240
2241#if !defined(_NTIGHT_)
2242                dp_T0[res] = dp_T0[arg];
2243#endif /* !_NTIGHT_ */
2244
2245#if defined(_INDO_)
2246#if defined(_INDOPRO_)
2247                copy_index_domain(res, arg, ind_dom);
2248#endif
2249#if defined(_NONLIND_)
2250                arg_index[res] = arg_index[arg];                           
2251#endif 
2252#else
2253#if !defined(_ZOS_) /* BREAK_ZOS */
2254                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
2255                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
2256
2257                FOR_0_LE_l_LT_pk
2258                TRES_INC = TARG_INC;
2259#endif
2260#endif /* ALL_TOGETHER_AGAIN */
2261                break;
2262
2263                /*--------------------------------------------------------------------------*/
2264            case neg_sign_a:                                        /* neg_sign_a */
2265                arg   = get_locint_f();
2266                res   = get_locint_f();
2267
2268                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2269
2270#if !defined(_NTIGHT_)
2271                dp_T0[res] = -dp_T0[arg];
2272#endif /* !_NTIGHT_ */
2273
2274#if defined(_INDO_)
2275#if defined(_INDOPRO_)
2276                copy_index_domain(res, arg, ind_dom);
2277#endif
2278#if defined(_NONLIND_)
2279                arg_index[res] = arg_index[arg];                           
2280#endif
2281#else
2282#if !defined(_ZOS_) /* BREAK_ZOS */
2283                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
2284                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
2285
2286#ifdef _INT_FOR_
2287                FOR_0_LE_l_LT_pk
2288                TRES_INC = TARG_INC;
2289#else
2290                FOR_0_LE_l_LT_pk
2291                TRES_INC = -TARG_INC;
2292#endif
2293#endif
2294#endif /* ALL_TOGETHER_AGAIN */
2295                break;
2296
2297
2298                /****************************************************************************/
2299                /*                                                         UNARY OPERATIONS */
2300
2301                /*--------------------------------------------------------------------------*/
2302            case exp_op:                          /* exponent operation    exp_op */
2303                arg = get_locint_f();
2304                res = get_locint_f();
2305
2306                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2307
2308#if !defined(_NTIGHT_)
2309                dp_T0[res] = exp(dp_T0[arg]);
2310#endif /* !_NTIGHT_ */
2311
2312                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
2313
2314#if defined(_INDO_)
2315#if defined(_INDOPRO_)
2316                copy_index_domain(res, arg, ind_dom);
2317#endif
2318#if defined(_NONLIND_)
2319                fod[opind].entry = maxopind+2;
2320                fod[opind].left = &fod[arg_index[arg]];
2321                fod[opind].right = NULL;
2322                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2323                arg_index[res] = opind++;               
2324#endif
2325#if defined(_NONLIND_OLD_)
2326                extend_nonlinearity_domain_unary(arg, ind_dom, nonl_dom);
2327#endif
2328#else
2329#if !defined(_ZOS_) /* BREAK_ZOS */
2330                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
2331                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
2332
2333#ifdef _INT_FOR_
2334                FOR_0_LE_l_LT_p
2335                TRES_FOINC = TARG_FOINC;
2336#else
2337                FOR_0_LE_l_LT_p
2338                FOR_0_LE_i_LT_k
2339                { /* olvo 980915 changed order to allow x = exp(x) */
2340#if defined(_HIGHER_ORDER_)
2341                    zOP      = dp_z+i;
2342                    (*zOP--) = (i+1) * (*Targ);
2343#endif /* _HIGHER_ORDER_ */
2344
2345                    TRES_FOINC = dp_T0[res] * TARG_INC;
2346
2347#if defined(_HIGHER_ORDER_)
2348                    TresOP = Tres-i;
2349
2350                    *Tres *= (i+1);
2351                    for (j=0;j<i;j++)
2352                    *Tres += (*TresOP++) * (*zOP--);
2353                    *Tres++ /= (i+1); /* important only for i>0 */
2354#endif /* _HIGHER_ORDER_ */
2355                }
2356
2357#endif
2358#endif
2359#endif /* ALL_TOGETHER_AGAIN */
2360                break;
2361
2362            /*--------------------------------------------------------------------------*/
2363        case sin_op:                              /* sine operation    sin_op */
2364                arg1 = get_locint_f();
2365                arg2 = get_locint_f();
2366                res  = get_locint_f();
2367
2368                IF_KEEP_WRITE_TAYLOR(arg2,keep,k,p) /* olvo 980710 covalue */
2369                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2370
2371#if !defined(_NTIGHT_)
2372                /* Note: always arg2 != arg1 */
2373                dp_T0[arg2] = cos(dp_T0[arg1]);
2374                dp_T0[res]  = sin(dp_T0[arg1]);
2375#endif /* !_NTIGHT_ */
2376
2377                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
2378
2379#if defined(_INDO_)
2380#if defined(_INDOPRO_)
2381                copy_index_domain(res, arg1, ind_dom);
2382#endif
2383#if defined(_NONLIND_)
2384                fod[opind].entry = maxopind+2;
2385                fod[opind].left = &fod[arg_index[arg1]];
2386                fod[opind].right = NULL;
2387                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2388                arg_index[res] = opind++;               
2389#endif
2390#if defined(_NONLIND_OLD_)
2391                extend_nonlinearity_domain_unary(arg1, ind_dom, nonl_dom);
2392#endif
2393#else
2394#if !defined(_ZOS_) /* BREAK_ZOS */
2395                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
2396                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
2397                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
2398
2399#ifdef _INT_FOR_
2400                FOR_0_LE_l_LT_p
2401                { /* olvo 980923 changed order to allow x = sin(x) */
2402                    TARG2_FOINC =  TARG1;
2403                    TRES_FOINC  =  TARG1_FOINC;
2404            }
2405#else
2406                FOR_0_LE_l_LT_p
2407                FOR_0_LE_i_LT_k
2408                { /* olvo 980921 changed order to allow x = sin(x) */
2409#if defined(_HIGHER_ORDER_)
2410                    zOP      = dp_z+i;
2411                    (*zOP--) = (i+1) * (*Targ1);
2412#endif /* _HIGHER_ORDER_ */
2413
2414                    /* Note: always arg2 != arg1 */
2415                    TARG2_FOINC = -dp_T0[res]  * TARG1;
2416                    TRES_FOINC  =  dp_T0[arg2] * TARG1_INC;
2417
2418#if defined(_HIGHER_ORDER_)
2419                    TresOP  = Tres-i;
2420                    Targ2OP = Targ2-i;
2421
2422                    *Tres  *= (i+1);
2423                    *Targ2 *= (i+1);
2424                    for (j=0;j<i;j++) {
2425                    *Tres  += (*Targ2OP++) * (*zOP);
2426                        *Targ2 -= (*TresOP++)  * (*zOP--);
2427                    }
2428                    *Targ2++ /= (i+1);
2429                    *Tres++  /= (i+1);
2430#endif /* _HIGHER_ORDER_ */
2431            }
2432#endif
2433#endif
2434#endif /* ALL_TOGETHER_AGAIN */
2435                break;
2436
2437                /*--------------------------------------------------------------------------*/
2438            case cos_op:                            /* cosine operation    cos_op */
2439                arg1 = get_locint_f();
2440                arg2 = get_locint_f();
2441                res  = get_locint_f();
2442
2443                IF_KEEP_WRITE_TAYLOR(arg2,keep,k,p) /* olvo 980710 covalue */
2444                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2445
2446#if !defined(_NTIGHT_)
2447                /* Note: always arg2 != arg1 */
2448                dp_T0[arg2] = sin(dp_T0[arg1]);
2449                dp_T0[res]  = cos(dp_T0[arg1]);
2450#endif /* !_NTIGHT_ */
2451
2452                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
2453
2454#if defined(_INDO_)
2455#if defined(_INDOPRO_)
2456                copy_index_domain(res, arg1, ind_dom);
2457#endif
2458#if defined(_NONLIND_)
2459                fod[opind].entry = maxopind+2;
2460                fod[opind].left = &fod[arg_index[arg1]];
2461                fod[opind].right = NULL;
2462                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2463                arg_index[res] = opind++;               
2464#endif
2465#if defined(_NONLIND_OLD_)
2466                extend_nonlinearity_domain_unary(arg1, ind_dom, nonl_dom);
2467#endif
2468
2469#else
2470#if !defined(_ZOS_) /* BREAK_ZOS */
2471                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
2472                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
2473                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
2474
2475#ifdef _INT_FOR_
2476                FOR_0_LE_l_LT_p
2477                { /* olvo 980923 changed order to allow x = cos(x) */
2478                    TARG2_FOINC = TARG1;
2479                    TRES_FOINC  = TARG1_FOINC;
2480            }
2481#else
2482                FOR_0_LE_l_LT_p
2483                FOR_0_LE_i_LT_k
2484                { /* olvo 980921 changed order to allow x = cos(x) */
2485#if defined(_HIGHER_ORDER_)
2486                    zOP      = dp_z+i;
2487                    (*zOP--) = (i+1) * (*Targ1);
2488#endif /* _HIGHER_ORDER_ */
2489
2490                    /* Note: always arg2 != arg1 */
2491                    TARG2_FOINC =  dp_T0[res]  * TARG1;
2492                    TRES_FOINC  = -dp_T0[arg2] * TARG1_INC;
2493
2494#if defined(_HIGHER_ORDER_)
2495                    TresOP  = Tres-i;
2496                    Targ2OP = Targ2-i;
2497
2498                    *Tres  *= (i+1);
2499                    *Targ2 *= (i+1);
2500                    for (j=0;j<i;j++) {
2501                    *Tres  -= (*Targ2OP++) * (*zOP);
2502                        *Targ2 += (*TresOP++)  * (*zOP--);
2503                    }
2504                    *Targ2++ /= (i+1);
2505                    *Tres++  /= (i+1);
2506#endif /* _HIGHER_ORDER_ */
2507            }
2508#endif
2509#endif
2510#endif /* ALL_TOGETHER_AGAIN */
2511                break;
2512
2513                /*--------------------------------------------------------------------------*/
2514            case atan_op:                                              /* atan_op */
2515                arg1 = get_locint_f();
2516                arg2 = get_locint_f();
2517                res  = get_locint_f();
2518
2519                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2520
2521#if !defined(_NTIGHT_)
2522                dp_T0[res]=atan(dp_T0[arg1]);
2523#endif /* !_NTIGHT_ */
2524
2525                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
2526
2527#if defined(_INDO_)
2528#if defined(_INDOPRO_)
2529                copy_index_domain(res, arg1, ind_dom);
2530#endif
2531#if defined(_NONLIND_)
2532                fod[opind].entry = maxopind+2;
2533                fod[opind].left = &fod[arg_index[arg1]];
2534                fod[opind].right = NULL;
2535                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2536                arg_index[res] = opind++;               
2537#endif
2538#if defined(_NONLIND_OLD_)
2539                extend_nonlinearity_domain_unary(arg1, ind_dom, nonl_dom);
2540#endif
2541#else
2542#if !defined(_ZOS_) /* BREAK_ZOS */
2543                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
2544                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
2545#ifdef _INT_FOR_
2546                FOR_0_LE_l_LT_p
2547                TRES_FOINC = TARG1_FOINC;
2548#else
2549                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
2550
2551                FOR_0_LE_l_LT_p
2552                { FOR_0_LE_i_LT_k
2553                  { /* olvo 980921 changed order to allow x = atan(x) */
2554#if defined(_HIGHER_ORDER_)
2555                      zOP      = dp_z+i;
2556                      (*zOP--) = (i+1) * (*Targ1);
2557#endif /* _HIGHER_ORDER_ */
2558
2559                      TRES_FOINC = dp_T0[arg2] * TARG1_INC;
2560
2561#if defined(_HIGHER_ORDER_)
2562                      Targ2OP = Targ2;
2563
2564                      *Tres *= (i+1);
2565                      for (j=0;j<i;j++)
2566                      *Tres  += (*Targ2OP++) * (*zOP--);
2567                      *Tres++ /= (i+1);
2568#endif /* _HIGHER_ORDER_ */
2569                  }
2570                  HOV_INC(Targ2, k)
2571                }
2572#endif
2573#endif
2574#endif /* ALL_TOGETHER_AGAIN */
2575                break;
2576
2577            /*--------------------------------------------------------------------------*/
2578        case asin_op:                                              /* asin_op */
2579            arg1 = get_locint_f();
2580                arg2 = get_locint_f();
2581                res  = get_locint_f();
2582
2583                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2584
2585#if !defined(_NTIGHT_)
2586                dp_T0[res] = asin(dp_T0[arg1]);
2587#endif /* !_NTIGHT_ */
2588
2589                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
2590
2591#if defined(_INDO_)
2592#if defined(_INDOPRO_)
2593                copy_index_domain(res, arg1, ind_dom);
2594#endif
2595#if defined(_NONLIND_)
2596                fod[opind].entry = maxopind+2;
2597                fod[opind].left = &fod[arg_index[arg1]];
2598                fod[opind].right = NULL;
2599                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2600                arg_index[res] = opind++;               
2601#endif
2602#if defined(_NONLIND_OLD_)
2603                extend_nonlinearity_domain_unary(arg1, ind_dom, nonl_dom);
2604#endif
2605#else
2606#if !defined(_ZOS_) /* BREAK_ZOS */
2607                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
2608                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
2609#ifdef _INT_FOR_
2610                FOR_0_LE_l_LT_p
2611                TRES_FOINC = TARG1_FOINC;
2612#else
2613                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
2614
2615                if (dp_T0[arg1] == 1.0)
2616                    FOR_0_LE_l_LT_p
2617                    { FOR_0_LE_i_LT_k
2618                      if (TARG1 > 0.0) {
2619                      r0 = make_nan();
2620                          VEC_INC(Targ1, k-i)
2621                          BREAK_FOR_I
2622                      } else
2623                          if (TARG1 < 0.0) {
2624                          r0 = make_inf();
2625                              VEC_INC(Targ1, k-i)
2626                              BREAK_FOR_I
2627                          } else {
2628                              r0 = 0.0;
2629                              Targ1++;
2630                          }
2631                  TRES = r0;
2632                  VEC_INC(Tres, k)
2633            } else
2634                    if (dp_T0[arg1] == -1.0)
2635                        FOR_0_LE_l_LT_p
2636                        { FOR_0_LE_i_LT_k
2637                          if (TARG1 > 0.0) {
2638                          r0 = make_inf();
2639                              VEC_INC(Targ1, k-i)
2640                              BREAK_FOR_I
2641                          } else
2642                              if (TARG1 < 0.0) {
2643                              r0 = make_nan();
2644                                  VEC_INC(Targ1, k-i)
2645                                  BREAK_FOR_I
2646                              } else {
2647                                  r0 = 0.0;
2648                                  Targ1++;
2649                              }
2650                  TRES = r0;
2651                  VEC_INC(Tres, k)
2652                } else
2653                        FOR_0_LE_l_LT_p {
2654                            FOR_0_LE_i_LT_k
2655                            { /* olvo 980921 changed order to allow x = asin(x) */
2656#if defined(_HIGHER_ORDER_)
2657                                zOP      = dp_z+i;
2658                                (*zOP--) = (i+1) * (*Targ1);
2659#endif /* _HIGHER_ORDER_ */
2660
2661                                TRES_FOINC = dp_T0[arg2] * TARG1_INC;
2662
2663#if defined(_HIGHER_ORDER_)
2664                                Targ2OP = Targ2;
2665
2666                                *Tres *= (i+1);
2667                                for (j=0;j<i;j++)
2668                                *Tres += (*Targ2OP++) * (*zOP--);
2669                                *Tres++ /= (i+1);
2670#endif /* _HIGHER_ORDER_ */
2671                            }
2672                            HOV_INC(Targ2, k)
2673                        }
2674#endif
2675#endif
2676#endif /* ALL_TOGETHER_AGAIN */
2677                        break;
2678
2679            /*--------------------------------------------------------------------------*/
2680        case acos_op:                                              /* acos_op */
2681                arg1 = get_locint_f();
2682                arg2 = get_locint_f();
2683                res  = get_locint_f();
2684
2685                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2686
2687#if !defined(_NTIGHT_)
2688                dp_T0[res] = acos(dp_T0[arg1]);
2689#endif /* !_NTIGHT_ */
2690
2691                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
2692
2693#if defined(_INDO_)
2694#if defined(_INDOPRO_)
2695                copy_index_domain(res, arg1, ind_dom);
2696#endif
2697#if defined(_NONLIND_)
2698                fod[opind].entry = maxopind+2;
2699                fod[opind].left = &fod[arg_index[arg1]];
2700                fod[opind].right = NULL;
2701                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2702                arg_index[res] = opind++;               
2703#endif
2704#if defined(_NONLIND_OLD_)
2705                extend_nonlinearity_domain_unary(arg1, ind_dom, nonl_dom);
2706#endif
2707#else
2708#if !defined(_ZOS_) /* BREAK_ZOS */
2709                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
2710                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
2711#ifdef _INT_FOR_
2712                FOR_0_LE_l_LT_p
2713                TRES_FOINC = TARG1_FOINC;
2714#else
2715                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
2716
2717                if (dp_T0[arg1] == 1.0)
2718                    FOR_0_LE_l_LT_p
2719                    { FOR_0_LE_i_LT_k
2720                      if (TARG1 > 0.0) {
2721                      r0 = make_nan();
2722                          VEC_INC(Targ1, k-i)
2723                          BREAK_FOR_I
2724                      } else
2725                          if (TARG1 < 0.0) {
2726                          r0 = -make_inf();
2727                              VEC_INC(Targ1, k-i)
2728                              BREAK_FOR_I
2729                          } else {
2730                              r0 = 0.0;
2731                              Targ1++;
2732                          }
2733                  TRES = r0;
2734                  VEC_INC(Tres, k)
2735            } else
2736                    if (dp_T0[arg1] == -1.0)
2737                        FOR_0_LE_l_LT_p
2738                        { FOR_0_LE_i_LT_k
2739                          if (TARG1 > 0.0) {
2740                          r0 = -make_inf();
2741                              VEC_INC(Targ1, k-i)
2742                              BREAK_FOR_I
2743                          } else
2744                              if (TARG1 < 0.0) {
2745                              r0 = make_nan();
2746                                  VEC_INC(Targ1, k-i)
2747                                  BREAK_FOR_I
2748                              } else {
2749                                  r0 = 0.0;
2750                                  Targ1++;
2751                              }
2752                  TRES = r0;
2753                  VEC_INC(Tres, k)
2754                } else
2755                        FOR_0_LE_l_LT_p {
2756                            FOR_0_LE_i_LT_k
2757                            { /* olvo 980921 changed order to allow x = acos(x) */
2758#if defined(_HIGHER_ORDER_)
2759                                zOP      = dp_z+i;
2760                                (*zOP--) = (i+1) * (*Targ1);
2761#endif /* _HIGHER_ORDER_ */
2762
2763                                TRES_FOINC = dp_T0[arg2] * TARG1_INC;
2764
2765#if defined(_HIGHER_ORDER_)
2766                                Targ2OP = Targ2;
2767
2768                                *Tres *= (i+1);
2769                                for (j=0;j<i;j++)
2770                                *Tres += (*Targ2OP++) * (*zOP--);
2771                                *Tres++ /= (i+1);
2772#endif /* _HIGHER_ORDER_ */
2773                            }
2774                            HOV_INC(Targ2, k)
2775                        }
2776#endif
2777#endif
2778#endif /* ALL_TOGETHER_AGAIN */
2779                        break;
2780
2781#ifdef ATRIG_ERF
2782
2783            /*--------------------------------------------------------------------------*/
2784        case asinh_op:                                            /* asinh_op */
2785                arg1 = get_locint_f();
2786                arg2 = get_locint_f();
2787                res  = get_locint_f();
2788
2789                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2790
2791#if !defined(_NTIGHT_)
2792                dp_T0[res] = asinh(dp_T0[arg1]);
2793#endif /* !_NTIGHT_ */
2794
2795                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
2796
2797#if defined(_INDO_)
2798#if defined(_INDOPRO_)
2799                copy_index_domain(res, arg1, ind_dom);
2800#endif
2801#if defined(_NONLIND_)
2802                fod[opind].entry = maxopind+2;
2803                fod[opind].left = &fod[arg_index[arg1]];
2804                fod[opind].right = NULL;
2805                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2806                arg_index[res] = opind++;               
2807#endif
2808#if defined(_NONLIND_OLD_)
2809                extend_nonlinearity_domain_unary(arg1, ind_dom, nonl_dom);
2810#endif
2811#else
2812#if !defined(_ZOS_) /* BREAK_ZOS */
2813                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
2814                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
2815#ifdef _INT_FOR_
2816                FOR_0_LE_l_LT_p
2817                TRES_FOINC = TARG1_FOINC;
2818#else
2819                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
2820
2821                FOR_0_LE_l_LT_p
2822                { FOR_0_LE_i_LT_k
2823                  { /* olvo 980921 changed order to allow x = asinh(x) */
2824#if defined(_HIGHER_ORDER_)
2825                      zOP      = dp_z+i;
2826                      (*zOP--) = (i+1) * (*Targ1);
2827#endif /* _HIGHER_ORDER_ */
2828
2829                      TRES_FOINC = dp_T0[arg2] * TARG1_INC;
2830
2831#if defined(_HIGHER_ORDER_)
2832                      Targ2OP = Targ2;
2833
2834                      *Tres *= (i+1);
2835                      for (j=0;j<i;j++)
2836                      *Tres += (*Targ2OP++) * (*zOP--);
2837                      *Tres++ /= (i+1);
2838#endif /* _HIGHER_ORDER_ */
2839                  }
2840                  HOV_INC(Targ2, k)
2841                }
2842#endif
2843#endif
2844#endif /* ALL_TOGETHER_AGAIN */
2845                break;
2846
2847            /*--------------------------------------------------------------------------*/
2848        case acosh_op:                                           /* acosh_op */
2849                arg1 = get_locint_f();
2850                arg2 = get_locint_f();
2851                res  = get_locint_f();
2852
2853                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2854
2855#if !defined(_NTIGHT_)
2856                dp_T0[res] = acosh(dp_T0[arg1]);
2857#endif /* !_NTIGHT_ */
2858
2859                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
2860
2861#if defined(_INDO_)
2862#if defined(_INDOPRO_)
2863                copy_index_domain(res, arg1, ind_dom);
2864#endif
2865#if defined(_NONLIND_)
2866                fod[opind].entry = maxopind+2;
2867                fod[opind].left = &fod[arg_index[arg1]];
2868                fod[opind].right = NULL;
2869                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2870                arg_index[res] = opind++;               
2871#endif
2872#if defined(_NONLIND_OLD_)
2873                extend_nonlinearity_domain_unary(arg1, ind_dom, nonl_dom);
2874#endif
2875#else
2876#if !defined(_ZOS_) /* BREAK_ZOS */
2877                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
2878                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
2879#ifdef _INT_FOR_
2880                FOR_0_LE_l_LT_p
2881                TRES_FOINC = TARG1_FOINC;
2882#else
2883                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
2884
2885                if (dp_T0[arg1] == 1.0)
2886                    FOR_0_LE_l_LT_p
2887                    { FOR_0_LE_i_LT_k
2888                      if (TARG1 > 0.0) {
2889                      r0 = make_inf();
2890                          VEC_INC(Targ1, k-i)
2891                          BREAK_FOR_I
2892                      } else
2893                          if (TARG1 < 0.0) {
2894                          r0 = make_nan();
2895                              VEC_INC(Targ1, k-i)
2896                              BREAK_FOR_I
2897                          } else {
2898                              r0 = 0.0;
2899                              Targ1++;
2900                          }
2901                  TRES_INC = r0;
2902#if defined(_HIGHER_ORDER_)
2903                  for (i=1;i<k;i++)
2904                  *Tres++ = make_nan();
2905#endif /* _HIGHER_ORDER_ */
2906                } else
2907                    FOR_0_LE_l_LT_p {
2908                        FOR_0_LE_i_LT_k
2909                        { /* olvo 980921 changed order to allow x = acosh(x) */
2910#if defined(_HIGHER_ORDER_)
2911                            zOP      = dp_z+i;
2912                            (*zOP--) = (i+1) * (*Targ1);
2913#endif /* _HIGHER_ORDER_ */
2914
2915                            TRES_FOINC = dp_T0[arg2] * TARG1_INC;
2916
2917#if defined(_HIGHER_ORDER_)
2918                            Targ2OP = Targ2;
2919
2920                            *Tres *= (i+1);
2921                            for (j=0;j<i;j++)
2922                                *Tres += (*Targ2OP++) * (*zOP--);
2923                                *Tres++ /= (i+1);
2924#endif /* _HIGHER_ORDER_ */
2925                            }
2926                            HOV_INC(Targ2, k)
2927                        }
2928#endif
2929#endif
2930#endif /* ALL_TOGETHER_AGAIN */
2931                        break;
2932
2933            /*--------------------------------------------------------------------------*/
2934        case atanh_op:                                            /* atanh_op */
2935                arg1 = get_locint_f();
2936                arg2 = get_locint_f();
2937                res  = get_locint_f();
2938
2939                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
2940
2941#if !defined(_NTIGHT_)
2942                dp_T0[res] = atanh(dp_T0[arg1]);
2943#endif /* !_NTIGHT_ */
2944
2945                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
2946
2947#if defined(_INDO_)
2948#if defined(_INDOPRO_)
2949                copy_index_domain(res, arg1, ind_dom);
2950#endif
2951#if defined(_NONLIND_)
2952                fod[opind].entry = maxopind+2;
2953                fod[opind].left = &fod[arg_index[arg1]];
2954                fod[opind].right = NULL;
2955                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
2956                arg_index[res] = opind++;               
2957#endif
2958#if defined(_NONLIND_OLD_)
2959                extend_nonlinearity_domain_unary(arg1, ind_dom, nonl_dom);
2960#endif
2961#else
2962#if !defined(_ZOS_) /* BREAK_ZOS */
2963                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
2964                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
2965#ifdef _INT_FOR_
2966                FOR_0_LE_l_LT_p
2967                TRES_FOINC = TARG1_FOINC;
2968#else
2969                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
2970
2971                if (dp_T0[arg1] == 1.0)
2972                    FOR_0_LE_l_LT_p
2973                    { FOR_0_LE_i_LT_k
2974                      if (TARG1 > 0.0) {
2975                      r0 = make_nan();
2976                          VEC_INC(Targ1, k-i)
2977                          BREAK_FOR_I
2978                      } else
2979                          if (TARG1 < 0.0) {
2980                          r0 = make_inf();
2981                              VEC_INC(Targ1, k-i)
2982                              BREAK_FOR_I
2983                          } else {
2984                              r0 = 0.0;
2985                              Targ1++;
2986                          }
2987                  TRES_INC = r0;
2988#if defined(_HIGHER_ORDER_)
2989                  for (i=1;i<k;i++)
2990                  *Tres++ = make_nan();
2991#endif /* _HIGHER_ORDER_ */
2992                } else
2993                    if (dp_T0[arg1] == -1.0)
2994                            FOR_0_LE_l_LT_p
2995                            { FOR_0_LE_i_LT_k
2996                              if (TARG1 > 0.0) {
2997                              r0 = make_inf();
2998                                  VEC_INC(Targ1, k-i)
2999                                  BREAK_FOR_I
3000                              } else
3001                                  if (TARG1 < 0.0) {
3002                                  r0 = make_nan();
3003                                      VEC_INC(Targ1, k-i)
3004                                      BREAK_FOR_I
3005                                  } else {
3006                                      r0 = 0.0;
3007                                      Targ1++;
3008                                  }
3009                  TRES_INC = r0;
3010#if defined(_HIGHER_ORDER_)
3011                  for (i=1;i<k;i++)
3012                  *Tres++ = make_nan();
3013#endif /* _HIGHER_ORDER_ */
3014                        } else
3015                            FOR_0_LE_l_LT_p {
3016                                FOR_0_LE_i_LT_k
3017                                { /* olvo 980921 changed order to allow x = atanh(x) */
3018#if defined(_HIGHER_ORDER_)
3019                                    zOP      = dp_z+i;
3020                                    (*zOP--) = (i+1) * (*Targ1);
3021#endif /* _HIGHER_ORDER_ */
3022
3023                                    TRES_FOINC = dp_T0[arg2] * TARG1_INC;
3024
3025#if defined(_HIGHER_ORDER_)
3026                                    Targ2OP = Targ2;
3027
3028                                    *Tres *= (i+1);
3029                                    for (j=0;j<i;j++)
3030                                        *Tres += (*Targ2OP++) * (*zOP--);
3031                                        *Tres++ /= (i+1);
3032#endif /* _HIGHER_ORDER_ */
3033                                    }
3034                                    HOV_INC(Targ2, k)
3035                                }
3036#endif
3037#endif
3038#endif /* ALL_TOGETHER_AGAIN */
3039                                break;
3040
3041            /*--------------------------------------------------------------------------*/
3042        case erf_op:                                                /* erf_op */
3043                arg1 = get_locint_f();
3044                arg2 = get_locint_f();
3045                res  = get_locint_f();
3046
3047                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
3048
3049#if !defined(_NTIGHT_)
3050                dp_T0[res] = erf(dp_T0[arg1]);
3051#endif /* !_NTIGHT_ */
3052
3053                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
3054#if defined(_INDO_)
3055#if defined(_INDOPRO_)
3056                copy_index_domain(res, arg1, ind_dom);
3057#endif
3058#if defined(_NONLIND_)
3059                fod[opind].entry = maxopind+2;
3060                fod[opind].left = &fod[arg_index[arg1]];
3061                fod[opind].right = NULL;
3062                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
3063                arg_index[res] = opind++;               
3064#endif       
3065#else
3066#if !defined(_ZOS_) /* BREAK_ZOS */
3067                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
3068                ASSIGN_T(Targ1,TAYLOR_BUFFER[arg1])
3069#ifdef _INT_FOR_
3070                FOR_0_LE_l_LT_p
3071                TRES_FOINC = TARG1_FOINC;
3072#else
3073                ASSIGN_T(Targ2,TAYLOR_BUFFER[arg2])
3074
3075                FOR_0_LE_l_LT_p
3076                { FOR_0_LE_i_LT_k
3077                  { /* olvo 980921 changed order to allow x = erf(x) */
3078#if defined(_HIGHER_ORDER_)
3079                      zOP      = dp_z+i;
3080                      (*zOP--) = (i+1) * (*Targ1);
3081#endif /* _HIGHER_ORDER_ */
3082
3083                      TRES_FOINC = dp_T0[arg2] * TARG1_INC;
3084
3085#if defined(_HIGHER_ORDER_)
3086                      Targ2OP = Targ2;
3087
3088                      *Tres *= (i+1);
3089                      for (j=0;j<i;j++)
3090                      *Tres += (*Targ2OP++) * (*zOP--);
3091                      *Tres++ /= (i+1);
3092#endif /* _HIGHER_ORDER_ */
3093                  }
3094                  HOV_INC(Targ2, k)
3095                }
3096#endif
3097#endif
3098#endif /* ALL_TOGETHER_AGAIN */
3099                break;
3100
3101#endif
3102
3103            /*--------------------------------------------------------------------------*/
3104        case log_op:                                                /* log_op */
3105                arg = get_locint_f();
3106                res = get_locint_f();
3107
3108                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
3109
3110#if defined(_INDO_)
3111#if defined(_INDOPRO_)
3112                copy_index_domain(res, arg, ind_dom);
3113#endif
3114#if defined(_NONLIND_)
3115                fod[opind].entry = maxopind+2;
3116                fod[opind].left = &fod[arg_index[arg]];
3117                fod[opind].right = NULL;
3118                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
3119                arg_index[res] = opind++;               
3120#endif
3121#if defined(_NONLIND_OLD_)
3122                extend_nonlinearity_domain_unary(arg, ind_dom, nonl_dom);
3123#endif
3124#else
3125#if !defined(_ZOS_) /* BREAK_ZOS */
3126                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
3127                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
3128
3129#ifdef _INT_FOR_
3130                FOR_0_LE_l_LT_p
3131                TRES_FOINC = TARG_INC;
3132#else
3133                divs = 1.0 / dp_T0[arg];
3134                FOR_0_LE_l_LT_p
3135                { if (dp_T0[arg] == 0.0) {
3136                  TargOP = Targ;
3137                  FOR_0_LE_i_LT_k
3138                  { if (*TargOP++ < 0.0) {
3139                        divs = make_nan();
3140                            BREAK_FOR_I
3141                        }
3142                      }
3143                  }
3144
3145                  /* olvo 980921 changed order to allow x = log(x) */
3146                  FOR_0_LE_i_LT_k
3147                  { TRES_FOINC = TARG_INC * divs;
3148#if defined(_HIGHER_ORDER_)
3149                    TresOP = Tres - i;
3150                    zOP    = dp_z+i;
3151
3152                    (*zOP--) = *Tres;
3153                    (*Tres) *= i+1;
3154                    for (j=0;j<i;j++)
3155                    (*Tres) -= (*zOP--) * (*TresOP++) * (j+1);
3156                    *Tres++ /= i+1;
3157#endif /* _HIGHER_ORDER_ */
3158                  }
3159                }
3160#endif
3161#endif
3162#endif /* ALL_TOGETHER_AGAIN */
3163#if !defined(_NTIGHT_)
3164                dp_T0[res] = log(dp_T0[arg]);
3165#endif /* !_NTIGHT_ */
3166
3167                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
3168                break;
3169
3170                /*--------------------------------------------------------------------------*/
3171            case pow_op:                                                /* pow_op */
3172                arg   = get_locint_f();
3173                res   = get_locint_f();
3174
3175                coval = get_val_f();
3176
3177                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
3178
3179                /* olvo 980921 necessary for reverse */
3180                if (arg == res) {
3181                    IF_KEEP_WRITE_TAYLOR(arg,keep,k,p)
3182                }
3183
3184#if !defined(_NTIGHT_)
3185
3186#ifndef _ZOS_ /* BREAK_ZOS */
3187#if !defined(_INT_FOR_)
3188                T0arg   = dp_T0[arg];
3189#endif
3190#endif /* ALL_TOGETHER_AGAIN */
3191
3192                dp_T0[res] =
3193                    pow(dp_T0[arg], coval);
3194#endif /* !_NTIGHT_ */
3195
3196                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
3197#if defined(_INDO_)
3198#if defined(_INDOPRO_)
3199                copy_index_domain(res, arg, ind_dom);
3200#endif
3201#if defined(_NONLIND_)
3202                fod[opind].entry = maxopind+2;
3203                fod[opind].left = &fod[arg_index[arg]];
3204                fod[opind].right = NULL;
3205                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
3206                arg_index[res] = opind++;               
3207#endif
3208#if defined(_NONLIND_OLD_)
3209                extend_nonlinearity_domain_unary(arg, ind_dom, nonl_dom);
3210#endif
3211#else
3212#ifndef _ZOS_ /* BREAK_ZOS */
3213                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
3214                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
3215
3216#ifdef _INT_FOR_
3217                FOR_0_LE_l_LT_p
3218                TRES_FOINC = TARG_INC;
3219#else
3220                if (T0arg == 0.0) {
3221                    if (coval <= 0.0)
3222                        FOR_0_LE_l_LT_pk
3223                        TRES_INC = make_nan();
3224                    else {
3225                        /* coval not a whole number */
3226                        if (coval - floor(coval) != 0) {
3227                            FOR_0_LE_l_LT_p
3228                            {
3229                                i = 0;
3230                                FOR_0_LE_i_LT_k
3231                                {
3232                                    if (coval - i > 1)
3233                                    TRES_INC = 0;
3234                                    if ((coval - i < 1) && (coval - i > 0))
3235                                        TRES_INC = make_inf();
3236                                        if (coval - i < 0)
3237                                            TRES_INC = make_nan();
3238                                        }
3239                                    }
3240                                } else {
3241                        if (coval == 1) {
3242                                FOR_0_LE_l_LT_pk
3243                                TRES_INC = TARG_INC;
3244                            } else
3245                                /* coval is an int > 1 */
3246                                /* the following is not efficient but at least it works */
3247                                /* it reformulates x^n into x* ... *x n times */
3248                            {
3249                                INC_pk_1(Targ)
3250                                INC_pk_1(Tres)
3251
3252                                FOR_p_GT_l_GE_0
3253                                {
3254                                    FOR_k_GT_i_GE_0
3255                                    {
3256                                        *Tres = 0;
3257                                        DEC_TRES_FO
3258#if defined(_HIGHER_ORDER_)
3259                                        if (i == k-1) {
3260                                        zOP = dp_z+k-1;
3261                                        for(j=k-1;j>=0;j--) {
3262                                                (*zOP--) = (*Targ--);
3263                                            }
3264                                        }
3265                                        for (j=0;j<i;j++) {
3266                                        *Tres += dp_z[j] *
3267                                                     dp_z[i-j-1];
3268                                        }
3269                                        Tres--;
3270#endif /* _HIGHER_ORDER_ */
3271                                    }
3272                                }
3273                                for(ii=3;ii<=coval;ii++) {
3274                                    ASSIGN_T(Tres, TAYLOR_BUFFER[res])
3275                                    ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
3276                                    INC_pk_1(Targ)
3277                                    INC_pk_1(Tres)
3278
3279                                    FOR_p_GT_l_GE_0
3280                                    {
3281                                        FOR_k_GT_i_GE_0
3282                                        {
3283                                            *Tres = 0;
3284                                            DEC_TRES_FO
3285#if defined(_HIGHER_ORDER_)
3286                                            TresOP = Tres-i;
3287                                            for (j=0;j<i;j++)
3288                                            *Tres += TresOP[j] * dp_z[i-j-1];
3289                                            Tres--;
3290#endif /* _HIGHER_ORDER_ */
3291                                        }
3292                                    }
3293                                }
3294                        }
3295                    }
3296                }
3297            } else {
3298                r0 = 1.0 / T0arg;
3299                FOR_0_LE_l_LT_p
3300                FOR_0_LE_i_LT_k
3301                { /* olvo 980921 changed order to allow x = pow(x,n) */
3302#ifdef _HIGHER_ORDER_
3303                    zOP      = dp_z+i;
3304                    (*zOP--) = (*Targ) * r0;
3305#endif /* _HIGHER_ORDER_ */
3306
3307                    TRES_FOINC = dp_T0[res] *
3308                                 TARG_INC * coval * r0;
3309
3310#ifdef _HIGHER_ORDER_
3311                    TresOP = Tres-i;
3312
3313                    (*Tres) *= i+1;
3314                    y = coval*i -1;
3315                    for (j=0;j<i;j++) {
3316                        *Tres += (*TresOP++) * (*zOP--) * y;
3317                            y -= coval + 1;
3318                        }
3319                        *Tres++ /= (i+1);
3320#endif /* _HIGHER_ORDER_ */
3321                    }
3322                }
3323#endif
3324#endif
3325#endif /* ALL_TOGETHER_AGAIN */
3326                break;
3327
3328                /*--------------------------------------------------------------------------*/
3329            case sqrt_op:                                              /* sqrt_op */
3330                arg = get_locint_f();
3331                res = get_locint_f();
3332
3333                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
3334
3335#if !defined(_NTIGHT_)
3336                dp_T0[res] = sqrt(dp_T0[arg]);
3337#endif /* !_NTIGHT_ */
3338
3339                ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
3340
3341#if defined(_INDO_)
3342#if defined(_INDOPRO_)
3343                copy_index_domain(res, arg, ind_dom);
3344#endif
3345#if defined(_NONLIND_)
3346                fod[opind].entry = maxopind+2;
3347                fod[opind].left = &fod[arg_index[arg]];
3348                fod[opind].right = NULL;
3349                traverse_unary(&fod[opind], nonl_dom, &fod[opind], indcheck+1,maxopind+2);
3350                arg_index[res] = opind++;               
3351#endif
3352#if defined(_NONLIND_OLD_)
3353                extend_nonlinearity_domain_unary(arg, ind_dom, nonl_dom);
3354#endif
3355#else
3356#if !defined(_ZOS_) /* BREAK_ZOS */
3357                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
3358                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
3359
3360#ifdef _INT_FOR_
3361                FOR_0_LE_l_LT_p
3362                TRES_FOINC = TARG_INC;
3363#else
3364                FOR_0_LE_l_LT_p
3365                { TargOP = Targ;
3366                  if (dp_T0[arg] == 0.0)
3367                  /* Note: <=> dp_T0[res] == 0.0 */
3368              { r0 = 0.0;
3369                  FOR_0_LE_i_LT_k
3370                  { if (TARG>0.0) {
3371                        r0 = make_inf();
3372                            VEC_INC(Targ, k-i)
3373                            BREAK_FOR_I
3374                        } else
3375                            if (TARG<0.0) {
3376                            r0 = make_nan();
3377                                VEC_INC(Targ, k-i)
3378                                BREAK_FOR_I
3379                            } else
3380                                Targ++;
3381                              }
3382                          }
3383                  else {
3384                      r0 = 0.5/dp_T0[res];
3385                  }
3386                  Targ = TargOP;
3387
3388                  even = 1;
3389                  FOR_0_LE_i_LT_k
3390                  { TRES_FOINC = r0 * TARG_INC;
3391#if defined(_HIGHER_ORDER_)
3392                    TresOP  = Tres-i;
3393                    TresOP2 = Tres-1;
3394
3395                    x = 0;
3396                    for (j=1;2*j-1<i;j++)
3397                    x += (*TresOP++) * (*TresOP2--);
3398                    x *= 2;
3399                    if (!even)
3400                        x += (*TresOP) * (*TresOP2); /* !!! */
3401                        even = !even;
3402                        *Tres++ -= r0*x;
3403#endif /* _HIGHER_ORDER_ */
3404                      }
3405                    }
3406#endif
3407#endif
3408#endif /* ALL_TOGETHER_AGAIN */
3409                    break;
3410
3411            /*--------------------------------------------------------------------------*/
3412        case gen_quad:                                            /* gen_quad */
3413            arg1 = get_locint_f();
3414                arg2 = get_locint_f();
3415                res  = get_locint_f();
3416
3417#if !defined(_NTIGHT_)
3418                if (get_val_f()!=dp_T0[arg1]) {
3419                    fprintf(DIAG_OUT,
3420                            "ADOL-C Warning: forward sweep aborted; tape invalid!\n");
3421                    IF_KEEP_TAYLOR_CLOSE
3422                    end_sweep();
3423                    return -2;
3424                }
3425#endif /* !_NTIGHT_ */
3426
3427                coval = get_val_f();
3428
3429                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
3430
3431#if !defined(_NTIGHT_)
3432                dp_T0[res] = coval;
3433#endif /* !_NTIGHT_ */
3434
3435#if defined(_INDO_)
3436               fprintf(DIAG_OUT,
3437                    "ADOL-C Warning: forward sweep aborted; sparse mode not available for gen_quad!\n");
3438               end_sweep();
3439               return -2;
3440#else
3441#if !defined(_ZOS_) /* BREAK_ZOS */
3442                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
3443                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
3444#ifdef _INT_FOR_
3445                FOR_0_LE_l_LT_p
3446                TRES_FOINC = TARG1_FOINC;
3447#else
3448                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
3449
3450                FOR_0_LE_l_LT_p
3451                { FOR_0_LE_i_LT_k
3452                  { /* olvo 980922 changed order to allow x = gen_quad(x) */
3453#if defined(_HIGHER_ORDER_)
3454                      zOP      = dp_z+i;
3455                      (*zOP--) = (i+1) * (*Targ1);
3456#endif /* _HIGHER_ORDER_ */
3457
3458                      TRES_FOINC = dp_T0[arg2] * TARG1_INC;
3459
3460#if defined(_HIGHER_ORDER_)
3461                      Targ2OP = Targ2;
3462
3463                      *Tres *= (i+1);
3464                      for (j=0;j<i;j++)
3465                      *Tres += (*Targ2OP++) * (*zOP--);
3466                      *Tres++ /= (i+1);
3467#endif /* _HIGHER_ORDER_ */
3468                  }
3469                  HOV_INC(Targ2, k)
3470                }
3471#endif
3472#endif
3473#endif /* ALL_TOGETHER_AGAIN */
3474                break;
3475
3476            /*--------------------------------------------------------------------------*/
3477        case min_op:                                                /* min_op */
3478            arg1  = get_locint_f();
3479                arg2  = get_locint_f();
3480                res   = get_locint_f();
3481                coval = get_val_f();
3482
3483                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
3484
3485#if !defined(_NTIGHT_)
3486                /* olvo 980923 changed order to allow x = min(x,y) etc. */
3487
3488                /* olvo/mitev 980721 return value (taken from below) */
3489                if (dp_T0[arg1] > dp_T0[arg2]) {
3490                    if (coval)
3491                        MINDEC(ret_c,2);
3492                } else
3493                    if (dp_T0[arg1] < dp_T0[arg2]) {
3494                        if (!coval)
3495                            MINDEC(ret_c,2);
3496                    } else
3497                        if (arg1 != arg2)
3498                            MINDEC(ret_c,1);
3499#endif /* !_NTIGHT_ */
3500
3501#if defined (_INDO_)
3502#if defined (_INDOPRO_)
3503#if defined (_TIGHT_)
3504                    if (dp_T0[arg1] < dp_T0[arg2])
3505                        copy_index_domain(res, arg1, ind_dom);
3506                    else {
3507                        if (dp_T0[arg1] > dp_T0[arg2])
3508                            copy_index_domain(res, arg2, ind_dom);
3509                        else
3510                            combine_2_index_domains(res, arg1, arg2, ind_dom);
3511                    }
3512#else
3513                    combine_2_index_domains(res, arg1, arg2, ind_dom);
3514#endif
3515#endif
3516#if defined(_NONLIND_)
3517#ifdef _TIGHT_
3518                    if (dp_T0[arg1] < dp_T0[arg2])
3519                      {
3520                        fod[opind].entry = maxopind+2;
3521                        fod[opind].left = &fod[arg_index[arg1]];
3522                        fod[opind].right = NULL;
3523                        arg_index[res] = opind++;               
3524                      }           
3525                    else {
3526                        if (dp_T0[arg1] > dp_T0[arg2])
3527                          {
3528                            fod[opind].entry = maxopind+2;
3529                            fod[opind].left = &fod[arg_index[arg2]];
3530                            fod[opind].right = NULL;
3531                            arg_index[res] = opind++;           
3532
3533                          }               
3534                        else
3535                          {
3536                            fod[opind].entry = maxopind+2;
3537                            fod[opind].left = &fod[arg_index[arg1]];
3538                            fod[opind].right = &fod[arg_index[arg2]];
3539                            arg_index[res] = opind++;           
3540                          }
3541                    }
3542#else
3543                    fod[opind].entry = maxopind+2;
3544                    fod[opind].left = &fod[arg_index[arg1]];
3545                    fod[opind].right = &fod[arg_index[arg2]];
3546                    arg_index[res] = opind++;           
3547                    arg_index[res] = opind++;           
3548#endif
3549#endif
3550#else
3551#if !defined(_ZOS_) /* BREAK_ZOS */
3552                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
3553                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
3554                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
3555
3556#ifdef _INT_FOR_
3557#ifdef _TIGHT_
3558                Tqo = NULL;
3559                if (dp_T0[arg1] > dp_T0[arg2])
3560                    Tqo = Targ2;
3561                else
3562                    if (dp_T0[arg1] < dp_T0[arg2])
3563                        Tqo = Targ1;
3564
3565                FOR_0_LE_l_LT_p
3566                { Targ = Tqo;
3567                  if (Targ == NULL) /* e.g. T0[arg1] == T0[arg2] */
3568                    { Targ1OP = Targ1;
3569                      Targ2OP = Targ2;
3570                      if (TARG1 > TARG2)
3571                          Targ = Targ2OP;
3572                      else
3573                          if (TARG1 < TARG2)
3574                              Targ = Targ1OP;
3575                      Targ1++;
3576                      Targ2++;
3577                      if (Targ == NULL) /* e.g. both are equal */
3578                          Targ = Targ1OP;
3579                  }
3580
3581                  TRES_INC = TARG_INC;
3582
3583                  if (Tqo)
3584                    Tqo++;
3585                }
3586
3587                dp_T0[res] = MIN_ADOLC(dp_T0[arg1], dp_T0[arg2]);
3588#endif /* _TIGHT_ */
3589#ifdef _NTIGHT_
3590                TRES_INC = TARG1_INC | TARG2_INC;
3591#endif /* _NTIGHT_ */
3592#else
3593                Tqo = NULL;
3594                if (dp_T0[arg1] > dp_T0[arg2])
3595                    Tqo = Targ2;
3596                else
3597                    if (dp_T0[arg1] < dp_T0[arg2])
3598                        Tqo = Targ1;
3599
3600                FOR_0_LE_l_LT_p
3601                { Targ = Tqo;
3602                  if (Targ == NULL) /* e.g. dp_T0[arg1] ==
3603                                                                                 dp_T0[arg2] */
3604              { Targ1OP = Targ1;
3605                  Targ2OP = Targ2;
3606                  FOR_0_LE_i_LT_k
3607                  { if (TARG1 > TARG2) {
3608                        Targ = Targ2OP;
3609                        VEC_INC(Targ1, k-i)
3610                            VEC_INC(Targ2, k-i)
3611                            BREAK_FOR_I
3612                        } else
3613                            if (TARG1 < TARG2) {
3614                            Targ = Targ1OP;
3615                            VEC_INC(Targ1, k-i)
3616                                VEC_INC(Targ2, k-i)
3617                                BREAK_FOR_I
3618                            }
3619                        Targ1++;
3620                        Targ2++;
3621                      }
3622                      if (Targ == NULL) /* e.g. both are equal */
3623                          Targ = Targ1OP;
3624                  }
3625
3626                  FOR_0_LE_i_LT_k
3627                  TRES_INC = TARG_INC;
3628
3629                  if (Tqo) {
3630                  VEC_INC(Tqo, k)
3631                  }
3632            }
3633#endif
3634#endif
3635#endif /* ALL_TOGETHER_AGAIN */
3636#if !defined(_NTIGHT_)
3637                dp_T0[res] =
3638                    MIN_ADOLC( dp_T0[arg1],
3639                               dp_T0[arg2] );
3640#endif /* !_NTIGHT_ */
3641                break;
3642
3643                /*--------------------------------------------------------------------------*/
3644            case abs_val:                                              /* abs_val */
3645                arg   = get_locint_f();
3646                res   = get_locint_f();
3647                coval = get_val_f();
3648
3649                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
3650
3651#if !defined(_NTIGHT_)
3652                /* olvo 980923 changed order to allow x = min(x,y) etc. */
3653
3654                /* olvo/mitev 980721 ec n3l (taken from below) */
3655                if (dp_T0[arg] < 0.0) {
3656                    if (coval)
3657                        MINDEC(ret_c,2);
3658                } else
3659                    if (dp_T0[arg] > 0.0) {
3660                        if (!coval)
3661                            MINDEC(ret_c,2);
3662                    }
3663#endif /* !_NTIGHT_ */
3664
3665#if defined(_INDO_)
3666#if defined(_INDOPRO_)
3667                copy_index_domain(res, arg, ind_dom);
3668#endif
3669#if defined(_NONLIND_)
3670                arg_index[res] = arg_index[arg];               
3671#endif
3672#else
3673#if !defined(_ZOS_) /* BREAK_ZOS */
3674                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
3675                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
3676
3677#ifdef _INT_FOR_
3678#ifdef _TIGHT_
3679                y = 0.0;
3680                if (dp_T0[arg] != 0.0) {
3681                    if (dp_T0[arg] < 0.0)
3682                        y = -1.0;
3683                    else
3684                        y = 1.0;
3685                }
3686                FOR_0_LE_l_LT_p
3687                { if ((y == 0.0) && (TARG != 0.0))
3688                  MINDEC(ret_c,1);
3689
3690                  TRES_INC = TARG_INC;
3691                }
3692
3693                dp_T0[res] = fabs(dp_T0[arg]);
3694#endif /* _TIGHT_ */
3695#ifdef _NTIGHT_
3696                FOR_0_LE_l_LT_p
3697                TRES_INC = TARG_INC;
3698#endif /* _NTIGHT_ */
3699#else
3700                y = 0.0;
3701                if (dp_T0[arg] != 0.0) {
3702                    if (dp_T0[arg] < 0.0)
3703                        y = -1.0;
3704                    else
3705                        y = 1.0;
3706                }
3707
3708                FOR_0_LE_l_LT_p
3709                { x = y;
3710                  FOR_0_LE_i_LT_k
3711                  { if ((x == 0.0) && (TARG != 0.0)) {
3712                    MINDEC(ret_c,1);
3713                        if (TARG < 0.0)
3714                            x = -1.0;
3715                        else
3716                            x = 1.0;
3717                    }
3718                    TRES_INC = x * TARG_INC;
3719              }
3720            }
3721#endif
3722#endif
3723#endif /* ALL_TOGETHER_AGAIN */
3724#if !defined(_NTIGHT_)
3725                dp_T0[res] = fabs(dp_T0[arg]);
3726#endif /* !_NTIGHT_ */
3727                break;
3728
3729                /*--------------------------------------------------------------------------*/
3730            case ceil_op:                                              /* ceil_op */
3731                arg   = get_locint_f();
3732                res   = get_locint_f();
3733                coval = get_val_f();
3734
3735                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
3736
3737#if !defined(_NTIGHT_)
3738                dp_T0[res]=ceil(dp_T0[arg]);
3739                /* olvo/mitev 980721 ec n2l (taken from below) */
3740                if (coval != dp_T0[res])
3741                    MINDEC(ret_c,2);
3742#endif /* !_NTIGHT_ */
3743
3744#if defined(_INDO_)
3745#if defined(_INDOPRO_)
3746#ifdef _TIGHT_
3747                ind_dom[res][0] = 0;
3748#else
3749                copy_index_domain(res, arg, ind_dom);
3750#endif /* _TIGHT_ */
3751#endif
3752#if defined(_NONLIND_)
3753                arg_index[res] = arg_index[arg];               
3754#endif
3755#else
3756#if !defined(_ZOS_) /* BREAK_ZOS */
3757                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
3758
3759                FOR_0_LE_l_LT_pk
3760                TRES_INC = 0.0;
3761#endif
3762#endif /* ALL_TOGETHER_AGAIN */
3763                break;
3764
3765                /*--------------------------------------------------------------------------*/
3766            case floor_op:                 /* Compute ceil of adouble    floor_op */
3767                arg   = get_locint_f();
3768                res   = get_locint_f();
3769                coval = get_val_f();
3770
3771                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
3772
3773#if !defined(_NTIGHT_)
3774                dp_T0[res] = floor(dp_T0[arg]);
3775                /* olvo/mitev 980721 ec n2l (taken from below) */
3776                if (coval != dp_T0[res])
3777                    MINDEC(ret_c,2);
3778#endif /* !_NTIGHT_ */
3779
3780#if defined(_INDO_)
3781#if defined(_INDOPRO_)
3782#ifdef _TIGHT_
3783                ind_dom[res][0] = 0;
3784#else
3785                copy_index_domain(res, arg, ind_dom);
3786#endif /* _TIGHT_ */
3787#endif
3788#if defined(_NONLIND_)
3789                arg_index[res] = arg_index[arg];               
3790#endif
3791#else
3792#if !defined(_ZOS_) /* BREAK_ZOS */
3793                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
3794
3795                FOR_0_LE_l_LT_pk
3796                TRES_INC = 0.0;
3797#endif
3798#endif /* ALL_TOGETHER_AGAIN */
3799                break;
3800
3801
3802                /****************************************************************************/
3803                /*                                                             CONDITIONALS */
3804
3805                /*--------------------------------------------------------------------------*/
3806            case cond_assign:                                      /* cond_assign */
3807                arg   = get_locint_f();
3808                arg1  = get_locint_f();
3809                arg2  = get_locint_f();
3810                res   = get_locint_f();
3811                coval = get_val_f();
3812
3813                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
3814
3815                /* olvo 980924 changed order to allow reflexive ops */
3816#if defined (_INDO_)
3817#if defined (_INDOPRO_)
3818#if defined (_TIGHT_)
3819                if (dp_T0[arg] > 0) {
3820                    if (coval <= 0.0)
3821                        MINDEC(ret_c,2);
3822                    dp_T0[res] = dp_T0[arg1];
3823
3824                    copy_index_domain(res, arg1, ind_dom);
3825
3826                } else {
3827                    if (coval > 0.0)
3828                        MINDEC(ret_c,2);
3829                    if (dp_T0[arg] == 0)
3830                        MINDEC(ret_c,0);
3831                    dp_T0[res] = dp_T0[arg2];
3832                    copy_index_domain(res, arg2, ind_dom);
3833                }
3834#else
3835                    combine_2_index_domains(res, arg1, arg2, ind_dom);
3836#endif
3837#endif
3838#if defined (_NONLIND_)
3839#ifdef _TIGHT_
3840                if (dp_T0[arg] > 0) {
3841                    if (coval <= 0.0)
3842                        MINDEC(ret_c,2);
3843                    dp_T0[res] = dp_T0[arg1];
3844
3845                    arg_index[res] = arg_index[arg1];           
3846                } else {
3847                    if (coval > 0.0)
3848                        MINDEC(ret_c,2);
3849                    if (dp_T0[arg] == 0)
3850                        MINDEC(ret_c,0);
3851                    dp_T0[res] = dp_T0[arg2];
3852
3853                    arg_index[res] = arg_index[arg2];           
3854                }
3855
3856#else
3857               arg_index[res] = opind++;               
3858#endif
3859#endif
3860#else
3861#if !defined(_ZOS_) /* BREAK_ZOS */
3862                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
3863                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
3864                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
3865#endif /* ALL_TOGETHER_AGAIN */
3866
3867#ifdef _INT_FOR_
3868#ifdef _TIGHT_
3869                coval = get_val_f();
3870
3871                if (dp_T0[arg] > 0)
3872                    FOR_0_LE_l_LT_pk
3873                    TRES_INC = TARG1_INC;
3874                else
3875                    FOR_0_LE_l_LT_pk
3876                    TRES_INC = TARG2_INC;
3877
3878                if (dp_T0[arg] > 0) {
3879                    if (coval <= 0.0)
3880                        MINDEC(ret_c,2);
3881                    dp_T0[res] = dp_T0[arg1];
3882                } else {
3883                    if (coval > 0.0)
3884                        MINDEC(ret_c,2);
3885                    if (dp_T0[arg] == 0)
3886                        MINDEC(ret_c,0);
3887                    dp_T0[res] = dp_T0[arg2];
3888                }
3889#endif /* _TIGHT_ */
3890#ifdef _NTIGHT_
3891                FOR_0_LE_l_LT_pk
3892                TRES_INC = TARG1_INC | TARG2_INC;
3893#endif /* _NTIGHT_ */
3894#else
3895#if !defined(_ZOS_) /* BREAK_ZOS */
3896                if (dp_T0[arg] > 0)
3897                    FOR_0_LE_l_LT_pk
3898                    TRES_INC = TARG1_INC;
3899                else
3900                    FOR_0_LE_l_LT_pk
3901                    TRES_INC = TARG2_INC;
3902#endif
3903
3904                if (dp_T0[arg] > 0) {
3905                    if (coval <= 0.0)
3906                        MINDEC(ret_c,2);
3907                    dp_T0[res] = dp_T0[arg1];
3908                } else {
3909                    if (coval > 0.0)
3910                        MINDEC(ret_c,2);
3911                    if (dp_T0[arg] == 0)
3912                        MINDEC(ret_c,0);
3913                    dp_T0[res] = dp_T0[arg2];
3914                }
3915#endif
3916#endif /* ALL_TOGETHER_AGAIN */
3917                break;
3918
3919                /*--------------------------------------------------------------------------*/
3920            case cond_assign_s:                                  /* cond_assign_s */
3921                arg   = get_locint_f();
3922                arg1  = get_locint_f();
3923                res   = get_locint_f();
3924                coval = get_val_f();
3925
3926                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
3927
3928                /* olvo 980924 changed order to allow reflexive ops */
3929#if defined(_INDO_)
3930#if defined(_INDOPRO_)
3931#ifdef _TIGHT_
3932                if (dp_T0[arg] > 0)
3933                    copy_index_domain(res, arg1, ind_dom);
3934#else
3935                merge_2_index_domains(res, arg1, ind_dom);
3936#endif
3937#endif
3938#if defined(_NONLIND_)
3939                arg_index[res] = arg_index[arg1];               
3940#endif
3941#else
3942#if !defined(_ZOS_) /* BREAK_ZOS */
3943                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
3944                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
3945#endif /* ALL_TOGETHER_AGAIN */
3946
3947#ifdef _INT_FOR_
3948#ifdef _TIGHT_
3949                coval = get_val_f();
3950
3951                if (dp_T0[arg] > 0)
3952#endif /* _TIGHT_ */
3953                    FOR_0_LE_l_LT_pk
3954                    TRES_INC = TARG1_INC;
3955
3956#ifdef _TIGHT_
3957                if (dp_T0[arg] > 0) {
3958                    if (coval <= 0.0)
3959                        MINDEC(ret_c,2);
3960                    dp_T0[res] = dp_T0[arg1];
3961                } else
3962                    if (dp_T0[arg] == 0)
3963                        MINDEC(ret_c,0);
3964#endif /* _TIGHT_ */
3965#else
3966#if !defined(_ZOS_) /* BREAK_ZOS */
3967                if (dp_T0[arg] > 0)
3968                    FOR_0_LE_l_LT_pk
3969                    TRES_INC = TARG1_INC;
3970#endif
3971                if (dp_T0[arg] > 0) {
3972                    if (coval <= 0.0)
3973                        MINDEC(ret_c,2);
3974                    dp_T0[res] = dp_T0[arg1];
3975                } else
3976                    if (dp_T0[arg] == 0)
3977                        MINDEC(ret_c,0);
3978#endif
3979#endif /* ALL_TOGETHER_AGAIN */
3980                break;
3981
3982
3983                /*--------------------------------------------------------------------------*/
3984                /* NEW CONDITIONALS */
3985                /*--------------------------------------------------------------------------*/
3986            case neq_a_a:
3987            case eq_a_a:
3988            case le_a_a:
3989            case ge_a_a:
3990            case lt_a_a:
3991            case gt_a_a:
3992                coval = get_val_f();
3993                arg = get_locint_f();
3994                arg1 = get_locint_f();
3995                res = get_locint_f();
3996#if !defined(_NTIGHT_)
3997                {
3998                    revreal retval = -1;
3999                    const char* opname = "";
4000                    switch (operation) {
4001                    case neq_a_a:
4002                        retval = (revreal)(dp_T0[arg] != dp_T0[arg1]);
4003                        opname = "neq_a_a";
4004                        break;
4005                    case eq_a_a:
4006                        retval = (revreal)(dp_T0[arg] == dp_T0[arg1]);
4007                        opname = "eq_a_a";
4008                        break;
4009                    case ge_a_a:
4010                        retval = (revreal)(dp_T0[arg] >= dp_T0[arg1]);
4011                        opname = "ge_a_a";
4012                        break;
4013                    case le_a_a:
4014                        retval = (revreal)(dp_T0[arg] <= dp_T0[arg1]);
4015                        opname = "le_a_a";
4016                        break;
4017                    case gt_a_a:
4018                        retval = (revreal)(dp_T0[arg] > dp_T0[arg1]);
4019                        opname = "gt_a_a";
4020                        break;
4021                    case lt_a_a:
4022                        retval = (revreal)(dp_T0[arg] < dp_T0[arg1]);
4023                        opname = "lt_a_a";
4024                        break;
4025                    }
4026                    if (retval != coval && ADOLC_GLOBAL_TAPE_VARS.branchSwitchWarning)
4027                        fprintf(DIAG_OUT,
4028                                "ADOL-C Warning: Branch switch detected in comparison "
4029                                "(operator %s).\n"
4030                                "Results may be unpredictable! Retaping recommended!\n",opname);
4031                    dp_T0[res] = retval;
4032                }
4033#endif
4034#if defined(_INDO_)
4035#if defined(_INDOPRO_)
4036                ind_dom[res][0]=0;
4037#endif
4038#if defined(_NONLIND_)
4039                fod[opind].entry = maxopind+2;
4040                fod[opind].left = NULL;
4041                fod[opind].right = NULL;
4042                arg_index[res] = opind++;
4043#endif
4044#else
4045#if !defined(_ZOS_) /* BREAK_ZOS */
4046                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
4047
4048                FOR_0_LE_l_LT_pk
4049                TRES_INC = 0;
4050#endif
4051#endif /* ALL_TOGETHER_AGAIN */
4052
4053                break;
4054
4055                /*--------------------------------------------------------------------------*/
4056            case subscript:
4057                coval = get_val_f();
4058                arg = get_locint_f();
4059                res = get_locint_f();
4060                {
4061                    size_t cnt, idx, numvar = (size_t)trunc(fabs(coval));
4062                    locint vectorloc;
4063                    vectorloc = get_locint_f();
4064#if !defined(_NTIGHT_)
4065                    idx = (size_t)trunc(fabs(dp_T0[arg]));
4066                    if (idx >= numvar)
4067                        fprintf(DIAG_OUT, "ADOL-C warning: index out of bounds while subscripting n=%z, idx=%z\n", numvar, idx);
4068                    arg1 = vectorloc+idx;
4069                    IF_KEEP_WRITE_TAYLOR(res,keep,k,p);
4070                    dp_T0[res] = dp_T0[arg1];
4071#if defined(_INDO_)
4072#if defined(_INDOPRO_)
4073                    copy_index_domain(res, arg1, ind_dom);
4074#endif
4075#if defined(_NONLIND_)
4076                    arg_index[res] = arg_index[arg1];
4077#endif
4078#else
4079#if !defined(_ZOS_) /* BREAK_ZOS */
4080                    ASSIGN_T(Targ1,TAYLOR_BUFFER[arg1])
4081                    ASSIGN_T(Tres,TAYLOR_BUFFER[res])
4082
4083                    FOR_0_LE_l_LT_pk
4084                    TRES_INC = TARG1_INC;
4085#endif
4086#endif
4087#else
4088                    fprintf(DIAG_OUT, "ADOL-C error: active subscripting does not work in safe mode, please use tight mode\n");
4089                    exit(-2);
4090#endif /* ALL_TOGETHER_AGAIN */
4091                }
4092                break;
4093
4094            case subscript_ref:
4095                coval = get_val_f();
4096                arg = get_locint_f();
4097                res = get_locint_f();
4098                {
4099                    size_t cnt, idx, numvar = (size_t)trunc(fabs(coval));
4100                    locint vectorloc;
4101                    vectorloc = get_locint_f();
4102#if !defined(_NTIGHT_)
4103                    idx = (size_t)trunc(fabs(dp_T0[arg]));
4104                    if (idx >= numvar)
4105                        fprintf(DIAG_OUT, "ADOL-C warning: index out of bounds while subscripting (ref) n=%z, idx=%z\n", numvar, idx);
4106                    arg1 = vectorloc+idx;
4107                    IF_KEEP_WRITE_TAYLOR(res,keep,k,p);
4108                    dp_T0[res] = arg1;
4109#else
4110                    fprintf(DIAG_OUT, "ADOL-C error: active subscripting does not work in safe mode, please use tight mode\n");
4111                    exit(-2);
4112#endif
4113                }
4114                break;
4115
4116            case ref_copyout:
4117                arg = get_locint_f();
4118                res = get_locint_f();
4119#if !defined(_NTIGHT_)
4120                arg1 = (size_t)trunc(fabs(dp_T0[arg]));
4121                IF_KEEP_WRITE_TAYLOR(res,keep,k,p);
4122                dp_T0[res] = dp_T0[arg1];
4123#if defined(_INDO_)
4124#if defined(_INDOPRO_)
4125                copy_index_domain(res, arg1, ind_dom);
4126#endif
4127#if defined(_NONLIND_)
4128                arg_index[res] = arg_index[arg1];
4129#endif
4130#else
4131#if !defined(_ZOS_) /* BREAK_ZOS */
4132                ASSIGN_T(Targ1,TAYLOR_BUFFER[arg1])
4133                ASSIGN_T(Tres,TAYLOR_BUFFER[res])
4134
4135                FOR_0_LE_l_LT_pk
4136                TRES_INC = TARG1_INC;
4137#endif
4138#endif
4139#else
4140                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4141                exit(-2);
4142#endif /* ALL_TOGETHER_AGAIN */
4143                break;
4144
4145            case ref_incr_a:
4146                arg = get_locint_f();
4147#if !defined(_NTIGHT_)
4148                arg1 = (size_t)trunc(fabs(dp_T0[arg]));
4149                IF_KEEP_WRITE_TAYLOR(arg1,keep,k,p);
4150                dp_T0[arg1]++;
4151#else
4152                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4153                exit(-2);
4154#endif
4155                break;
4156
4157            case ref_decr_a:
4158                arg = get_locint_f();
4159#if !defined(_NTIGHT_)
4160                arg1 = (size_t)trunc(fabs(dp_T0[arg]));
4161                IF_KEEP_WRITE_TAYLOR(arg1,keep,k,p);
4162                dp_T0[arg1]--;
4163#else
4164                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4165                exit(-2);
4166#endif
4167                break;
4168
4169            case ref_assign_d:
4170                arg = get_locint_f();
4171                coval = get_val_f();
4172               
4173#if !defined(_NTIGHT_)
4174                arg1 = (size_t)trunc(fabs(dp_T0[arg]));
4175                IF_KEEP_WRITE_TAYLOR(arg1,keep,k,p)
4176                dp_T0[arg1] = coval;
4177#if defined(_INDO_)
4178#if defined(_INDOPRO_)
4179                ind_dom[arg1][0] = 0;
4180#endif
4181#if defined(_NONLIND_)
4182                fod[opind].entry = maxopind+2;
4183                fod[opind].left = NULL;
4184                fod[opind].right = NULL;
4185                arg_index[arg1] = opind++;
4186#endif
4187#else
4188#if !defined(_ZOS_)
4189                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
4190
4191                FOR_0_LE_l_LT_pk
4192                TARG1_INC = 0;
4193#endif
4194#endif
4195#else
4196                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4197                exit(-2);
4198#endif
4199                break;
4200
4201            case ref_assign_d_zero:
4202                arg = get_locint_f();
4203
4204#if !defined(_NTIGHT_)
4205                arg1 = (size_t)trunc(fabs(dp_T0[arg]));
4206                IF_KEEP_WRITE_TAYLOR(arg1,keep,k,p)
4207                dp_T0[arg1] = 0.0;
4208#if defined(_INDO_)
4209#if defined(_INDOPRO_)
4210                ind_dom[arg1][0] = 0;
4211#endif
4212#if defined(_NONLIND_)
4213                fod[opind].entry = maxopind+2;
4214                fod[opind].left = NULL;
4215                fod[opind].right = NULL;
4216                arg_index[arg1] = opind++;
4217#endif
4218#else
4219#if !defined(_ZOS_)
4220                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
4221
4222                FOR_0_LE_l_LT_pk
4223                TARG1_INC = 0;
4224#endif
4225#endif
4226#else
4227                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4228                exit(-2);
4229#endif
4230                break;
4231
4232            case ref_assign_d_one:
4233                arg = get_locint_f();
4234
4235#if !defined(_NTIGHT_)
4236                arg1 = (size_t)trunc(fabs(dp_T0[arg]));
4237                IF_KEEP_WRITE_TAYLOR(arg1,keep,k,p)
4238                dp_T0[arg1] = 1.0;
4239#if defined(_INDO_)
4240#if defined(_INDOPRO_)
4241                ind_dom[arg1][0] = 0;
4242#endif
4243#if defined(_NONLIND_)
4244                fod[opind].entry = maxopind+2;
4245                fod[opind].left = NULL;
4246                fod[opind].right = NULL;
4247                arg_index[arg1] = opind++;
4248#endif
4249#else
4250#if !defined(_ZOS_)
4251                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
4252
4253                FOR_0_LE_l_LT_pk
4254                TARG1_INC = 0;
4255#endif
4256#endif
4257#else
4258                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4259                exit(-2);
4260#endif
4261                break;
4262
4263            case ref_assign_a:           /* assign an adouble variable an    assign_a */
4264                /* adouble value. (=) */
4265                arg = get_locint_f();
4266                res = get_locint_f();
4267
4268#if !defined(_NTIGHT_)
4269                arg1 = (size_t)trunc(fabs(dp_T0[res]));
4270                IF_KEEP_WRITE_TAYLOR(arg1,keep,k,p)
4271                dp_T0[arg1] = dp_T0[arg];
4272#if defined(_INDO_)
4273#if defined(_INDOPRO_)
4274                copy_index_domain(arg1, arg, ind_dom);
4275#endif
4276#if defined(_NONLIND_)
4277                arg_index[arg1] = arg_index[arg];
4278#endif
4279#else
4280#if !defined(_ZOS_) /* BREAK_ZOS */
4281                ASSIGN_T(Targ,TAYLOR_BUFFER[arg])
4282                ASSIGN_T(Targ1,TAYLOR_BUFFER[arg1])
4283
4284                FOR_0_LE_l_LT_pk
4285                TARG1_INC = TARG_INC;
4286#endif
4287#endif
4288#else
4289                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4290                exit(-2);
4291#endif /* ALL_TOGETHER_AGAIN */
4292                break;
4293
4294            case ref_assign_ind:       /* assign an adouble variable an    assign_ind */
4295                /* independent double value (<<=) */
4296                arg = get_locint_f();
4297
4298
4299#if !defined(_NTIGHT_)
4300                res = (size_t)trunc(fabs(dp_T0[arg]));
4301                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4302                dp_T0[res] = basepoint[indexi];
4303#if defined(_INDO_)
4304#if defined(_INDOPRO_)
4305                ind_dom[res][0] = 1;
4306                ind_dom[res][2] = indexi;
4307#endif
4308#if defined(_NONLIND_)
4309                fod[opind].entry = indexi;
4310                fod[opind].left = NULL;
4311                fod[opind].right = NULL;
4312                arg_index[res] = opind++;
4313#endif
4314#else
4315#if !defined(_ZOS_) /* BREAK_ZOS */
4316                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
4317
4318#ifdef _INT_FOR_
4319                FOR_0_LE_l_LT_p
4320                TRES_INC = ARGUMENT(indexi,l,i);
4321#else
4322                FOR_0_LE_l_LT_p
4323                FOR_0_LE_i_LT_k
4324                TRES_INC = ARGUMENT(indexi,l,i);
4325#endif
4326#endif
4327#endif
4328#else
4329                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4330                exit(-2);
4331#endif /* ALL_TOGETHER_AGAIN */
4332                ++indexi;
4333                break;
4334
4335            case ref_eq_plus_d:            /* Add a floating point to an    eq_plus_d */
4336                /* adouble. (+=) */
4337                arg  = get_locint_f();
4338                coval = get_val_f();
4339
4340
4341#if !defined(_NTIGHT_)
4342                res = (size_t)trunc(fabs(dp_T0[arg]));
4343                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4344                dp_T0[res] += coval;
4345#else
4346                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4347                exit(-2);
4348#endif /* !_NTIGHT_ */
4349                break;
4350
4351                /*--------------------------------------------------------------------------*/
4352            case ref_eq_plus_a:             /* Add an adouble to another    eq_plus_a */
4353                /* adouble. (+=) */
4354                arg = get_locint_f();
4355                arg1 = get_locint_f();
4356
4357#if !defined(_NTIGHT_)
4358                res = (size_t)trunc(fabs(dp_T0[arg1]));
4359                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4360                dp_T0[res] += dp_T0[arg];
4361#if defined(_INDO_)
4362#if defined(_INDOPRO_)
4363                merge_2_index_domains(res, arg, ind_dom);
4364#endif
4365#if defined(_NONLIND_)
4366                fod[opind].entry = maxopind+2;
4367                fod[opind].left = &fod[arg_index[res]];
4368                fod[opind].right = &fod[arg_index[arg]];
4369                arg_index[res] = opind++;
4370#endif
4371#else
4372#if !defined(_ZOS_) /* BREAK_ZOS */
4373                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
4374                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
4375
4376#ifdef _INT_FOR_
4377                FOR_0_LE_l_LT_pk
4378                TRES_INC |= TARG_INC;
4379#else
4380                FOR_0_LE_l_LT_pk
4381                TRES_INC += TARG_INC;
4382#endif
4383#endif
4384#endif
4385#else
4386                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4387                exit(-2);
4388#endif /* ALL_TOGETHER_AGAIN */
4389                break;
4390
4391            case ref_eq_min_d:       /* Subtract a floating point from an    eq_min_d */
4392                /* adouble. (-=) */
4393                arg = get_locint_f();
4394                coval = get_val_f();
4395
4396#if !defined(_NTIGHT_)
4397                res = (size_t)trunc(fabs(dp_T0[arg]));
4398                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4399                dp_T0[res] -= coval;
4400#else
4401                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4402                exit(-2);
4403#endif /* !_NTIGHT_ */
4404                break;
4405
4406                /*--------------------------------------------------------------------------*/
4407            case ref_eq_min_a:        /* Subtract an adouble from another    eq_min_a */
4408                /* adouble. (-=) */
4409                arg = get_locint_f();
4410                arg1 = get_locint_f();
4411
4412#if !defined(_NTIGHT_)
4413                res = (size_t)trunc(fabs(dp_T0[arg1]));
4414                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4415                dp_T0[res] -= dp_T0[arg];
4416#if defined(_INDO_)
4417#if defined(_INDOPRO_)
4418                merge_2_index_domains(res, arg, ind_dom);
4419#endif
4420#if defined(_NONLIND_)
4421                fod[opind].entry = maxopind+2;
4422                fod[opind].left = &fod[arg_index[res]];
4423                fod[opind].right = &fod[arg_index[arg]];
4424                arg_index[res] = opind++;
4425#endif
4426#else
4427#if !defined(_ZOS_) /* BREAK_ZOS */
4428                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
4429                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
4430
4431#ifdef _INT_FOR_
4432                FOR_0_LE_l_LT_pk
4433                TRES_INC |= TARG_INC;
4434#else
4435                FOR_0_LE_l_LT_pk
4436                TRES_INC -= TARG_INC;
4437#endif
4438#endif
4439#endif
4440#else
4441                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4442                exit(-2);
4443#endif /* ALL_TOGETHER_AGAIN */
4444                break;
4445
4446            case ref_eq_mult_d:              /* Multiply an adouble by a    eq_mult_d */
4447                /* flaoting point. (*=) */
4448                arg = get_locint_f();
4449                coval = get_val_f();
4450
4451#if !defined(_NTIGHT_)
4452                res = (size_t)trunc(fabs(dp_T0[arg]));
4453                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4454                dp_T0[res] *= coval;
4455#if !defined(_INDO_)
4456#if !defined(_ZOS_) /* BREAK_ZOS */
4457#if !defined( _INT_FOR_)
4458
4459                FOR_0_LE_l_LT_pk
4460                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
4461
4462                FOR_0_LE_l_LT_pk
4463                TRES_INC *= coval;
4464#endif
4465#endif
4466#endif
4467#else
4468                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4469                exit(-2);
4470#endif /* ALL_TOGETHER_AGAIN */
4471                break;
4472
4473            case ref_eq_mult_a:       /* Multiply one adouble by another    eq_mult_a */
4474                /* (*=) */
4475                arg = get_locint_f();
4476                arg1 = get_locint_f();
4477
4478#if !defined(_NTIGHT_)
4479                res = (size_t)trunc(fabs(dp_T0[arg1]));
4480                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4481#if defined(_INDO_)
4482#if defined(_INDOPRO_)
4483                merge_2_index_domains(res, arg, ind_dom);
4484#endif
4485#if defined(_NONLIND_)
4486                fod[opind].entry = maxopind+2;
4487                fod[opind].left = &fod[arg_index[res]];
4488                fod[opind].right = &fod[arg_index[arg]];
4489                traverse_unary(&fod[arg_index[res]], nonl_dom, &fod[arg_index[arg]], indcheck+1,maxopind+2);
4490                traverse_unary(&fod[arg_index[arg]], nonl_dom, &fod[arg_index[res]], indcheck+1,maxopind+2);
4491                arg_index[res] = opind++;
4492#endif
4493#if defined(_NONLIND_OLD_)
4494                extend_nonlinearity_domain_binary(res, arg, ind_dom, nonl_dom);
4495#endif
4496#else
4497#if !defined(_ZOS_) /* BREAK_ZOS */
4498                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
4499                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
4500
4501                INC_pk_1(Tres)
4502                INC_pk_1(Targ)
4503
4504#ifdef _INT_FOR_
4505                FOR_p_GT_l_GE_0
4506                TRES_FODEC |= TARG_DEC;
4507#else
4508                FOR_p_GT_l_GE_0
4509                FOR_k_GT_i_GE_0
4510                { TRES_FODEC = dp_T0[res]*TARG_DEC +
4511                               TRES*dp_T0[arg];
4512                  DEC_TRES_FO
4513#ifdef _HIGHER_ORDER_
4514                  TresOP = Tres-i;
4515                  TargOP = Targ;
4516
4517                  for (j=0;j<i;j++)
4518                  *Tres += (*TresOP++) * (*TargOP--);
4519                  Tres--;
4520#endif /* _HIGHER_ORDER_ */
4521                }
4522#endif
4523#endif
4524#endif
4525                dp_T0[res] *= dp_T0[arg];
4526#else
4527                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4528                exit(-2);
4529#endif /* !_NTIGHT_ */
4530                break;
4531
4532            case ref_cond_assign:                                      /* cond_assign */
4533                arg   = get_locint_f();
4534                arg1  = get_locint_f();
4535                arg2  = get_locint_f();
4536                { 
4537                    locint ref = get_locint_f();
4538                    coval = get_val_f();
4539#if !defined(_NTIGHT_)
4540                    res   = (size_t)trunc(fabs(dp_T0[ref]));
4541
4542                    IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4543
4544                /* olvo 980924 changed order to allow reflexive ops */
4545#if defined(_INDO_)
4546                if (dp_T0[arg] > 0) {
4547                    if (coval <= 0.0)
4548                        MINDEC(ret_c,2);
4549                    dp_T0[res] = dp_T0[arg1];
4550
4551#if defined(_INDOPRO_)
4552                    copy_index_domain(res, arg1, ind_dom);
4553#endif
4554#if defined(_NONLIND_)
4555                    arg_index[res] = arg_index[arg1];
4556#endif
4557                } else {
4558                    if (coval > 0.0)
4559                        MINDEC(ret_c,2);
4560                    if (dp_T0[arg] == 0)
4561                        MINDEC(ret_c,0);
4562                    dp_T0[res] = dp_T0[arg2];
4563
4564#if defined(_INDOPRO_)
4565                    copy_index_domain(res, arg2, ind_dom);
4566#endif
4567#if defined(_NONLIND_)
4568                    arg_index[res] = arg_index[arg2];
4569#endif
4570                }
4571#else
4572#if !defined(_ZOS_) /* BREAK_ZOS */
4573                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
4574                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
4575                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
4576#endif /* ALL_TOGETHER_AGAIN */
4577
4578#ifdef _INT_FOR_
4579                coval = get_val_f();
4580
4581                if (dp_T0[arg] > 0)
4582                    FOR_0_LE_l_LT_pk
4583                    TRES_INC = TARG1_INC;
4584                else
4585                    FOR_0_LE_l_LT_pk
4586                    TRES_INC = TARG2_INC;
4587
4588                if (dp_T0[arg] > 0) {
4589                    if (coval <= 0.0)
4590                        MINDEC(ret_c,2);
4591                    dp_T0[res] = dp_T0[arg1];
4592                } else {
4593                    if (coval > 0.0)
4594                        MINDEC(ret_c,2);
4595                    if (dp_T0[arg] == 0)
4596                        MINDEC(ret_c,0);
4597                    dp_T0[res] = dp_T0[arg2];
4598                }
4599                FOR_0_LE_l_LT_pk
4600                TRES_INC = TARG1_INC | TARG2_INC;
4601#else
4602#if !defined(_ZOS_) /* BREAK_ZOS */
4603                if (dp_T0[arg] > 0)
4604                    FOR_0_LE_l_LT_pk
4605                    TRES_INC = TARG1_INC;
4606                else
4607                    FOR_0_LE_l_LT_pk
4608                    TRES_INC = TARG2_INC;
4609#endif
4610
4611                if (dp_T0[arg] > 0) {
4612                    if (coval <= 0.0)
4613                        MINDEC(ret_c,2);
4614                    dp_T0[res] = dp_T0[arg1];
4615                } else {
4616                    if (coval > 0.0)
4617                        MINDEC(ret_c,2);
4618                    if (dp_T0[arg] == 0)
4619                        MINDEC(ret_c,0);
4620                    dp_T0[res] = dp_T0[arg2];
4621                }
4622#endif
4623#endif
4624#else
4625                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4626                exit(-2);
4627#endif /* ALL_TOGETHER_AGAIN */
4628                }
4629                break;
4630
4631            case ref_cond_assign_s:                                  /* cond_assign_s */
4632                arg   = get_locint_f();
4633                arg1  = get_locint_f();
4634                arg2   = get_locint_f();
4635                coval = get_val_f();
4636
4637#if !defined(_NTIGHT_)
4638                res = (size_t)trunc(fabs(dp_T0[arg2]));
4639                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4640
4641                /* olvo 980924 changed order to allow reflexive ops */
4642#if defined(_INDO_)
4643                if (dp_T0[arg] > 0) {
4644#if defined(_INDOPRO_)
4645                    copy_index_domain(res, arg1, ind_dom);
4646#endif
4647#if defined(_NONLIND_)
4648                    arg_index[res] = arg_index[arg1];
4649#endif
4650                }
4651#else
4652#if !defined(_ZOS_) /* BREAK_ZOS */
4653                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
4654                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
4655#endif /* ALL_TOGETHER_AGAIN */
4656
4657#ifdef _INT_FOR_
4658                coval = get_val_f();
4659
4660                if (dp_T0[arg] > 0)
4661                    FOR_0_LE_l_LT_pk
4662                    TRES_INC = TARG1_INC;
4663
4664                if (dp_T0[arg] > 0) {
4665                    if (coval <= 0.0)
4666                        MINDEC(ret_c,2);
4667                    dp_T0[res] = dp_T0[arg1];
4668                } else
4669                    if (dp_T0[arg] == 0)
4670                        MINDEC(ret_c,0);
4671#else
4672#if !defined(_ZOS_) /* BREAK_ZOS */
4673                if (dp_T0[arg] > 0)
4674                    FOR_0_LE_l_LT_pk
4675                    TRES_INC = TARG1_INC;
4676#endif
4677                if (dp_T0[arg] > 0) {
4678                    if (coval <= 0.0)
4679                        MINDEC(ret_c,2);
4680                    dp_T0[res] = dp_T0[arg1];
4681                } else
4682                    if (dp_T0[arg] == 0)
4683                        MINDEC(ret_c,0);
4684#endif
4685#endif
4686#else
4687                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4688                exit(-2);
4689#endif /* ALL_TOGETHER_AGAIN */
4690                break;
4691
4692                /****************************************************************************/
4693                /*                                                          REMAINING STUFF */
4694
4695                /*--------------------------------------------------------------------------*/
4696            case take_stock_op:                                  /* take_stock_op */
4697                size = get_locint_f();
4698                res  = get_locint_f();
4699                d    = get_val_v_f(size);
4700
4701                for (ls=0;ls<size;ls++) {
4702#if !defined(_NTIGHT_)
4703                    dp_T0[res]=*d;
4704#endif /* !_NTIGHT_ */
4705#if !defined(_INDO_)
4706#if !defined(_ZOS_) /* BREAK_ZOS */
4707                    ASSIGN_T(Tres,TAYLOR_BUFFER[res])
4708
4709                    FOR_0_LE_l_LT_pk
4710                    TRES_INC = 0;
4711
4712#endif /* ALL_TOGETHER_AGAIN */
4713                    res++;
4714#if !defined(_NTIGHT_)
4715                    d++;
4716#endif /* !_NTIGHT_ */
4717#endif
4718                }
4719                break;
4720
4721                /*--------------------------------------------------------------------------*/
4722            case death_not:                                          /* death_not */
4723                arg1=get_locint_f();
4724                arg2=get_locint_f();
4725
4726#ifdef _KEEP_
4727                if (keep) {
4728                    do {
4729                        IF_KEEP_WRITE_TAYLOR(arg2,keep,k,p)
4730                    } while(arg1 < arg2-- );
4731                }
4732#endif
4733                break;
4734
4735                /*--------------------------------------------------------------------------*/
4736#if defined(_EXTERN_) /* ZOS,  FOS, FOV up to now */
4737            case ext_diff:                       /* extern differntiated function */
4738                ADOLC_CURRENT_TAPE_INFOS.ext_diff_fct_index=get_locint_f();
4739                n=get_locint_f();
4740                m=get_locint_f();
4741                ADOLC_CURRENT_TAPE_INFOS.lowestXLoc_for = get_locint_f();
4742                ADOLC_CURRENT_TAPE_INFOS.lowestYLoc_for = get_locint_f();
4743                ADOLC_CURRENT_TAPE_INFOS.cpIndex = get_locint_f();
4744                edfct=get_ext_diff_fct(ADOLC_CURRENT_TAPE_INFOS.ext_diff_fct_index);
4745
4746                if (edfct->ADOLC_EXT_FCT_POINTER==NULL)
4747                    fail(ADOLC_EXT_DIFF_NULLPOINTER_DIFFFUNC);
4748                if (n>0) {
4749                    if (edfct->dp_x==NULL) fail(ADOLC_EXT_DIFF_NULLPOINTER_ARGUMENT);
4750#if !defined(_ZOS_)
4751                    if (ADOLC_EXT_POINTER_X==NULL) fail(ADOLC_EXT_DIFF_NULLPOINTER_ARGUMENT);
4752#endif
4753                }
4754                if (m>0) {
4755                    if (edfct->dp_y==NULL) fail(ADOLC_EXT_DIFF_NULLPOINTER_ARGUMENT);
4756#if !defined(_ZOS_)
4757                    if (ADOLC_EXT_POINTER_Y==NULL) fail(ADOLC_EXT_DIFF_NULLPOINTER_ARGUMENT);
4758#endif
4759                }
4760
4761                arg = ADOLC_CURRENT_TAPE_INFOS.lowestXLoc_for;
4762                for (loop=0; loop<n; ++loop) {
4763                    if (edfct->dp_x_changes) {
4764                      IF_KEEP_WRITE_TAYLOR(arg, keep, k, p);
4765                    }
4766                    edfct->dp_x[loop]=dp_T0[arg];
4767#if !defined(_ZOS_)
4768                    ADOLC_EXT_LOOP
4769                        ADOLC_EXT_POINTER_X[loop]ADOLC_EXT_SUBSCRIPT =
4770                        TAYLOR_BUFFER[arg]ADOLC_EXT_SUBSCRIPT;
4771#endif
4772                    ++arg;
4773                }
4774                arg = ADOLC_CURRENT_TAPE_INFOS.lowestYLoc_for;
4775                for (loop=0; loop<m; ++loop) {
4776                    if (edfct->dp_y_priorRequired) {
4777                      IF_KEEP_WRITE_TAYLOR(arg, keep, k, p);
4778                    }
4779                    edfct->dp_y[loop]=dp_T0[arg];
4780#if !defined(_ZOS_)
4781                    ADOLC_EXT_LOOP
4782                        ADOLC_EXT_POINTER_Y[loop]ADOLC_EXT_SUBSCRIPT =
4783                        TAYLOR_BUFFER[arg]ADOLC_EXT_SUBSCRIPT;
4784#endif
4785                    ++arg;
4786                }
4787
4788                ext_retc = edfct->ADOLC_EXT_FCT_COMPLETE;
4789                MINDEC(ret_c, ext_retc);
4790
4791                res = ADOLC_CURRENT_TAPE_INFOS.lowestXLoc_for;
4792                for (loop=0; loop<n; ++loop) {
4793                    dp_T0[res]=edfct->dp_x[loop];
4794#if !defined(_ZOS_)
4795                    ADOLC_EXT_LOOP
4796                        TAYLOR_BUFFER[res]ADOLC_EXT_SUBSCRIPT =
4797                        ADOLC_EXT_POINTER_X[loop]ADOLC_EXT_SUBSCRIPT;
4798#endif
4799                    ++res;
4800                }
4801                res = ADOLC_CURRENT_TAPE_INFOS.lowestYLoc_for;
4802                for (loop=0; loop<m; ++loop) {
4803                    dp_T0[res]=edfct->dp_y[loop];
4804#if !defined(_ZOS_)
4805                    ADOLC_EXT_LOOP
4806                        TAYLOR_BUFFER[res]ADOLC_EXT_SUBSCRIPT =
4807                        ADOLC_EXT_POINTER_Y[loop]ADOLC_EXT_SUBSCRIPT;
4808#endif
4809                    ++res;
4810                }
4811
4812                break;
4813#endif
4814
4815                /*--------------------------------------------------------------------------*/
4816
4817            default:                                                   /* default */
4818                /* Die here, we screwed up */
4819
4820                fprintf(DIAG_OUT,"ADOL-C fatal error in " GENERATED_FILENAME " ("
4821                        __FILE__
4822                        ") : no such operation %d\n", operation);
4823                exit(-1);
4824                break;
4825
4826        } /* endswitch */
4827
4828        /* Read the next operation */
4829        operation=get_op_f();
4830#if defined(ADOLC_DEBUG)
4831        ++countPerOperation[operation];
4832#endif /* ADOLC_DEBUG */
4833    }  /* endwhile */
4834
4835
4836#if defined(ADOLC_DEBUG)
4837    printf("\nTape contains:\n");
4838    for (v = 0; v < 256; ++v)
4839        if (countPerOperation[v] > 0)
4840            printf("operation %3d: %6d time(s) - %6d taylors written (%10.2f per operation)\n", v, countPerOperation[v], taylorPerOperation[v], (double)taylorPerOperation[v] / (double)countPerOperation[v]);
4841    printf("\n");
4842#endif /* ADOLC_DEBUG */
4843
4844#if defined(_KEEP_)
4845    if (keep) taylor_close(taylbuf);
4846#endif
4847
4848    /* clean up */
4849#if !defined (_NTIGHT_)
4850    free(dp_T0);
4851#endif /* !_NTIGHT_ */
4852#if !defined(_INDO_)
4853#if !defined(_ZOS_)
4854#   if defined(_FOS_)
4855    free(dp_T);
4856#   else
4857#if !defined (_INT_FOR_)
4858    myfree2(dpp_T);
4859    free(dp_Ttemp);
4860#endif /* !_NTIGHT_ */
4861#endif
4862#endif
4863#endif
4864#if defined(_HIGHER_ORDER_)
4865    free(dp_z);
4866#endif
4867
4868    end_sweep();
4869
4870
4871#if defined(_INDO_)
4872#if defined(_INDOPRO_)
4873    for(i=0;i<max_ind_dom;i++)
4874      {
4875        free(ind_dom[i]);
4876      }
4877    free(ind_dom);
4878#endif
4879#if defined(_NONLIND_)
4880    for( i=0; i < indcheck; i++) {
4881      traverse_crs(&nonl_dom[i],&sod[i],indcheck+1);
4882      free_tree(&nonl_dom[i],indcheck+1);
4883      crs[i] = (unsigned int*) malloc(sizeof(unsigned int) * (sod[i].entry+1));
4884      crs[i][0] = sod[i].entry;
4885      temp = sod[i].left;
4886      for( ii=1; ii <=sod[i].entry; ii++)
4887        {
4888          crs[i][ii] = temp->entry;
4889          temp1 = temp->left;
4890          free(temp);
4891          temp = temp1;
4892        }
4893    }
4894
4895    free(sod);
4896    free(nonl_dom);
4897    free(fod);
4898    free(arg_index);
4899
4900#endif
4901#if defined(_NONLIND_OLD_)
4902
4903    for( i=0; i < indcheck; i++) {
4904       crs[i] = (unsigned int*) malloc(sizeof(unsigned int) * (nonl_dom[i][0]+1));
4905       crs[i][0] = nonl_dom[i][0];
4906       for(l=1; l < crs[i][0]+1; l++)
4907          crs[i][l] = nonl_dom[i][l+1];
4908       free(nonl_dom[i]);
4909    }
4910    free(nonl_dom);
4911
4912#endif
4913#endif
4914    return ret_c;
4915}
4916
4917/****************************************************************************/
4918
4919#if defined(_INDOPRO_) && !defined(_NONLIND_OLD_)
4920
4921/****************************************************************************/
4922/* set operations for propagation of index domains                          */
4923
4924/*--------------------------------------------------------------------------*/
4925/* operations on index domains                                              */
4926
4927#if defined(_TIGHT_)
4928void copy_index_domain(int res, int arg, locint **ind_dom) {
4929
4930   int i;
4931
4932   if (ind_dom[arg][0] > ind_dom[res][1])
4933     {
4934       free(ind_dom[res]);
4935       ind_dom[res] = (locint *)  malloc(sizeof(locint) * 2*(ind_dom[arg][0]+1));
4936       ind_dom[res][1] = 2*ind_dom[arg][0];
4937     }
4938
4939
4940    for(i=2;i<ind_dom[arg][0]+2;i++)
4941       ind_dom[res][i] = ind_dom[arg][i];
4942    ind_dom[res][0] = ind_dom[arg][0];
4943}
4944
4945
4946void merge_2_index_domains(int res, int arg, locint **ind_dom)
4947{
4948
4949  int num,num1,num2, i,j,k,l;
4950  locint *temp_array, *arg_ind_dom, *res_ind_dom;
4951
4952  if (ind_dom[res][0] == 0)
4953    copy_index_domain(res,arg,ind_dom);
4954  else
4955    {
4956      if (res != arg)
4957     {
4958       arg_ind_dom = ind_dom[arg];
4959       res_ind_dom = ind_dom[res];
4960
4961       num  = ind_dom[res][0];
4962       num1 = arg_ind_dom[0];
4963       num2 = ind_dom[res][1];
4964
4965       if (num2 < num1+num)
4966         num2 = num1+num;
4967
4968       temp_array = (locint *)  malloc(sizeof(locint)* (num2+2));
4969       temp_array[1] = num2;
4970
4971       i = 2;
4972       j = 2;
4973       k = 2;
4974       num += 2;
4975       num1 += 2;
4976       while ((i< num) && (j < num1))
4977         {
4978           if (res_ind_dom[i] < arg_ind_dom[j])
4979          {
4980            temp_array[k] = res_ind_dom[i];
4981            i++; k++;
4982          }
4983           else
4984          {
4985            if (res_ind_dom[i] == arg_ind_dom[j])
4986              {
4987                temp_array[k] = arg_ind_dom[j];
4988                i++;j++;k++;
4989              }
4990            else
4991              {
4992                temp_array[k] = arg_ind_dom[j];
4993                j++;k++;
4994              }
4995          }
4996         }
4997       for(l = i;l<num;l++)
4998         {
4999           temp_array[k] = res_ind_dom[l];
5000           k++;
5001         }
5002       for(l = j;l<num1;l++)
5003         {
5004           temp_array[k] = arg_ind_dom[l];
5005           k++;
5006         }
5007       temp_array[0] = k-2;
5008       free(ind_dom[res]);
5009       ind_dom[res]=temp_array;
5010     }
5011    }
5012
5013
5014}
5015
5016void combine_2_index_domains(int res, int arg1, int arg2, locint **ind_dom) {
5017
5018    if (res != arg1)
5019       copy_index_domain(res, arg1, ind_dom);
5020
5021    merge_2_index_domains(res, arg2, ind_dom);
5022}
5023
5024void merge_3_index_domains(int res, int arg1, int arg2, locint **ind_dom) {
5025    merge_2_index_domains(res, arg1, ind_dom);
5026    merge_2_index_domains(res, arg2, ind_dom);
5027}
5028
5029
5030
5031#endif
5032#endif
5033
5034
5035#if defined(_NONLIND_)
5036#if defined(_TIGHT_)
5037
5038void free_tree(IndexElement* tree, int num)
5039{
5040
5041  if (tree->left != NULL)
5042    {
5043      free_tree(tree->left,num);
5044    }
5045  if (tree->right != NULL)
5046    {
5047      free_tree(tree->right,num);
5048     }
5049    {
5050      if (tree->entry == num)
5051        free(tree);
5052
5053    }
5054 
5055}
5056void traverse_crs(IndexElement* tree,  IndexElement_sod* sod, int num)
5057{
5058
5059  IndexElement_sod *temp, *temp1;
5060  int ii;
5061
5062  if (tree->left != NULL)
5063    {
5064      traverse_crs(tree->left, sod, num);
5065    }
5066  if (tree->right != NULL)
5067    {
5068      traverse_crs(tree->right, sod, num);
5069    }
5070  if (tree->entry < num)
5071    {
5072      temp = sod->left;
5073      if (temp == NULL)
5074        {
5075          temp = (struct IndexElement_sod*) malloc(sizeof(struct IndexElement_sod));
5076          temp->left = NULL;
5077          temp->entry = tree->entry;
5078          sod->entry++;
5079          sod->left=temp;
5080        }
5081      else
5082        {
5083          while ((temp->entry < tree->entry) && (temp->left != NULL))
5084            {
5085              temp1 = temp;
5086              temp = temp->left;
5087            }
5088          if (temp->left == NULL)
5089            {
5090              if(temp->entry < tree->entry)
5091                {
5092                  temp->left = (struct IndexElement_sod*) malloc(sizeof(struct IndexElement_sod));
5093                  temp = temp->left;
5094                  temp->left = NULL;
5095                  temp->entry = tree->entry;
5096                  sod->entry++;
5097                }
5098              if(temp->entry > tree->entry)
5099                {
5100                  temp->left = (struct IndexElement_sod*) malloc(sizeof(struct IndexElement_sod));
5101                  temp->left->entry = temp->entry;
5102                  temp->left->left = NULL;
5103                  temp->entry = tree->entry;
5104                  sod->entry++;
5105                }
5106            }
5107          else
5108            {
5109              if (temp->entry > tree->entry)
5110                {
5111                  temp1 = (struct IndexElement_sod*) malloc(sizeof(struct IndexElement_sod));
5112                  temp1->left = temp->left;
5113                  temp1->entry = temp->entry;
5114                  temp->entry = tree->entry;
5115                  temp->left=temp1;
5116                  sod->entry++;
5117                }
5118             
5119            }
5120        }
5121    }
5122}
5123
5124void traverse_unary(IndexElement* tree,  IndexElement* nonl_dom,  IndexElement* fodi, int num, int maxopind)
5125{
5126  IndexElement *temp;
5127
5128  if (tree->left != NULL)
5129    {
5130      traverse_unary(tree->left, nonl_dom, fodi, num, maxopind);
5131      if (tree->right != NULL)
5132        {
5133          traverse_unary(tree->right, nonl_dom, fodi, num, maxopind);
5134        }
5135     }
5136  else
5137    {
5138      if(tree->entry<maxopind)
5139        {
5140          temp = (struct IndexElement*) malloc(sizeof(struct IndexElement));
5141          temp->right = fodi;
5142          temp->left = nonl_dom[tree->entry].left;
5143          temp->entry= num;
5144          nonl_dom[tree->entry].left = temp;
5145        }
5146    }
5147}
5148
5149#endif
5150#endif
5151
5152#if defined(_NONLIND_OLD_)
5153#if defined(_TIGHT_)
5154
5155void extend_nonlinearity_domain_binary_step
5156(int arg1, int arg2, locint **ind_dom, locint **nonl_dom) 
5157{
5158  int index,num,num1, num2, i,j,k,l,m;
5159  locint *temp_nonl, *index_nonl_dom, *arg1_ind_dom, *arg2_ind_dom;
5160
5161  num = ind_dom[arg2][0];
5162
5163  for(m=2;m<ind_dom[arg1][0]+2;m++) 
5164    {
5165      index = ind_dom[arg1][m];
5166      index_nonl_dom = nonl_dom[index];
5167
5168      if (index_nonl_dom[0] == 0)  /* empty list */
5169        {
5170          if ( index_nonl_dom[1] < num)
5171            {
5172              free(index_nonl_dom);
5173              index_nonl_dom = (locint*) malloc(sizeof(locint)*2*(num+1) );
5174              index_nonl_dom[1] = 2*num;
5175            }
5176          for(i=2;i<num+2;i++)      /* append index domain list of "arg" */
5177            index_nonl_dom[i] = ind_dom[arg2][i];
5178          index_nonl_dom[0] = num;
5179        } 
5180      else 
5181        { /* merge lists */
5182          num1 = index_nonl_dom[0];
5183          num2 = index_nonl_dom[1];
5184         
5185          if (num1+num > num2)
5186            num2 = num1+num;
5187         
5188          temp_nonl = (locint*) malloc(sizeof(locint)*(num2+2));
5189          temp_nonl[1] = num2;
5190         
5191          i = 2;
5192          k = 2;
5193          j = 2;
5194          num1 +=2;
5195          num2 = num+2;
5196          while ((i<num1) && (j < num2)){
5197            if (ind_dom[arg2][j] < index_nonl_dom[i]) /* < */ {
5198              temp_nonl[k] = ind_dom[arg2][j];
5199              j++; k++;
5200            } else {
5201              if (ind_dom[arg2][j] == index_nonl_dom[i])  /* == */ {
5202                temp_nonl[k] = ind_dom[arg2][j];
5203                j++; k++; i++;
5204              } else {
5205                temp_nonl[k] = index_nonl_dom[i];
5206                i++; k++;
5207              }
5208            }
5209          }
5210          for(l = j;l<num2;l++) {
5211            temp_nonl[k] = ind_dom[arg2][l];
5212            k++;
5213          }
5214          for(l = i;l<num1;l++) {
5215            temp_nonl[k] = index_nonl_dom[l];
5216            k++;
5217          }
5218          temp_nonl[0] = k-2; 
5219          free((char*) nonl_dom[index]);
5220          nonl_dom[index] = temp_nonl;
5221        }
5222    }
5223}
5224
5225void extend_nonlinearity_domain_unary
5226(int arg, locint **ind_dom, locint **nonl_dom) {
5227    extend_nonlinearity_domain_binary_step(arg, arg, ind_dom, nonl_dom);
5228}
5229
5230void extend_nonlinearity_domain_binary
5231(int arg1, int arg2, locint **ind_dom, locint **nonl_dom) {
5232    extend_nonlinearity_domain_binary_step(arg1, arg2, ind_dom, nonl_dom);
5233    extend_nonlinearity_domain_binary_step(arg2, arg1, ind_dom, nonl_dom);
5234}
5235
5236
5237#endif
5238#endif
5239END_C_DECLS
Note: See TracBrowser for help on using the repository browser.