source: stable/2.4/ADOL-C/src/uni5_for.c @ 412

Last change on this file since 412 was 376, checked in by kulshres, 7 years ago

use the option --enable-advanced-branching for new comparison operations

the configure option --enable-advanced-branching sets ADOLC_ADVANCED_BRANCHING
and this is unset by default. The new comparison operators that can
handle branching in conjunction with condassign and advectors will only
be compiled if this option is used.

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 376 2012-12-14 11:22:56Z 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#if defined(ADOLC_ADVANCED_BRANCHING)
3987            case neq_a_a:
3988            case eq_a_a:
3989            case le_a_a:
3990            case ge_a_a:
3991            case lt_a_a:
3992            case gt_a_a:
3993                coval = get_val_f();
3994                arg = get_locint_f();
3995                arg1 = get_locint_f();
3996                res = get_locint_f();
3997#if !defined(_NTIGHT_)
3998                {
3999                    revreal retval = -1;
4000                    const char* opname = "";
4001                    switch (operation) {
4002                    case neq_a_a:
4003                        retval = (revreal)(dp_T0[arg] != dp_T0[arg1]);
4004                        opname = "neq_a_a";
4005                        break;
4006                    case eq_a_a:
4007                        retval = (revreal)(dp_T0[arg] == dp_T0[arg1]);
4008                        opname = "eq_a_a";
4009                        break;
4010                    case ge_a_a:
4011                        retval = (revreal)(dp_T0[arg] >= dp_T0[arg1]);
4012                        opname = "ge_a_a";
4013                        break;
4014                    case le_a_a:
4015                        retval = (revreal)(dp_T0[arg] <= dp_T0[arg1]);
4016                        opname = "le_a_a";
4017                        break;
4018                    case gt_a_a:
4019                        retval = (revreal)(dp_T0[arg] > dp_T0[arg1]);
4020                        opname = "gt_a_a";
4021                        break;
4022                    case lt_a_a:
4023                        retval = (revreal)(dp_T0[arg] < dp_T0[arg1]);
4024                        opname = "lt_a_a";
4025                        break;
4026                    }
4027                    if (retval != coval && ADOLC_GLOBAL_TAPE_VARS.branchSwitchWarning)
4028                        fprintf(DIAG_OUT,
4029                                "ADOL-C Warning: Branch switch detected in comparison "
4030                                "(operator %s).\n"
4031                                "Results may be unpredictable! Retaping recommended!\n",opname);
4032                    dp_T0[res] = retval;
4033                }
4034#endif
4035#if defined(_INDO_)
4036#if defined(_INDOPRO_)
4037                ind_dom[res][0]=0;
4038#endif
4039#if defined(_NONLIND_)
4040                fod[opind].entry = maxopind+2;
4041                fod[opind].left = NULL;
4042                fod[opind].right = NULL;
4043                arg_index[res] = opind++;
4044#endif
4045#else
4046#if !defined(_ZOS_) /* BREAK_ZOS */
4047                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
4048
4049                FOR_0_LE_l_LT_pk
4050                TRES_INC = 0;
4051#endif
4052#endif /* ALL_TOGETHER_AGAIN */
4053
4054                break;
4055#endif /* ADVANCED_BRANCHING */
4056
4057                /*--------------------------------------------------------------------------*/
4058            case subscript:
4059                coval = get_val_f();
4060                arg = get_locint_f();
4061                res = get_locint_f();
4062                {
4063                    size_t cnt, idx, numvar = (size_t)trunc(fabs(coval));
4064                    locint vectorloc;
4065                    vectorloc = get_locint_f();
4066#if !defined(_NTIGHT_)
4067                    idx = (size_t)trunc(fabs(dp_T0[arg]));
4068                    if (idx >= numvar)
4069                        fprintf(DIAG_OUT, "ADOL-C warning: index out of bounds while subscripting n=%z, idx=%z\n", numvar, idx);
4070                    arg1 = vectorloc+idx;
4071                    IF_KEEP_WRITE_TAYLOR(res,keep,k,p);
4072                    dp_T0[res] = dp_T0[arg1];
4073#if defined(_INDO_)
4074#if defined(_INDOPRO_)
4075                    copy_index_domain(res, arg1, ind_dom);
4076#endif
4077#if defined(_NONLIND_)
4078                    arg_index[res] = arg_index[arg1];
4079#endif
4080#else
4081#if !defined(_ZOS_) /* BREAK_ZOS */
4082                    ASSIGN_T(Targ1,TAYLOR_BUFFER[arg1])
4083                    ASSIGN_T(Tres,TAYLOR_BUFFER[res])
4084
4085                    FOR_0_LE_l_LT_pk
4086                    TRES_INC = TARG1_INC;
4087#endif
4088#endif
4089#else
4090                    fprintf(DIAG_OUT, "ADOL-C error: active subscripting does not work in safe mode, please use tight mode\n");
4091                    exit(-2);
4092#endif /* ALL_TOGETHER_AGAIN */
4093                }
4094                break;
4095
4096            case subscript_ref:
4097                coval = get_val_f();
4098                arg = get_locint_f();
4099                res = get_locint_f();
4100                {
4101                    size_t cnt, idx, numvar = (size_t)trunc(fabs(coval));
4102                    locint vectorloc;
4103                    vectorloc = get_locint_f();
4104#if !defined(_NTIGHT_)
4105                    idx = (size_t)trunc(fabs(dp_T0[arg]));
4106                    if (idx >= numvar)
4107                        fprintf(DIAG_OUT, "ADOL-C warning: index out of bounds while subscripting (ref) n=%z, idx=%z\n", numvar, idx);
4108                    arg1 = vectorloc+idx;
4109                    IF_KEEP_WRITE_TAYLOR(res,keep,k,p);
4110                    dp_T0[res] = arg1;
4111#else
4112                    fprintf(DIAG_OUT, "ADOL-C error: active subscripting does not work in safe mode, please use tight mode\n");
4113                    exit(-2);
4114#endif
4115                }
4116                break;
4117
4118            case ref_copyout:
4119                arg = get_locint_f();
4120                res = get_locint_f();
4121#if !defined(_NTIGHT_)
4122                arg1 = (size_t)trunc(fabs(dp_T0[arg]));
4123                IF_KEEP_WRITE_TAYLOR(res,keep,k,p);
4124                dp_T0[res] = dp_T0[arg1];
4125#if defined(_INDO_)
4126#if defined(_INDOPRO_)
4127                copy_index_domain(res, arg1, ind_dom);
4128#endif
4129#if defined(_NONLIND_)
4130                arg_index[res] = arg_index[arg1];
4131#endif
4132#else
4133#if !defined(_ZOS_) /* BREAK_ZOS */
4134                ASSIGN_T(Targ1,TAYLOR_BUFFER[arg1])
4135                ASSIGN_T(Tres,TAYLOR_BUFFER[res])
4136
4137                FOR_0_LE_l_LT_pk
4138                TRES_INC = TARG1_INC;
4139#endif
4140#endif
4141#else
4142                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4143                exit(-2);
4144#endif /* ALL_TOGETHER_AGAIN */
4145                break;
4146
4147            case ref_incr_a:
4148                arg = get_locint_f();
4149#if !defined(_NTIGHT_)
4150                arg1 = (size_t)trunc(fabs(dp_T0[arg]));
4151                IF_KEEP_WRITE_TAYLOR(arg1,keep,k,p);
4152                dp_T0[arg1]++;
4153#else
4154                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4155                exit(-2);
4156#endif
4157                break;
4158
4159            case ref_decr_a:
4160                arg = get_locint_f();
4161#if !defined(_NTIGHT_)
4162                arg1 = (size_t)trunc(fabs(dp_T0[arg]));
4163                IF_KEEP_WRITE_TAYLOR(arg1,keep,k,p);
4164                dp_T0[arg1]--;
4165#else
4166                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4167                exit(-2);
4168#endif
4169                break;
4170
4171            case ref_assign_d:
4172                arg = get_locint_f();
4173                coval = get_val_f();
4174               
4175#if !defined(_NTIGHT_)
4176                arg1 = (size_t)trunc(fabs(dp_T0[arg]));
4177                IF_KEEP_WRITE_TAYLOR(arg1,keep,k,p)
4178                dp_T0[arg1] = coval;
4179#if defined(_INDO_)
4180#if defined(_INDOPRO_)
4181                ind_dom[arg1][0] = 0;
4182#endif
4183#if defined(_NONLIND_)
4184                fod[opind].entry = maxopind+2;
4185                fod[opind].left = NULL;
4186                fod[opind].right = NULL;
4187                arg_index[arg1] = opind++;
4188#endif
4189#else
4190#if !defined(_ZOS_)
4191                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
4192
4193                FOR_0_LE_l_LT_pk
4194                TARG1_INC = 0;
4195#endif
4196#endif
4197#else
4198                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4199                exit(-2);
4200#endif
4201                break;
4202
4203            case ref_assign_d_zero:
4204                arg = get_locint_f();
4205
4206#if !defined(_NTIGHT_)
4207                arg1 = (size_t)trunc(fabs(dp_T0[arg]));
4208                IF_KEEP_WRITE_TAYLOR(arg1,keep,k,p)
4209                dp_T0[arg1] = 0.0;
4210#if defined(_INDO_)
4211#if defined(_INDOPRO_)
4212                ind_dom[arg1][0] = 0;
4213#endif
4214#if defined(_NONLIND_)
4215                fod[opind].entry = maxopind+2;
4216                fod[opind].left = NULL;
4217                fod[opind].right = NULL;
4218                arg_index[arg1] = opind++;
4219#endif
4220#else
4221#if !defined(_ZOS_)
4222                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
4223
4224                FOR_0_LE_l_LT_pk
4225                TARG1_INC = 0;
4226#endif
4227#endif
4228#else
4229                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4230                exit(-2);
4231#endif
4232                break;
4233
4234            case ref_assign_d_one:
4235                arg = get_locint_f();
4236
4237#if !defined(_NTIGHT_)
4238                arg1 = (size_t)trunc(fabs(dp_T0[arg]));
4239                IF_KEEP_WRITE_TAYLOR(arg1,keep,k,p)
4240                dp_T0[arg1] = 1.0;
4241#if defined(_INDO_)
4242#if defined(_INDOPRO_)
4243                ind_dom[arg1][0] = 0;
4244#endif
4245#if defined(_NONLIND_)
4246                fod[opind].entry = maxopind+2;
4247                fod[opind].left = NULL;
4248                fod[opind].right = NULL;
4249                arg_index[arg1] = opind++;
4250#endif
4251#else
4252#if !defined(_ZOS_)
4253                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
4254
4255                FOR_0_LE_l_LT_pk
4256                TARG1_INC = 0;
4257#endif
4258#endif
4259#else
4260                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4261                exit(-2);
4262#endif
4263                break;
4264
4265            case ref_assign_a:           /* assign an adouble variable an    assign_a */
4266                /* adouble value. (=) */
4267                arg = get_locint_f();
4268                res = get_locint_f();
4269
4270#if !defined(_NTIGHT_)
4271                arg1 = (size_t)trunc(fabs(dp_T0[res]));
4272                IF_KEEP_WRITE_TAYLOR(arg1,keep,k,p)
4273                dp_T0[arg1] = dp_T0[arg];
4274#if defined(_INDO_)
4275#if defined(_INDOPRO_)
4276                copy_index_domain(arg1, arg, ind_dom);
4277#endif
4278#if defined(_NONLIND_)
4279                arg_index[arg1] = arg_index[arg];
4280#endif
4281#else
4282#if !defined(_ZOS_) /* BREAK_ZOS */
4283                ASSIGN_T(Targ,TAYLOR_BUFFER[arg])
4284                ASSIGN_T(Targ1,TAYLOR_BUFFER[arg1])
4285
4286                FOR_0_LE_l_LT_pk
4287                TARG1_INC = TARG_INC;
4288#endif
4289#endif
4290#else
4291                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4292                exit(-2);
4293#endif /* ALL_TOGETHER_AGAIN */
4294                break;
4295
4296            case ref_assign_ind:       /* assign an adouble variable an    assign_ind */
4297                /* independent double value (<<=) */
4298                arg = get_locint_f();
4299
4300
4301#if !defined(_NTIGHT_)
4302                res = (size_t)trunc(fabs(dp_T0[arg]));
4303                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4304                dp_T0[res] = basepoint[indexi];
4305#if defined(_INDO_)
4306#if defined(_INDOPRO_)
4307                ind_dom[res][0] = 1;
4308                ind_dom[res][2] = indexi;
4309#endif
4310#if defined(_NONLIND_)
4311                fod[opind].entry = indexi;
4312                fod[opind].left = NULL;
4313                fod[opind].right = NULL;
4314                arg_index[res] = opind++;
4315#endif
4316#else
4317#if !defined(_ZOS_) /* BREAK_ZOS */
4318                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
4319
4320#ifdef _INT_FOR_
4321                FOR_0_LE_l_LT_p
4322                TRES_INC = ARGUMENT(indexi,l,i);
4323#else
4324                FOR_0_LE_l_LT_p
4325                FOR_0_LE_i_LT_k
4326                TRES_INC = ARGUMENT(indexi,l,i);
4327#endif
4328#endif
4329#endif
4330#else
4331                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4332                exit(-2);
4333#endif /* ALL_TOGETHER_AGAIN */
4334                ++indexi;
4335                break;
4336
4337            case ref_eq_plus_d:            /* Add a floating point to an    eq_plus_d */
4338                /* adouble. (+=) */
4339                arg  = get_locint_f();
4340                coval = get_val_f();
4341
4342
4343#if !defined(_NTIGHT_)
4344                res = (size_t)trunc(fabs(dp_T0[arg]));
4345                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4346                dp_T0[res] += coval;
4347#else
4348                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4349                exit(-2);
4350#endif /* !_NTIGHT_ */
4351                break;
4352
4353                /*--------------------------------------------------------------------------*/
4354            case ref_eq_plus_a:             /* Add an adouble to another    eq_plus_a */
4355                /* adouble. (+=) */
4356                arg = get_locint_f();
4357                arg1 = get_locint_f();
4358
4359#if !defined(_NTIGHT_)
4360                res = (size_t)trunc(fabs(dp_T0[arg1]));
4361                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4362                dp_T0[res] += dp_T0[arg];
4363#if defined(_INDO_)
4364#if defined(_INDOPRO_)
4365                merge_2_index_domains(res, arg, ind_dom);
4366#endif
4367#if defined(_NONLIND_)
4368                fod[opind].entry = maxopind+2;
4369                fod[opind].left = &fod[arg_index[res]];
4370                fod[opind].right = &fod[arg_index[arg]];
4371                arg_index[res] = opind++;
4372#endif
4373#else
4374#if !defined(_ZOS_) /* BREAK_ZOS */
4375                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
4376                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
4377
4378#ifdef _INT_FOR_
4379                FOR_0_LE_l_LT_pk
4380                TRES_INC |= TARG_INC;
4381#else
4382                FOR_0_LE_l_LT_pk
4383                TRES_INC += TARG_INC;
4384#endif
4385#endif
4386#endif
4387#else
4388                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4389                exit(-2);
4390#endif /* ALL_TOGETHER_AGAIN */
4391                break;
4392
4393            case ref_eq_min_d:       /* Subtract a floating point from an    eq_min_d */
4394                /* adouble. (-=) */
4395                arg = get_locint_f();
4396                coval = get_val_f();
4397
4398#if !defined(_NTIGHT_)
4399                res = (size_t)trunc(fabs(dp_T0[arg]));
4400                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4401                dp_T0[res] -= coval;
4402#else
4403                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4404                exit(-2);
4405#endif /* !_NTIGHT_ */
4406                break;
4407
4408                /*--------------------------------------------------------------------------*/
4409            case ref_eq_min_a:        /* Subtract an adouble from another    eq_min_a */
4410                /* adouble. (-=) */
4411                arg = get_locint_f();
4412                arg1 = get_locint_f();
4413
4414#if !defined(_NTIGHT_)
4415                res = (size_t)trunc(fabs(dp_T0[arg1]));
4416                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4417                dp_T0[res] -= dp_T0[arg];
4418#if defined(_INDO_)
4419#if defined(_INDOPRO_)
4420                merge_2_index_domains(res, arg, ind_dom);
4421#endif
4422#if defined(_NONLIND_)
4423                fod[opind].entry = maxopind+2;
4424                fod[opind].left = &fod[arg_index[res]];
4425                fod[opind].right = &fod[arg_index[arg]];
4426                arg_index[res] = opind++;
4427#endif
4428#else
4429#if !defined(_ZOS_) /* BREAK_ZOS */
4430                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
4431                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
4432
4433#ifdef _INT_FOR_
4434                FOR_0_LE_l_LT_pk
4435                TRES_INC |= TARG_INC;
4436#else
4437                FOR_0_LE_l_LT_pk
4438                TRES_INC -= TARG_INC;
4439#endif
4440#endif
4441#endif
4442#else
4443                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4444                exit(-2);
4445#endif /* ALL_TOGETHER_AGAIN */
4446                break;
4447
4448            case ref_eq_mult_d:              /* Multiply an adouble by a    eq_mult_d */
4449                /* flaoting point. (*=) */
4450                arg = get_locint_f();
4451                coval = get_val_f();
4452
4453#if !defined(_NTIGHT_)
4454                res = (size_t)trunc(fabs(dp_T0[arg]));
4455                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4456                dp_T0[res] *= coval;
4457#if !defined(_INDO_)
4458#if !defined(_ZOS_) /* BREAK_ZOS */
4459#if !defined( _INT_FOR_)
4460
4461                FOR_0_LE_l_LT_pk
4462                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
4463
4464                FOR_0_LE_l_LT_pk
4465                TRES_INC *= coval;
4466#endif
4467#endif
4468#endif
4469#else
4470                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4471                exit(-2);
4472#endif /* ALL_TOGETHER_AGAIN */
4473                break;
4474
4475            case ref_eq_mult_a:       /* Multiply one adouble by another    eq_mult_a */
4476                /* (*=) */
4477                arg = get_locint_f();
4478                arg1 = get_locint_f();
4479
4480#if !defined(_NTIGHT_)
4481                res = (size_t)trunc(fabs(dp_T0[arg1]));
4482                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4483#if defined(_INDO_)
4484#if defined(_INDOPRO_)
4485                merge_2_index_domains(res, arg, ind_dom);
4486#endif
4487#if defined(_NONLIND_)
4488                fod[opind].entry = maxopind+2;
4489                fod[opind].left = &fod[arg_index[res]];
4490                fod[opind].right = &fod[arg_index[arg]];
4491                traverse_unary(&fod[arg_index[res]], nonl_dom, &fod[arg_index[arg]], indcheck+1,maxopind+2);
4492                traverse_unary(&fod[arg_index[arg]], nonl_dom, &fod[arg_index[res]], indcheck+1,maxopind+2);
4493                arg_index[res] = opind++;
4494#endif
4495#if defined(_NONLIND_OLD_)
4496                extend_nonlinearity_domain_binary(res, arg, ind_dom, nonl_dom);
4497#endif
4498#else
4499#if !defined(_ZOS_) /* BREAK_ZOS */
4500                ASSIGN_T(Tres, TAYLOR_BUFFER[res])
4501                ASSIGN_T(Targ, TAYLOR_BUFFER[arg])
4502
4503                INC_pk_1(Tres)
4504                INC_pk_1(Targ)
4505
4506#ifdef _INT_FOR_
4507                FOR_p_GT_l_GE_0
4508                TRES_FODEC |= TARG_DEC;
4509#else
4510                FOR_p_GT_l_GE_0
4511                FOR_k_GT_i_GE_0
4512                { TRES_FODEC = dp_T0[res]*TARG_DEC +
4513                               TRES*dp_T0[arg];
4514                  DEC_TRES_FO
4515#ifdef _HIGHER_ORDER_
4516                  TresOP = Tres-i;
4517                  TargOP = Targ;
4518
4519                  for (j=0;j<i;j++)
4520                  *Tres += (*TresOP++) * (*TargOP--);
4521                  Tres--;
4522#endif /* _HIGHER_ORDER_ */
4523                }
4524#endif
4525#endif
4526#endif
4527                dp_T0[res] *= dp_T0[arg];
4528#else
4529                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4530                exit(-2);
4531#endif /* !_NTIGHT_ */
4532                break;
4533
4534            case ref_cond_assign:                                      /* cond_assign */
4535                arg   = get_locint_f();
4536                arg1  = get_locint_f();
4537                arg2  = get_locint_f();
4538                { 
4539                    locint ref = get_locint_f();
4540                    coval = get_val_f();
4541#if !defined(_NTIGHT_)
4542                    res   = (size_t)trunc(fabs(dp_T0[ref]));
4543
4544                    IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4545
4546                /* olvo 980924 changed order to allow reflexive ops */
4547#if defined(_INDO_)
4548                if (dp_T0[arg] > 0) {
4549                    if (coval <= 0.0)
4550                        MINDEC(ret_c,2);
4551                    dp_T0[res] = dp_T0[arg1];
4552
4553#if defined(_INDOPRO_)
4554                    copy_index_domain(res, arg1, ind_dom);
4555#endif
4556#if defined(_NONLIND_)
4557                    arg_index[res] = arg_index[arg1];
4558#endif
4559                } else {
4560                    if (coval > 0.0)
4561                        MINDEC(ret_c,2);
4562                    if (dp_T0[arg] == 0)
4563                        MINDEC(ret_c,0);
4564                    dp_T0[res] = dp_T0[arg2];
4565
4566#if defined(_INDOPRO_)
4567                    copy_index_domain(res, arg2, ind_dom);
4568#endif
4569#if defined(_NONLIND_)
4570                    arg_index[res] = arg_index[arg2];
4571#endif
4572                }
4573#else
4574#if !defined(_ZOS_) /* BREAK_ZOS */
4575                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
4576                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
4577                ASSIGN_T(Targ2, TAYLOR_BUFFER[arg2])
4578#endif /* ALL_TOGETHER_AGAIN */
4579
4580#ifdef _INT_FOR_
4581                coval = get_val_f();
4582
4583                if (dp_T0[arg] > 0)
4584                    FOR_0_LE_l_LT_pk
4585                    TRES_INC = TARG1_INC;
4586                else
4587                    FOR_0_LE_l_LT_pk
4588                    TRES_INC = TARG2_INC;
4589
4590                if (dp_T0[arg] > 0) {
4591                    if (coval <= 0.0)
4592                        MINDEC(ret_c,2);
4593                    dp_T0[res] = dp_T0[arg1];
4594                } else {
4595                    if (coval > 0.0)
4596                        MINDEC(ret_c,2);
4597                    if (dp_T0[arg] == 0)
4598                        MINDEC(ret_c,0);
4599                    dp_T0[res] = dp_T0[arg2];
4600                }
4601                FOR_0_LE_l_LT_pk
4602                TRES_INC = TARG1_INC | TARG2_INC;
4603#else
4604#if !defined(_ZOS_) /* BREAK_ZOS */
4605                if (dp_T0[arg] > 0)
4606                    FOR_0_LE_l_LT_pk
4607                    TRES_INC = TARG1_INC;
4608                else
4609                    FOR_0_LE_l_LT_pk
4610                    TRES_INC = TARG2_INC;
4611#endif
4612
4613                if (dp_T0[arg] > 0) {
4614                    if (coval <= 0.0)
4615                        MINDEC(ret_c,2);
4616                    dp_T0[res] = dp_T0[arg1];
4617                } else {
4618                    if (coval > 0.0)
4619                        MINDEC(ret_c,2);
4620                    if (dp_T0[arg] == 0)
4621                        MINDEC(ret_c,0);
4622                    dp_T0[res] = dp_T0[arg2];
4623                }
4624#endif
4625#endif
4626#else
4627                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4628                exit(-2);
4629#endif /* ALL_TOGETHER_AGAIN */
4630                }
4631                break;
4632
4633            case ref_cond_assign_s:                                  /* cond_assign_s */
4634                arg   = get_locint_f();
4635                arg1  = get_locint_f();
4636                arg2   = get_locint_f();
4637                coval = get_val_f();
4638
4639#if !defined(_NTIGHT_)
4640                res = (size_t)trunc(fabs(dp_T0[arg2]));
4641                IF_KEEP_WRITE_TAYLOR(res,keep,k,p)
4642
4643                /* olvo 980924 changed order to allow reflexive ops */
4644#if defined(_INDO_)
4645                if (dp_T0[arg] > 0) {
4646#if defined(_INDOPRO_)
4647                    copy_index_domain(res, arg1, ind_dom);
4648#endif
4649#if defined(_NONLIND_)
4650                    arg_index[res] = arg_index[arg1];
4651#endif
4652                }
4653#else
4654#if !defined(_ZOS_) /* BREAK_ZOS */
4655                ASSIGN_T(Tres,  TAYLOR_BUFFER[res])
4656                ASSIGN_T(Targ1, TAYLOR_BUFFER[arg1])
4657#endif /* ALL_TOGETHER_AGAIN */
4658
4659#ifdef _INT_FOR_
4660                coval = get_val_f();
4661
4662                if (dp_T0[arg] > 0)
4663                    FOR_0_LE_l_LT_pk
4664                    TRES_INC = TARG1_INC;
4665
4666                if (dp_T0[arg] > 0) {
4667                    if (coval <= 0.0)
4668                        MINDEC(ret_c,2);
4669                    dp_T0[res] = dp_T0[arg1];
4670                } else
4671                    if (dp_T0[arg] == 0)
4672                        MINDEC(ret_c,0);
4673#else
4674#if !defined(_ZOS_) /* BREAK_ZOS */
4675                if (dp_T0[arg] > 0)
4676                    FOR_0_LE_l_LT_pk
4677                    TRES_INC = TARG1_INC;
4678#endif
4679                if (dp_T0[arg] > 0) {
4680                    if (coval <= 0.0)
4681                        MINDEC(ret_c,2);
4682                    dp_T0[res] = dp_T0[arg1];
4683                } else
4684                    if (dp_T0[arg] == 0)
4685                        MINDEC(ret_c,0);
4686#endif
4687#endif
4688#else
4689                fprintf(DIAG_OUT, "ADOL-C error: active vector element referencing does not work in safe mode, please use tight mode\n");
4690                exit(-2);
4691#endif /* ALL_TOGETHER_AGAIN */
4692                break;
4693
4694                /****************************************************************************/
4695                /*                                                          REMAINING STUFF */
4696
4697                /*--------------------------------------------------------------------------*/
4698            case take_stock_op:                                  /* take_stock_op */
4699                size = get_locint_f();
4700                res  = get_locint_f();
4701                d    = get_val_v_f(size);
4702
4703                for (ls=0;ls<size;ls++) {
4704#if !defined(_NTIGHT_)
4705                    dp_T0[res]=*d;
4706#endif /* !_NTIGHT_ */
4707#if !defined(_INDO_)
4708#if !defined(_ZOS_) /* BREAK_ZOS */
4709                    ASSIGN_T(Tres,TAYLOR_BUFFER[res])
4710
4711                    FOR_0_LE_l_LT_pk
4712                    TRES_INC = 0;
4713
4714#endif /* ALL_TOGETHER_AGAIN */
4715                    res++;
4716#if !defined(_NTIGHT_)
4717                    d++;
4718#endif /* !_NTIGHT_ */
4719#endif
4720                }
4721                break;
4722
4723                /*--------------------------------------------------------------------------*/
4724            case death_not:                                          /* death_not */
4725                arg1=get_locint_f();
4726                arg2=get_locint_f();
4727
4728#ifdef _KEEP_
4729                if (keep) {
4730                    do {
4731                        IF_KEEP_WRITE_TAYLOR(arg2,keep,k,p)
4732                    } while(arg1 < arg2-- );
4733                }
4734#endif
4735                break;
4736
4737                /*--------------------------------------------------------------------------*/
4738#if defined(_EXTERN_) /* ZOS,  FOS, FOV up to now */
4739            case ext_diff:                       /* extern differntiated function */
4740                ADOLC_CURRENT_TAPE_INFOS.ext_diff_fct_index=get_locint_f();
4741                n=get_locint_f();
4742                m=get_locint_f();
4743                ADOLC_CURRENT_TAPE_INFOS.lowestXLoc_for = get_locint_f();
4744                ADOLC_CURRENT_TAPE_INFOS.lowestYLoc_for = get_locint_f();
4745                ADOLC_CURRENT_TAPE_INFOS.cpIndex = get_locint_f();
4746                edfct=get_ext_diff_fct(ADOLC_CURRENT_TAPE_INFOS.ext_diff_fct_index);
4747
4748                if (edfct->ADOLC_EXT_FCT_POINTER==NULL)
4749                    fail(ADOLC_EXT_DIFF_NULLPOINTER_DIFFFUNC);
4750                if (n>0) {
4751                    if (edfct->dp_x==NULL) fail(ADOLC_EXT_DIFF_NULLPOINTER_ARGUMENT);
4752#if !defined(_ZOS_)
4753                    if (ADOLC_EXT_POINTER_X==NULL) fail(ADOLC_EXT_DIFF_NULLPOINTER_ARGUMENT);
4754#endif
4755                }
4756                if (m>0) {
4757                    if (edfct->dp_y==NULL) fail(ADOLC_EXT_DIFF_NULLPOINTER_ARGUMENT);
4758#if !defined(_ZOS_)
4759                    if (ADOLC_EXT_POINTER_Y==NULL) fail(ADOLC_EXT_DIFF_NULLPOINTER_ARGUMENT);
4760#endif
4761                }
4762
4763                arg = ADOLC_CURRENT_TAPE_INFOS.lowestXLoc_for;
4764                for (loop=0; loop<n; ++loop) {
4765                    if (edfct->dp_x_changes) {
4766                      IF_KEEP_WRITE_TAYLOR(arg, keep, k, p);
4767                    }
4768                    edfct->dp_x[loop]=dp_T0[arg];
4769#if !defined(_ZOS_)
4770                    ADOLC_EXT_LOOP
4771                        ADOLC_EXT_POINTER_X[loop]ADOLC_EXT_SUBSCRIPT =
4772                        TAYLOR_BUFFER[arg]ADOLC_EXT_SUBSCRIPT;
4773#endif
4774                    ++arg;
4775                }
4776                arg = ADOLC_CURRENT_TAPE_INFOS.lowestYLoc_for;
4777                for (loop=0; loop<m; ++loop) {
4778                    if (edfct->dp_y_priorRequired) {
4779                      IF_KEEP_WRITE_TAYLOR(arg, keep, k, p);
4780                    }
4781                    edfct->dp_y[loop]=dp_T0[arg];
4782#if !defined(_ZOS_)
4783                    ADOLC_EXT_LOOP
4784                        ADOLC_EXT_POINTER_Y[loop]ADOLC_EXT_SUBSCRIPT =
4785                        TAYLOR_BUFFER[arg]ADOLC_EXT_SUBSCRIPT;
4786#endif
4787                    ++arg;
4788                }
4789
4790                ext_retc = edfct->ADOLC_EXT_FCT_COMPLETE;
4791                MINDEC(ret_c, ext_retc);
4792
4793                res = ADOLC_CURRENT_TAPE_INFOS.lowestXLoc_for;
4794                for (loop=0; loop<n; ++loop) {
4795                    dp_T0[res]=edfct->dp_x[loop];
4796#if !defined(_ZOS_)
4797                    ADOLC_EXT_LOOP
4798                        TAYLOR_BUFFER[res]ADOLC_EXT_SUBSCRIPT =
4799                        ADOLC_EXT_POINTER_X[loop]ADOLC_EXT_SUBSCRIPT;
4800#endif
4801                    ++res;
4802                }
4803                res = ADOLC_CURRENT_TAPE_INFOS.lowestYLoc_for;
4804                for (loop=0; loop<m; ++loop) {
4805                    dp_T0[res]=edfct->dp_y[loop];
4806#if !defined(_ZOS_)
4807                    ADOLC_EXT_LOOP
4808                        TAYLOR_BUFFER[res]ADOLC_EXT_SUBSCRIPT =
4809                        ADOLC_EXT_POINTER_Y[loop]ADOLC_EXT_SUBSCRIPT;
4810#endif
4811                    ++res;
4812                }
4813
4814                break;
4815#endif
4816
4817                /*--------------------------------------------------------------------------*/
4818
4819            default:                                                   /* default */
4820                /* Die here, we screwed up */
4821
4822                fprintf(DIAG_OUT,"ADOL-C fatal error in " GENERATED_FILENAME " ("
4823                        __FILE__
4824                        ") : no such operation %d\n", operation);
4825                exit(-1);
4826                break;
4827
4828        } /* endswitch */
4829
4830        /* Read the next operation */
4831        operation=get_op_f();
4832#if defined(ADOLC_DEBUG)
4833        ++countPerOperation[operation];
4834#endif /* ADOLC_DEBUG */
4835    }  /* endwhile */
4836
4837
4838#if defined(ADOLC_DEBUG)
4839    printf("\nTape contains:\n");
4840    for (v = 0; v < 256; ++v)
4841        if (countPerOperation[v] > 0)
4842            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]);
4843    printf("\n");
4844#endif /* ADOLC_DEBUG */
4845
4846#if defined(_KEEP_)
4847    if (keep) taylor_close(taylbuf);
4848#endif
4849
4850    /* clean up */
4851#if !defined (_NTIGHT_)
4852    free(dp_T0);
4853#endif /* !_NTIGHT_ */
4854#if !defined(_INDO_)
4855#if !defined(_ZOS_)
4856#   if defined(_FOS_)
4857    free(dp_T);
4858#   else
4859#if !defined (_INT_FOR_)
4860    myfree2(dpp_T);
4861    free(dp_Ttemp);
4862#endif /* !_NTIGHT_ */
4863#endif
4864#endif
4865#endif
4866#if defined(_HIGHER_ORDER_)
4867    free(dp_z);
4868#endif
4869
4870    end_sweep();
4871
4872
4873#if defined(_INDO_)
4874#if defined(_INDOPRO_)
4875    for(i=0;i<max_ind_dom;i++)
4876      {
4877        free(ind_dom[i]);
4878      }
4879    free(ind_dom);
4880#endif
4881#if defined(_NONLIND_)
4882    for( i=0; i < indcheck; i++) {
4883      traverse_crs(&nonl_dom[i],&sod[i],indcheck+1);
4884      free_tree(&nonl_dom[i],indcheck+1);
4885      crs[i] = (unsigned int*) malloc(sizeof(unsigned int) * (sod[i].entry+1));
4886      crs[i][0] = sod[i].entry;
4887      temp = sod[i].left;
4888      for( ii=1; ii <=sod[i].entry; ii++)
4889        {
4890          crs[i][ii] = temp->entry;
4891          temp1 = temp->left;
4892          free(temp);
4893          temp = temp1;
4894        }
4895    }
4896
4897    free(sod);
4898    free(nonl_dom);
4899    free(fod);
4900    free(arg_index);
4901
4902#endif
4903#if defined(_NONLIND_OLD_)
4904
4905    for( i=0; i < indcheck; i++) {
4906       crs[i] = (unsigned int*) malloc(sizeof(unsigned int) * (nonl_dom[i][0]+1));
4907       crs[i][0] = nonl_dom[i][0];
4908       for(l=1; l < crs[i][0]+1; l++)
4909          crs[i][l] = nonl_dom[i][l+1];
4910       free(nonl_dom[i]);
4911    }
4912    free(nonl_dom);
4913
4914#endif
4915#endif
4916    return ret_c;
4917}
4918
4919/****************************************************************************/
4920
4921#if defined(_INDOPRO_) && !defined(_NONLIND_OLD_)
4922
4923/****************************************************************************/
4924/* set operations for propagation of index domains                          */
4925
4926/*--------------------------------------------------------------------------*/
4927/* operations on index domains                                              */
4928
4929#if defined(_TIGHT_)
4930void copy_index_domain(int res, int arg, locint **ind_dom) {
4931
4932   int i;
4933
4934   if (ind_dom[arg][0] > ind_dom[res][1])
4935     {
4936       free(ind_dom[res]);
4937       ind_dom[res] = (locint *)  malloc(sizeof(locint) * 2*(ind_dom[arg][0]+1));
4938       ind_dom[res][1] = 2*ind_dom[arg][0];
4939     }
4940
4941
4942    for(i=2;i<ind_dom[arg][0]+2;i++)
4943       ind_dom[res][i] = ind_dom[arg][i];
4944    ind_dom[res][0] = ind_dom[arg][0];
4945}
4946
4947
4948void merge_2_index_domains(int res, int arg, locint **ind_dom)
4949{
4950
4951  int num,num1,num2, i,j,k,l;
4952  locint *temp_array, *arg_ind_dom, *res_ind_dom;
4953
4954  if (ind_dom[res][0] == 0)
4955    copy_index_domain(res,arg,ind_dom);
4956  else
4957    {
4958      if (res != arg)
4959     {
4960       arg_ind_dom = ind_dom[arg];
4961       res_ind_dom = ind_dom[res];
4962
4963       num  = ind_dom[res][0];
4964       num1 = arg_ind_dom[0];
4965       num2 = ind_dom[res][1];
4966
4967       if (num2 < num1+num)
4968         num2 = num1+num;
4969
4970       temp_array = (locint *)  malloc(sizeof(locint)* (num2+2));
4971       temp_array[1] = num2;
4972
4973       i = 2;
4974       j = 2;
4975       k = 2;
4976       num += 2;
4977       num1 += 2;
4978       while ((i< num) && (j < num1))
4979         {
4980           if (res_ind_dom[i] < arg_ind_dom[j])
4981          {
4982            temp_array[k] = res_ind_dom[i];
4983            i++; k++;
4984          }
4985           else
4986          {
4987            if (res_ind_dom[i] == arg_ind_dom[j])
4988              {
4989                temp_array[k] = arg_ind_dom[j];
4990                i++;j++;k++;
4991              }
4992            else
4993              {
4994                temp_array[k] = arg_ind_dom[j];
4995                j++;k++;
4996              }
4997          }
4998         }
4999       for(l = i;l<num;l++)
5000         {
5001           temp_array[k] = res_ind_dom[l];
5002           k++;
5003         }
5004       for(l = j;l<num1;l++)
5005         {
5006           temp_array[k] = arg_ind_dom[l];
5007           k++;
5008         }
5009       temp_array[0] = k-2;
5010       free(ind_dom[res]);
5011       ind_dom[res]=temp_array;
5012     }
5013    }
5014
5015
5016}
5017
5018void combine_2_index_domains(int res, int arg1, int arg2, locint **ind_dom) {
5019
5020    if (res != arg1)
5021       copy_index_domain(res, arg1, ind_dom);
5022
5023    merge_2_index_domains(res, arg2, ind_dom);
5024}
5025
5026void merge_3_index_domains(int res, int arg1, int arg2, locint **ind_dom) {
5027    merge_2_index_domains(res, arg1, ind_dom);
5028    merge_2_index_domains(res, arg2, ind_dom);
5029}
5030
5031
5032
5033#endif
5034#endif
5035
5036
5037#if defined(_NONLIND_)
5038#if defined(_TIGHT_)
5039
5040void free_tree(IndexElement* tree, int num)
5041{
5042
5043  if (tree->left != NULL)
5044    {
5045      free_tree(tree->left,num);
5046    }
5047  if (tree->right != NULL)
5048    {
5049      free_tree(tree->right,num);
5050     }
5051    {
5052      if (tree->entry == num)
5053        free(tree);
5054
5055    }
5056 
5057}
5058void traverse_crs(IndexElement* tree,  IndexElement_sod* sod, int num)
5059{
5060
5061  IndexElement_sod *temp, *temp1;
5062  int ii;
5063
5064  if (tree->left != NULL)
5065    {
5066      traverse_crs(tree->left, sod, num);
5067    }
5068  if (tree->right != NULL)
5069    {
5070      traverse_crs(tree->right, sod, num);
5071    }
5072  if (tree->entry < num)
5073    {
5074      temp = sod->left;
5075      if (temp == NULL)
5076        {
5077          temp = (struct IndexElement_sod*) malloc(sizeof(struct IndexElement_sod));
5078          temp->left = NULL;
5079          temp->entry = tree->entry;
5080          sod->entry++;
5081          sod->left=temp;
5082        }
5083      else
5084        {
5085          while ((temp->entry < tree->entry) && (temp->left != NULL))
5086            {
5087              temp1 = temp;
5088              temp = temp->left;
5089            }
5090          if (temp->left == NULL)
5091            {
5092              if(temp->entry < tree->entry)
5093                {
5094                  temp->left = (struct IndexElement_sod*) malloc(sizeof(struct IndexElement_sod));
5095                  temp = temp->left;
5096                  temp->left = NULL;
5097                  temp->entry = tree->entry;
5098                  sod->entry++;
5099                }
5100              if(temp->entry > tree->entry)
5101                {
5102                  temp->left = (struct IndexElement_sod*) malloc(sizeof(struct IndexElement_sod));
5103                  temp->left->entry = temp->entry;
5104                  temp->left->left = NULL;
5105                  temp->entry = tree->entry;
5106                  sod->entry++;
5107                }
5108            }
5109          else
5110            {
5111              if (temp->entry > tree->entry)
5112                {
5113                  temp1 = (struct IndexElement_sod*) malloc(sizeof(struct IndexElement_sod));
5114                  temp1->left = temp->left;
5115                  temp1->entry = temp->entry;
5116                  temp->entry = tree->entry;
5117                  temp->left=temp1;
5118                  sod->entry++;
5119                }
5120             
5121            }
5122        }
5123    }
5124}
5125
5126void traverse_unary(IndexElement* tree,  IndexElement* nonl_dom,  IndexElement* fodi, int num, int maxopind)
5127{
5128  IndexElement *temp;
5129
5130  if (tree->left != NULL)
5131    {
5132      traverse_unary(tree->left, nonl_dom, fodi, num, maxopind);
5133      if (tree->right != NULL)
5134        {
5135          traverse_unary(tree->right, nonl_dom, fodi, num, maxopind);
5136        }
5137     }
5138  else
5139    {
5140      if(tree->entry<maxopind)
5141        {
5142          temp = (struct IndexElement*) malloc(sizeof(struct IndexElement));
5143          temp->right = fodi;
5144          temp->left = nonl_dom[tree->entry].left;
5145          temp->entry= num;
5146          nonl_dom[tree->entry].left = temp;
5147        }
5148    }
5149}
5150
5151#endif
5152#endif
5153
5154#if defined(_NONLIND_OLD_)
5155#if defined(_TIGHT_)
5156
5157void extend_nonlinearity_domain_binary_step
5158(int arg1, int arg2, locint **ind_dom, locint **nonl_dom) 
5159{
5160  int index,num,num1, num2, i,j,k,l,m;
5161  locint *temp_nonl, *index_nonl_dom, *arg1_ind_dom, *arg2_ind_dom;
5162
5163  num = ind_dom[arg2][0];
5164
5165  for(m=2;m<ind_dom[arg1][0]+2;m++) 
5166    {
5167      index = ind_dom[arg1][m];
5168      index_nonl_dom = nonl_dom[index];
5169
5170      if (index_nonl_dom[0] == 0)  /* empty list */
5171        {
5172          if ( index_nonl_dom[1] < num)
5173            {
5174              free(index_nonl_dom);
5175              index_nonl_dom = (locint*) malloc(sizeof(locint)*2*(num+1) );
5176              index_nonl_dom[1] = 2*num;
5177            }
5178          for(i=2;i<num+2;i++)      /* append index domain list of "arg" */
5179            index_nonl_dom[i] = ind_dom[arg2][i];
5180          index_nonl_dom[0] = num;
5181        } 
5182      else 
5183        { /* merge lists */
5184          num1 = index_nonl_dom[0];
5185          num2 = index_nonl_dom[1];
5186         
5187          if (num1+num > num2)
5188            num2 = num1+num;
5189         
5190          temp_nonl = (locint*) malloc(sizeof(locint)*(num2+2));
5191          temp_nonl[1] = num2;
5192         
5193          i = 2;
5194          k = 2;
5195          j = 2;
5196          num1 +=2;
5197          num2 = num+2;
5198          while ((i<num1) && (j < num2)){
5199            if (ind_dom[arg2][j] < index_nonl_dom[i]) /* < */ {
5200              temp_nonl[k] = ind_dom[arg2][j];
5201              j++; k++;
5202            } else {
5203              if (ind_dom[arg2][j] == index_nonl_dom[i])  /* == */ {
5204                temp_nonl[k] = ind_dom[arg2][j];
5205                j++; k++; i++;
5206              } else {
5207                temp_nonl[k] = index_nonl_dom[i];
5208                i++; k++;
5209              }
5210            }
5211          }
5212          for(l = j;l<num2;l++) {
5213            temp_nonl[k] = ind_dom[arg2][l];
5214            k++;
5215          }
5216          for(l = i;l<num1;l++) {
5217            temp_nonl[k] = index_nonl_dom[l];
5218            k++;
5219          }
5220          temp_nonl[0] = k-2; 
5221          free((char*) nonl_dom[index]);
5222          nonl_dom[index] = temp_nonl;
5223        }
5224    }
5225}
5226
5227void extend_nonlinearity_domain_unary
5228(int arg, locint **ind_dom, locint **nonl_dom) {
5229    extend_nonlinearity_domain_binary_step(arg, arg, ind_dom, nonl_dom);
5230}
5231
5232void extend_nonlinearity_domain_binary
5233(int arg1, int arg2, locint **ind_dom, locint **nonl_dom) {
5234    extend_nonlinearity_domain_binary_step(arg1, arg2, ind_dom, nonl_dom);
5235    extend_nonlinearity_domain_binary_step(arg2, arg1, ind_dom, nonl_dom);
5236}
5237
5238
5239#endif
5240#endif
5241END_C_DECLS
Note: See TracBrowser for help on using the repository browser.