source: trunk/ADOL-C/src/ho_rev.c @ 708

Last change on this file since 708 was 708, checked in by kulshres, 3 years ago

Merge branch 'master' of 'gitclone' into 'svn'

The following changes have been merged:

commit e2291bde44a282a133894b0db350aeb0b92a87db
Author: Mladen Banovic <mladenbanovic2705@…>
Date: Fri Jul 8 10:15:51 2016 +0200

Add methods getNumLiveVar and getNumDir in adtl.h, change counter type in FOR_I_EQ_0_LT_NUMDIR macro to size_t (instead of int). Update chunk size of BOOST pool in adouble_tl.cpp according to adouble::numDir.

commit 2ffb294465b973bfd4bf1f73d84478f8233c0d2f
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Thu Jun 23 12:32:14 2016 +0200

implement missing ref_eq_mult_p und ref_eq_min_p in ho_rev.c

somehow these were left out when parameters were being implemented.

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

commit 8cf0e5c1bd36f1dcf3be72cd67de631b2e1d0ee6
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Thu Jun 23 12:31:04 2016 +0200

make sure the result is the last locint written in trace for each operation

since we're trying to generate ascii traces in the future, we'll need this
convention that the last location is the result, and previous locations
are arguments. This has been the case for almost all operations anyway
except for a few new one's that I wrote without keeping this in mind.

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

commit 9ae0ff220f37463f2ed85cafc8a626c24e472f2f
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Tue Jun 21 14:16:27 2016 +0200

on some compilers newer boost interferes with AC_FUNC_MALLOC test

so do AC_FUNC_MALLOC and AC_FUNC_REALLOC as usual and check for boost
library later.

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

commit b746f620772cc8cce53e8f350adc6281279caf72
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Mon Jun 20 15:32:22 2016 +0200

make Klaus Röbenack's name UTF-8 instead of ISO-8859-1

These are the only places where we're not simple ASCII or UTF-8 already

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

commit 1171aa3961b5eb46a5d2ee64751c02a393a8a6f5
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Fri Jun 17 10:42:39 2016 +0200

correct short_ref document about include file

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

commit 2c6b2aac2ef04431ece2c6ff80e574aa2e58814b
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Fri Jun 17 10:40:34 2016 +0200

correct error message to new semantics

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

commit 506cde73451740bf0a15eff7d4abb158ee719ab0
Author: mflehmig <martin.flehmig@…>
Date: Fri Jun 17 10:14:26 2016 +0200

Fixed include of ADOL-C header.

ADOL-C header was included in old fashion (without adolc directory) for this example.

commit 2a023d3281d3d6d9824bad724a5768e3ee2fff94
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Thu Jun 16 13:50:39 2016 +0200

Try to use boost::pool for allocating advals in traceless vector mode

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

commit 80f1e2019ac1faab96fe06f3e9da47efcc1bcd23
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Mon May 23 15:13:22 2016 +0200

correct a reference in doc and rebuild

commit d7ab5283afe58bacb2e8739d72ede4e17f4c8081
Author: Mladen Banovic <mladenbanovic2705@…>
Date: Fri May 20 16:42:13 2016 +0200

Update section 7 of adolc-manual related to the Traceless forward differentiation.

commit bedb8e36f959c5272e4610fe504acc83208e5e9d
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Tue May 17 16:09:36 2016 +0200

macro name correction

commit 92ff596a0331776901df7f172ca347572e3daafd
Author: Kshitij Kulshreshtha <kshitij@…>
Date: Tue May 17 15:56:17 2016 +0200

Add a warning about using static build of ADOL-C

static build of ADOL-C does not call constructors
for internal global objects, thereby causing
segmentation faults.

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

  • Property svn:keywords set to Author Date Id Revision
File size: 115.5 KB
Line 
1/*----------------------------------------------------------------------------
2 ADOL-C -- Automatic Differentiation by Overloading in C++
3 File:     ho_rev.c
4 Revision: $Id: ho_rev.c 708 2016-07-12 08:18:44Z kulshres $
5 Contents: Contains the routines :
6           hos_reverse (higher-order-scalar reverse mode):
7              define _HOS_
8           hos_ov_reverse (higher-order-scalar reverse mode on vectors):
9              define _HOS_OV_
10           hov_reverse (higher-order-vector reverse mode):
11              define _HOV_
12 
13 Copyright (c) Andrea Walther, Andreas Griewank, Andreas Kowarz,
14               Hristo Mitev, Sebastian Schlenkrich, Jean Utke, Olaf Vogel,
15               Kshitij Kulshreshtha
16
17 This file is part of ADOL-C. This software is provided as open source.
18 Any use, reproduction, or distribution of the software constitutes
19 recipient's acceptance of the terms of the accompanying license file.
20 
21----------------------------------------------------------------------------*/
22
23/*****************************************************************************
24 
25  There are four basic versions of the procedure `reverse', which
26  are optimized for the cases of scalar or vector reverse sweeps
27  with first or higher derivatives, respectively. In the calling
28  sequence this distinction is apparent from the type of the
29  parameters `lagrange' and `results'. The former may be left out
30  and the integer parameters `depen', `indep', `degre', and `nrows'
31  must be set or default according to the following matrix of
32  calling cases.
33 
34           no lagrange         double* lagrange     double** lagrange
35 
36double*   gradient of scalar   weight vector times    infeasible
37results   valued function      Jacobian product       combination
38 
39          ( depen = 1 ,         ( depen > 0 ,         
40            degre = 0 ,           degre = 0 ,              ------
41            nrows = 1 )           nrows = 1 )
42 
43double**  Jacobian of vector   weight vector times     weight matrix
44results   valued function      Taylor-Jacobians        times Jacobian
45           
46          ( 0 < depen           ( depen > 0 ,          ( depen > 0 ,
47              = nrows ,           degre > 0 ,            degre = 0 ,
48            degre = 0 )           nrows = 1 )            nrows > 0 )
49 
50double*** full family of         ------------          weigth matrix x
51results   Taylor-Jacobians       ------------          Taylor Jacobians
52 
53*****************************************************************************/
54
55/****************************************************************************/
56/*                                                                   MACROS */
57#undef _ADOLC_VECTOR_
58#undef _HIGHER_ORDER_
59
60/*--------------------------------------------------------------------------*/
61#ifdef _HOS_
62#define GENERATED_FILENAME "hos_reverse"
63
64#define _HIGHER_ORDER_
65
66#define RESULTS(l,indexi,k) results[indexi][k]
67#define LAGRANGE(l,indexd,k)  lagrange[indexd][k]
68
69#define HOV_INC(T,degree) {}
70#define HOS_OV_INC(T,degree) {}
71
72#define GET_TAYL(loc,depth,p) \
73    { \
74        UPDATE_TAYLORREAD(depth) \
75        get_taylors(loc,depth); \
76    }
77
78/*--------------------------------------------------------------------------*/
79#elif _HOS_OV_
80#define GENERATED_FILENAME "hos_ov_reverse"
81
82#define _HIGHER_ORDER_
83
84#define RESULTS(l,indexi,k) results[l][indexi][k]
85#define LAGRANGE(l,indexd,k)  lagrange[indexd][k]
86
87#define HOV_INC(T,degree) T += degree;
88#define HOS_OV_INC(T,degree) T += degree;
89
90#define GET_TAYL(loc,depth,p) \
91    { \
92        UPDATE_TAYLORREAD(depth * p) \
93        get_taylors_p(loc,depth,p); \
94    }
95
96/*--------------------------------------------------------------------------*/
97#elif _HOV_
98#define GENERATED_FILENAME "hov_reverse"
99
100#define _ADOLC_VECTOR_
101#define _HIGHER_ORDER_
102
103#define RESULTS(l,indexi,k) results[l][indexi][k]
104#define LAGRANGE(l,indexd,k)  lagrange[l][indexd][k]
105
106#define IF_HOV_
107#define ENDIF_HOV_
108
109#define HOV_INC(T,degree) T += degree;
110#define HOS_OV_INC(T,degree)
111
112#define GET_TAYL(loc,depth,p) \
113    { \
114        UPDATE_TAYLORREAD(depth) \
115        get_taylors(loc,depth); \
116    }
117
118#else
119#error Error ! Define [_HOS_ | _HOS_OV_ | _HOV_]
120#endif
121
122/*--------------------------------------------------------------------------*/
123/*                                                     access to variables  */
124
125#ifdef _FOS_                                     /* why?, not in fo_rev.c ? */
126#define ARES       *Ares
127#define AARG       *Aarg
128#define AARG1      *Aarg1
129#define AARG2      *Aarg2
130#define AQO        *Aqo
131
132#define ARES_INC   *Ares
133#define AARG_INC   *Aarg
134#define AARG1_INC  *Aarg1
135#define AARG2_INC  *Aarg2
136#define AQO_INC    *Aqo
137
138#define ARES_INC_O  Ares
139#define AARG_INC_O  Aarg
140#define AARG1_INC_O Aarg1
141#define AARG2_INC_O Aarg2
142#define AQO_INC_O   Aqo
143
144#define ASSIGN_A(a,b)  a = &b;
145#define HOS_OV_ASSIGN_A(Aqo,  rp_Atemp)
146#define FOR_0_LE_l_LT_q l = 0;
147
148#elif _HOS_OV_
149#define ARES       *Ares
150#define AARG       *Aarg
151#define AARG1      *Aarg1
152#define AARG2      *Aarg2
153#define AQO        *Aqo
154
155#define ARES_INC   *Ares++
156#define AARG_INC   *Aarg++
157#define AARG1_INC  *Aarg1++
158#define AARG2_INC  *Aarg2++
159#define AQO_INC    *Aqo++
160
161#define ARES_INC_O  Ares++
162#define AARG_INC_O  Aarg++
163#define AARG1_INC_O Aarg1++
164#define AARG2_INC_O Aarg2++
165#define AQO_INC_O   Aqo++
166
167#define ASSIGN_A(a,b)  a = b;
168#define HOS_OV_ASSIGN_A(a, b) a = b;
169#define FOR_0_LE_l_LT_q for(l=0;l<q;l++)
170
171#else  /* _FOV_, _HOS_, _HOV_ */
172#define ARES       *Ares
173#define AARG       *Aarg
174#define AARG1      *Aarg1
175#define AARG2      *Aarg2
176#define AQO        *Aqo
177
178#define ARES_INC   *Ares++
179#define AARG_INC   *Aarg++
180#define AARG1_INC  *Aarg1++
181#define AARG2_INC  *Aarg2++
182#define AQO_INC    *Aqo++
183
184#define ARES_INC_O  Ares++
185#define AARG_INC_O  Aarg++
186#define AARG1_INC_O Aarg1++
187#define AARG2_INC_O Aarg2++
188#define AQO_INC_O   Aqo++
189
190#define ASSIGN_A(a,b)  a = b;
191#define HOS_OV_ASSIGN_A(Aqo, Atemp)
192#define FOR_0_LE_l_LT_q l = 0;
193#endif
194
195#ifdef _HIGHER_ORDER_
196
197#define TRES      *Tres                  /* why ? not used here */
198#define TARG      *Targ
199#define TARG1     *Targ1
200#define TARG2     *Targ2
201
202#define ASSIGN_T(a,b)  a = b;
203#else
204
205#define TRES       rpp_T[res]
206#define TARG       rpp_T[arg]
207#define TARG1      rpp_T[arg1]
208#define TARG2      rpp_T[arg2]
209
210#define ASSIGN_T(a,b)
211#endif
212
213/*--------------------------------------------------------------------------*/
214/*                                                              loop stuff  */
215#ifdef _ADOLC_VECTOR_
216#define FOR_0_LE_l_LT_p for (l=0; l<p; l++)
217#define FOR_p_GT_l_GE_0 for (l=p-1; l>=0; l--)  /* why ? not used here */
218#elif _HOS_OV_
219#define FOR_0_LE_l_LT_p for (l=0; l<p; l++)
220#define FOR_p_GT_l_GE_0                         /* why ? not used here */
221#else
222#define FOR_0_LE_l_LT_p
223#define FOR_p_GT_l_GE_0                         /* why ? not used here */
224#endif
225
226#ifdef _HIGHER_ORDER_
227#define FOR_0_LE_i_LT_k for (i=0; i<k; i++)
228#define FOR_k_GT_i_GE_0 for (i=k-1; i>=0; i--)
229#else
230#define FOR_0_LE_i_LT_k
231#define FOR_k_GT_i_GE_0
232#endif
233
234#ifdef _HOV_
235#define FOR_0_LE_l_LT_pk1 for (l=0; l<pk1; l++)
236#define FOR_0_LE_l_LT_pk for (l=0; l<k; l++)
237#elif _FOV_
238#define FOR_0_LE_l_LT_pk1 for (l=0; l<p; l++)
239#define FOR_0_LE_l_LT_pk for (l=0; l<k; l++)
240#elif _HOS_
241#define FOR_0_LE_l_LT_pk1 for (l=0; l<k1; l++)
242#define FOR_0_LE_l_LT_pk for (l=0; l<k; l++)
243#elif _HOS_OV_
244#define FOR_0_LE_l_LT_pk1 for (l=0; l<pk1; l++)
245#define FOR_0_LE_l_LT_pk for (l=0; l<p*k; l++)
246#else
247#define FOR_0_LE_l_LT_pk1
248#define FOR_0_LE_l_LT_pk
249#endif
250
251/*--------------------------------------------------------------------------*/
252/*                                                         VEC_COMPUTED_* */
253#ifdef _ADOLC_VECTOR
254#define VEC_COMPUTED_INIT   computed = 0;
255#define VEC_COMPUTED_CHECK  if (computed == 0) { computed = 1;
256#define VEC_COMPUTED_END    }
257#else
258#define VEC_COMPUTED_INIT
259#define VEC_COMPUTED_CHECK
260#define VEC_COMPUTED_END
261#endif
262
263/* END Macros */
264
265
266/****************************************************************************/
267/*                                                       NECESSARY INCLUDES */
268#include <adolc/interfaces.h>
269#include <adolc/adalloc.h>
270#include "oplate.h"
271#include "taping_p.h"
272#include <adolc/convolut.h>
273#include "dvlparms.h"
274
275#include <math.h>
276
277#if defined(ADOLC_DEBUG)
278#include <string.h>
279#endif /* ADOLC_DEBUG */
280
281BEGIN_C_DECLS
282
283/*--------------------------------------------------------------------------*/
284/*                                                   Local Static Variables */
285
286#ifdef _HOS_
287/***************************************************************************/
288/* Higher Order Scalar Reverse Pass.                                       */
289/***************************************************************************/
290int hos_reverse(short   tnum,        /* tape id */
291                int     depen,       /* consistency chk on # of deps */
292                int     indep,       /* consistency chk on # of indeps */
293                int     degre,       /* highest derivative degre  */
294                double  *lagrange,   /* range weight vector       */
295                double  **results)   /* matrix of coefficient vectors */
296{ int i, j, rc;
297    double** L = myalloc2(depen,degre+1);
298    for ( i = 0; i < depen; ++i ) {
299        L[i][0] = lagrange[i];
300        for ( j = 1; j <= degre; ++j )
301            L[i][j] = 0.0;
302    }
303    rc = hos_ti_reverse(tnum,depen,indep,degre,L,results);
304    myfree2(L);
305    return rc;
306}
307
308int hos_ti_reverse(
309    short   tnum,        /* tape id */
310    int     depen,       /* consistency chk on # of deps */
311    int     indep,       /* consistency chk on # of indeps */
312    int     degre,       /* highest derivative degre  */
313    double  **lagrange,  /* range weight vectors       */
314    double  **results)   /* matrix of coefficient vectors */
315
316#elif _HOS_OV_
317
318/***************************************************************************/
319/* Higher Order Scalar Reverse Pass, Vector Keep.                          */
320/***************************************************************************/
321int hos_ov_reverse(short   tnum,       /* tape id */
322                   int     depen,       /* consistency chk on # of deps */
323                   int     indep,       /* consistency chk on # of indeps */
324                   int     degre,       /* highest derivative degre  */
325                   int     nrows,       /* # of Jacobian rows calculated */
326                   double  **lagrange,  /* range weight vector       */
327                   double  ***results)  /* matrix of coefficient vectors */
328
329#elif _HOV_
330/***************************************************************************/
331/* Higher Order Vector Reverse Pass.                                       */
332/***************************************************************************/
333int hov_reverse(short   tnum,        /* tape id */
334                int     depen,       /* consistency chk on # of deps */
335                int     indep,       /* consistency chk on # of indeps */
336                int     degre,       /* highest derivative degre */
337                int     nrows,       /* # of Jacobian rows calculated */
338                double  **lagrange,  /* domain weight vector */
339                double  ***results,  /* matrix of coefficient vectors */
340                short   **nonzero )  /* structural sparsity  pattern  */
341{ int i, j, k, rc;
342    double*** L = myalloc3(nrows,depen,degre+1);
343    for ( k = 0; k < nrows; ++k )
344        for ( i = 0; i < depen; ++i ) {
345            L[k][i][0] = lagrange[k][i];
346            for ( j = 1; j <= degre; ++j )
347                L[k][i][j] = 0.0;
348        }
349    rc = hov_ti_reverse(tnum,depen,indep,degre,nrows,L,results,nonzero);
350    myfree3(L);
351    return rc;
352}
353
354int hov_ti_reverse(
355    short   tnum,        /* tape id */
356    int     depen,       /* consistency chk on # of deps */
357    int     indep,       /* consistency chk on # of indeps */
358    int     degre,       /* highest derivative degre */
359    int     nrows,       /* # of Jacobian rows calculated */
360    double  ***lagrange, /* domain weight vectors */
361    double  ***results,  /* matrix of coefficient vectors */
362    short   **nonzero )  /* structural sparsity  pattern  */
363
364#endif
365
366{
367    /************************************************************************/
368    /*                                                       ALL VARIABLES  */
369    unsigned char operation;   /* operation code */
370    int dc, ret_c=3;
371
372    locint size = 0;
373    locint res  = 0;
374    locint arg  = 0;
375    locint arg1 = 0;
376    locint arg2 = 0;
377
378    double coval = 0;
379
380    int indexi = 0,  indexd = 0;
381
382    /* loop indices */
383    int i, j, l, ls;
384
385    /* other necessary variables */
386    double *x;
387    int *jj;
388    int taycheck;
389    int numdep,numind;
390    double aTmp;
391
392    /*----------------------------------------------------------------------*/
393    /* Taylor stuff */
394    revreal *Tres, *Targ, *Targ1, *Targ2, *Tqo, *rp_Ttemp, *rp_Ttemp2;
395    revreal **rpp_T;
396
397    /*----------------------------------------------------------------------*/
398    /* Adjoint stuff */
399#ifdef _FOS_
400    double Atemp;
401# define A_TEMP Atemp
402#endif
403    revreal *Ares, *Aarg=NULL, *Aarg1, *Aarg2, *Aqo, *rp_Atemp, *rp_Atemp2;
404    revreal **rpp_A, *AP1, *AP2;
405
406    /*----------------------------------------------------------------------*/
407    int k = degre + 1;
408    int k1 = k + 1;
409    revreal comp;
410
411#ifdef _ADOLC_VECTOR_
412    int p = nrows;
413#endif
414
415#ifdef _HOV_
416    int pk1 = p*k1;
417    int q = 1;
418#elif _HOS_OV_
419    int p = nrows;
420    int pk1 = p*k1;
421    int q = p;
422#else
423    int q = 1;
424#endif
425        locint qq;
426
427    ADOLC_OPENMP_THREAD_NUMBER;
428    ADOLC_OPENMP_GET_THREAD_NUMBER;
429
430#if defined(ADOLC_DEBUG)
431    /************************************************************************/
432    /*                                                       DEBUG MESSAGES */
433    fprintf(DIAG_OUT,"Call of %s(..) with tag: %d, n: %d, m %d,\n",
434            GENERATED_FILENAME, tnum, indep, depen);
435
436#ifdef _HIGHER_ORDER_
437    fprintf(DIAG_OUT,"                    degree: %d\n",degre);
438#endif
439#ifdef _ADOLC_VECTOR_
440    fprintf(DIAG_OUT,"                    p: %d\n\n",nrows);
441#endif
442
443#endif
444
445
446    /************************************************************************/
447    /*                                                                INITs */
448
449    /*----------------------------------------------------------------------*/
450    /* Set up stuff for the tape */
451
452    /* Initialize the Reverse Sweep */
453    init_rev_sweep(tnum);
454
455    if ( (depen != ADOLC_CURRENT_TAPE_INFOS.stats[NUM_DEPENDENTS]) ||
456            (indep != ADOLC_CURRENT_TAPE_INFOS.stats[NUM_INDEPENDENTS]) )
457        fail(ADOLC_REVERSE_COUNTS_MISMATCH);
458
459    indexi = ADOLC_CURRENT_TAPE_INFOS.stats[NUM_INDEPENDENTS] - 1;
460    indexd = ADOLC_CURRENT_TAPE_INFOS.stats[NUM_DEPENDENTS] - 1;
461
462
463    /************************************************************************/
464    /*                                              MEMORY ALLOCATION STUFF */
465
466    /*----------------------------------------------------------------------*/
467#ifdef _HOS_                                                         /* HOS */
468    rpp_A = (revreal**)malloc(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES] *
469            sizeof(revreal*));
470    if (rpp_A == NULL) fail(ADOLC_MALLOC_FAILED);
471    Aqo = (revreal*)malloc(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES] * k1 *
472            sizeof(revreal));
473    if (Aqo == NULL) fail(ADOLC_MALLOC_FAILED);
474    for (i=0; i<ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES]; i++) {
475        rpp_A[i] = Aqo;
476        Aqo += k1;
477    }
478    rpp_T = (revreal**)malloc(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES] *
479            sizeof(revreal*));
480    if (rpp_T == NULL) fail(ADOLC_MALLOC_FAILED);
481    Tqo = (revreal *)malloc(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES] *
482            k * sizeof(revreal));
483    if (Tqo ==NULL) fail(ADOLC_MALLOC_FAILED);
484    for (i=0; i<ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES]; i++) {
485        rpp_T[i] = Tqo;
486        Tqo += k;
487    }
488    rp_Atemp  = (revreal *)malloc(k1 * sizeof(revreal));
489    rp_Atemp2 = (revreal *)malloc(k1 * sizeof(revreal));
490    rp_Ttemp2 = (revreal *)malloc(k * sizeof(revreal));
491    ADOLC_CURRENT_TAPE_INFOS.workMode = ADOLC_HOS_REVERSE;
492    /*----------------------------------------------------------------------*/
493#elif _HOV_                                                          /* HOV */
494    rpp_A = (revreal**)malloc(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES] *
495            sizeof(revreal*));
496    if (rpp_A == NULL) fail(ADOLC_MALLOC_FAILED);
497    Aqo = (revreal*)malloc(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES] * pk1 *
498            sizeof(revreal));
499    if (Aqo == NULL) fail(ADOLC_MALLOC_FAILED);
500    for (i=0; i<ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES]; i++) {
501        rpp_A[i] = Aqo;
502        Aqo += pk1;
503    }
504    rpp_T = (revreal**)malloc(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES] *
505            sizeof(revreal*));
506    if (rpp_T == NULL) fail(ADOLC_MALLOC_FAILED);
507    Tqo = (revreal*)malloc(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES] *
508            k * sizeof(revreal));
509    if (Tqo == NULL) fail(ADOLC_MALLOC_FAILED);
510    for (i=0; i<ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES]; i++) {
511        rpp_T[i] = Tqo;
512        Tqo += k;
513    }
514    rp_Atemp  = (revreal*) malloc(pk1 * sizeof(revreal));
515    rp_Atemp2 = (revreal*) malloc(pk1 * sizeof(revreal));
516    rp_Ttemp2 = (revreal*) malloc(k * sizeof(revreal));
517    ADOLC_CURRENT_TAPE_INFOS.workMode = ADOLC_HOV_REVERSE;
518    /*----------------------------------------------------------------------*/
519#elif _HOS_OV_                                                    /* HOS_OV */
520    rpp_A = (revreal**)malloc(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES] *
521            sizeof(revreal*));
522    if (rpp_A == NULL) fail(ADOLC_MALLOC_FAILED);
523    Aqo = (revreal*)malloc(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES] * pk1 *
524            sizeof(revreal));
525    if (Aqo == NULL) fail(ADOLC_MALLOC_FAILED);
526    for (i=0; i<ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES]; i++) {
527        rpp_A[i] = Aqo;
528        Aqo += pk1;
529    }
530    rpp_T = (revreal**)malloc(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES] *
531            sizeof(revreal*));
532    if (rpp_T == NULL) fail(ADOLC_MALLOC_FAILED);
533    Tqo = (revreal*)malloc(ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES] *
534            p * k * sizeof(revreal));
535    if (Tqo == NULL) fail(ADOLC_MALLOC_FAILED);
536    for (i=0; i<ADOLC_CURRENT_TAPE_INFOS.stats[NUM_MAX_LIVES]; i++) {
537        rpp_T[i] = Tqo;
538        Tqo += p*k;
539    }
540    rp_Atemp  = (revreal*) malloc(pk1 * sizeof(revreal));
541    rp_Atemp2 = (revreal*) malloc(pk1 * sizeof(revreal));
542    rp_Ttemp2 = (revreal*) malloc(p*k*sizeof(revreal));
543    ADOLC_CURRENT_TAPE_INFOS.workMode = ADOLC_HOV_REVERSE;
544#endif
545    rp_Ttemp  = (revreal*) malloc(k*sizeof(revreal));
546    if (rp_Ttemp == NULL) fail(ADOLC_MALLOC_FAILED);
547    if (rp_Ttemp2 == NULL) fail(ADOLC_MALLOC_FAILED);
548    x = myalloc1(q);
549    jj = (int*)malloc(q*sizeof(int));
550    if (jj == NULL) fail(ADOLC_MALLOC_FAILED);
551
552    /************************************************************************/
553    /*                                                TAYLOR INITIALIZATION */
554    ADOLC_CURRENT_TAPE_INFOS.rpp_A = rpp_A;
555    ADOLC_CURRENT_TAPE_INFOS.rpp_T = rpp_T;
556    taylor_back(tnum,&numdep,&numind,&taycheck);
557
558    if(taycheck != degre) {
559        fprintf(DIAG_OUT,"\n ADOL-C error: reverse fails because it was not"
560                " preceded\nby a forward sweep with degree>%i,"
561                " keep=%i!\n",degre,degre+1);
562        adolc_exit(-2,"",__func__,__FILE__,__LINE__);
563    };
564
565    if((numdep != depen)||(numind != indep)) {
566        fprintf(DIAG_OUT,"\n ADOL-C error: reverse fails on tape %d because "
567                "the number of\nindependent and/or dependent variables"
568                " given to reverse are\ninconsistent with that of the"
569                "  internal taylor array.\n",tnum);
570        adolc_exit(-2,"",__func__,__FILE__,__LINE__);
571    }
572
573
574    /************************************************************************/
575    /*                                                        REVERSE SWEEP */
576
577#if defined(ADOLC_DEBUG)
578    int v = 0;
579    unsigned int countPerOperation[256], taylorPerOperation[256];
580    memset(countPerOperation, 0, 1024);
581    memset(taylorPerOperation, 0, 1024);
582#   define UPDATE_TAYLORREAD(X) taylorPerOperation[operation] += X;
583#else
584#   define UPDATE_TAYLORREAD(X)
585#endif /* ADOLC_DEBUG */
586
587    operation=get_op_r();
588#if defined(ADOLC_DEBUG)
589    ++countPerOperation[operation];
590#endif /* ADOLC_DEBUG */
591
592    while (operation != start_of_tape) { 
593        /* Switch statement to execute the operations in Reverse */
594        switch (operation) {
595
596
597                /************************************************************/
598                /*                                                  MARKERS */
599
600                /*----------------------------------------------------------*/
601            case end_of_op:                                    /* end_of_op */
602                get_op_block_r();
603                operation = get_op_r();
604                /* Skip next operation, it's another end_of_op */
605                break;
606
607                /*----------------------------------------------------------*/
608            case end_of_int:                                  /* end_of_int */
609                get_loc_block_r(); /* Get the next int block */
610                break;
611
612                /*----------------------------------------------------------*/
613            case end_of_val:                                  /* end_of_val */
614                get_val_block_r(); /* Get the next val block */
615                break;
616
617                /*----------------------------------------------------------*/
618            case start_of_tape:                            /* start_of_tape */
619                break;
620            case end_of_tape:                                /* end_of_tape */
621                discard_params_r();
622                break;
623
624
625                /************************************************************/
626                /*                                               COMPARISON */
627
628                /*----------------------------------------------------------*/
629            case eq_zero  :                                      /* eq_zero */
630                arg   = get_locint_r();
631
632                ret_c = 0;
633                break;
634
635                /*----------------------------------------------------------*/
636            case neq_zero :                                     /* neq_zero */
637            case gt_zero  :                                      /* gt_zero */
638            case lt_zero :                                       /* lt_zero */
639                arg   = get_locint_r();
640                break;
641
642                /*----------------------------------------------------------*/
643            case ge_zero :                                       /* ge_zero */
644            case le_zero :                                       /* le_zero */
645                arg   = get_locint_r();
646
647                if (*rpp_T[arg] == 0)
648                    ret_c = 0;
649                break;
650
651
652                /************************************************************/
653                /*                                              ASSIGNMENTS */
654
655                /*----------------------------------------------------------*/
656            case assign_a:     /* assign an adouble variable an    assign_a */
657                /* adouble value. (=) */
658                res = get_locint_r();
659                arg = get_locint_r();
660
661                ASSIGN_A(Aarg, rpp_A[arg])
662                ASSIGN_A(Ares, rpp_A[res])
663
664                FOR_0_LE_l_LT_p
665                if  (0 == ARES) {
666                    HOV_INC(Aarg, k1)
667                    HOV_INC(Ares, k1)
668                } else {
669                    MAXDEC(AARG,ARES);
670                    AARG_INC_O;
671                    ARES_INC = 0.0;
672                    FOR_0_LE_i_LT_k
673                    { /* ! no tempory */
674                        AARG_INC += ARES;
675                        ARES_INC = 0.0;
676                    }
677                }
678
679                GET_TAYL(res,k,p)
680                break;
681
682                /*----------------------------------------------------------*/
683            case assign_d:      /* assign an adouble variable a    assign_d */
684                /* double value. (=) */
685                res   = get_locint_r();
686                coval = get_val_r();
687
688                ASSIGN_A(Ares, rpp_A[res])
689
690                FOR_0_LE_l_LT_pk1
691                ARES_INC = 0.0;
692
693                GET_TAYL(res,k,p)
694                break;
695
696                /*----------------------------------------------------------*/
697            case neg_sign_p:
698            case recipr_p:
699            case assign_p:      /* assign an adouble variable a    assign_d */
700                /* double value. (=) */
701                res   = get_locint_r();
702                arg   = get_locint_r();
703                coval = ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.paramstore[arg];
704
705                ASSIGN_A(Ares, rpp_A[res])
706
707                FOR_0_LE_l_LT_pk1
708                ARES_INC = 0.0;
709
710                GET_TAYL(res,k,p)
711                break;
712
713                /*----------------------------------------------------------*/
714            case assign_d_zero: /* assign an adouble a        assign_d_zero */
715            case assign_d_one:  /* double value. (=)           assign_d_one */
716                res   = get_locint_r();
717
718                ASSIGN_A(Ares, rpp_A[res])
719
720                FOR_0_LE_l_LT_pk1
721                ARES_INC = 0.0;
722
723                GET_TAYL(res,k,p)
724                break;
725
726                /*--------------------------------------------------------------------------*/
727            case assign_ind:       /* assign an adouble variable an    assign_ind */
728                /* independent double value (<<=) */
729                res = get_locint_r();
730
731                ASSIGN_A(Ares, rpp_A[res])
732
733                FOR_0_LE_l_LT_p
734                {
735#ifdef _HOV_
736                    if (nonzero) /* ??? question: why here? */
737                    nonzero[l][indexi] = (int)ARES;
738#endif /* _HOV_ */
739                    ARES_INC_O;
740                    FOR_0_LE_i_LT_k
741                        RESULTS(l,indexi,i) = ARES_INC;
742                }
743
744                GET_TAYL(res,k,p)
745                    indexi--;
746                break;
747
748                /*--------------------------------------------------------------------------*/
749            case assign_dep:           /* assign a float variable a    assign_dep */
750                /* dependent adouble value. (>>=) */
751                res = get_locint_r();
752
753                ASSIGN_A(Ares, rpp_A[res])
754                ASSIGN_A(Aarg, rpp_A[res])   /* just a helpful pointers */
755
756                FOR_0_LE_l_LT_p
757                { ARES_INC_O;
758                  dc = -1;
759                  FOR_0_LE_i_LT_k
760                  { ARES_INC = LAGRANGE(l,indexd,i);
761                    if (LAGRANGE(l,indexd,i)) dc = i;
762                  }
763                  AARG = (dc < 0)? 0.0 : (dc > 0)? 2.0 : 1.0;
764                  HOV_INC(Aarg, k1)
765                }
766                indexd--;
767            break;
768
769
770            /****************************************************************************/
771            /*                                                   OPERATION + ASSIGNMENT */
772
773            /*--------------------------------------------------------------------------*/
774        case eq_plus_d:            /* Add a floating point to an    eq_plus_d */
775            /* adouble. (+=) */
776            res   = get_locint_r();
777                coval = get_val_r();
778
779                GET_TAYL(res,k,p)
780                break;
781
782                /*--------------------------------------------------------------------------*/
783        case eq_plus_p:            /* Add a floating point to an    eq_plus_p */
784            /* adouble. (+=) */
785           res   = get_locint_r();
786           arg   = get_locint_r();
787                coval = ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.paramstore[arg];
788
789                GET_TAYL(res,k,p)
790                break;
791
792                /*--------------------------------------------------------------------------*/
793            case eq_plus_a:             /* Add an adouble to another    eq_plus_a */
794                /* adouble. (+=) */
795                res = get_locint_r();
796                arg = get_locint_r();
797
798                ASSIGN_A(Ares, rpp_A[res])
799                ASSIGN_A(Aarg, rpp_A[arg]);
800
801                FOR_0_LE_l_LT_p
802                if  (0 == ARES) {
803                    HOV_INC(Ares, k1)
804                    HOV_INC(Aarg, k1)
805                } else {
806                    MAXDEC(AARG,ARES);
807                    AARG_INC_O;
808                    ARES_INC_O;
809                    FOR_0_LE_i_LT_k
810                    AARG_INC += ARES_INC;
811                }
812
813                GET_TAYL(res,k,p)
814                break;
815
816                /*--------------------------------------------------------------------------*/
817            case eq_min_d:       /* Subtract a floating point from an    eq_min_d */
818                /* adouble. (-=) */
819                res   = get_locint_r();
820                coval = get_val_r();
821
822                GET_TAYL(res,k,p)
823                break;
824
825                /*--------------------------------------------------------------------------*/
826            case eq_min_p:       /* Subtract a floating point from an    eq_min_p */
827                /* adouble. (-=) */
828                res   = get_locint_r();
829                arg   = get_locint_r();
830                coval = ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.paramstore[arg];
831
832                GET_TAYL(res,k,p)
833                break;
834
835                /*--------------------------------------------------------------------------*/
836            case eq_min_a:        /* Subtract an adouble from another    eq_min_a */
837                /* adouble. (-=) */
838                res = get_locint_r();
839                arg = get_locint_r();
840
841                ASSIGN_A(Ares, rpp_A[res])
842                ASSIGN_A(Aarg, rpp_A[arg])
843
844                FOR_0_LE_l_LT_p
845                if  (0==ARES) {
846                    HOV_INC(Ares, k1)
847                    HOV_INC(Aarg, k1)
848                } else {
849                    MAXDEC(AARG,ARES);
850                    AARG_INC_O;
851                    ARES_INC_O;
852                    FOR_0_LE_i_LT_k
853                    AARG_INC -= ARES_INC;
854                }
855
856                GET_TAYL(res,k,p)
857                break;
858
859                /*--------------------------------------------------------------------------*/
860            case eq_mult_d:              /* Multiply an adouble by a    eq_mult_d */
861                /* flaoting point. (*=) */
862                res   = get_locint_r();
863                coval = get_val_r();
864
865                ASSIGN_A(Ares, rpp_A[res])
866
867                FOR_0_LE_l_LT_p
868                if ( 0 == ARES_INC )
869                    HOV_INC(Ares, k)
870                    else
871                        FOR_0_LE_i_LT_k
872                        ARES_INC *= coval;
873
874                GET_TAYL(res,k,p)
875                break;
876
877                /*--------------------------------------------------------------------------*/
878            case eq_mult_p:              /* Multiply an adouble by a    eq_mult_p */
879                /* flaoting point. (*=) */
880                res   = get_locint_r();
881                arg   = get_locint_r();
882                coval = ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.paramstore[arg];
883
884                ASSIGN_A(Ares, rpp_A[res])
885
886                FOR_0_LE_l_LT_p
887                if ( 0 == ARES_INC )
888                    HOV_INC(Ares, k)
889                    else
890                        FOR_0_LE_i_LT_k
891                        ARES_INC *= coval;
892
893                GET_TAYL(res,k,p)
894                break;
895
896                /*--------------------------------------------------------------------------*/
897            case eq_mult_a:       /* Multiply one adouble by another    eq_mult_a */
898                /* (*=) */
899                res = get_locint_r();
900                arg = get_locint_r();
901
902                GET_TAYL(res,k,p)
903
904                ASSIGN_A(Ares, rpp_A[res])
905                ASSIGN_A(Aarg, rpp_A[arg])
906                ASSIGN_A(Aqo,  rp_Atemp)
907                ASSIGN_T(Tres, rpp_T[res])
908                ASSIGN_T(Targ, rpp_T[arg])
909
910                FOR_0_LE_l_LT_p {
911                    if (0 == ARES) {
912                    HOV_INC(Aarg, k1)
913                        HOV_INC(Ares, k1)
914                    } else {
915                        MAXDEC(ARES,2.0);
916                        MAXDEC(AARG,ARES);
917                        AARG_INC_O;
918                        ARES_INC_O;
919                        conv(k,Ares,Targ,rp_Atemp);
920                        if(arg != res) {
921                            inconv(k,Ares,Tres,Aarg);
922                            FOR_0_LE_i_LT_k
923                            ARES_INC = AQO_INC;
924                        } else
925                            FOR_0_LE_i_LT_k
926                            ARES_INC = 2.0 * AQO_INC;
927                        HOV_INC(Aarg,k)
928                        HOS_OV_INC(Tres,k)
929                        HOS_OV_INC(Targ,k)
930                        HOS_OV_ASSIGN_A(Aqo,  rp_Atemp)
931                    }
932            }
933                break;
934
935                /*--------------------------------------------------------------------------*/
936            case incr_a:                        /* Increment an adouble    incr_a */
937            case decr_a:                        /* Increment an adouble    decr_a */
938                res   = get_locint_r();
939
940                GET_TAYL(res,k,p)
941                break;
942
943
944                /****************************************************************************/
945                /*                                                        BINARY OPERATIONS */
946
947                /*--------------------------------------------------------------------------*/
948            case plus_a_a:                 /* : Add two adoubles. (+)    plus a_a */
949                res  = get_locint_r();
950                arg2 = get_locint_r();
951                arg1 = get_locint_r();
952
953                ASSIGN_A(Ares,  rpp_A[res])
954                ASSIGN_A(Aarg1, rpp_A[arg1])
955                ASSIGN_A(Aarg2, rpp_A[arg2])
956
957                FOR_0_LE_l_LT_p
958                if  (0 == ARES) {
959                    HOV_INC(Ares,  k1)
960                    HOV_INC(Aarg1, k1)
961                    HOV_INC(Aarg2, k1)
962                } else {
963                    aTmp = ARES;
964                    ARES_INC = 0.0;
965                    MAXDEC(AARG1,aTmp);
966                    MAXDEC(AARG2,aTmp);
967                    AARG2_INC_O;
968                    AARG1_INC_O;
969                    FOR_0_LE_i_LT_k
970                    { aTmp = ARES;
971                      ARES_INC = 0.0;
972                      AARG1_INC += aTmp;
973                      AARG2_INC += aTmp;
974                    }
975                }
976
977                GET_TAYL(res,k,p)
978                break;
979
980                /*--------------------------------------------------------------------------*/
981            case plus_a_p:             /* Add an adouble and a double    plus_a_p */
982            case min_a_p:                /* Subtract an adouble from a    min_a_p */
983                /* (+) */
984                res   = get_locint_r();
985                arg1  = get_locint_r();
986                arg   = get_locint_r();
987                coval = ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.paramstore[arg1];
988
989                ASSIGN_A(Ares, rpp_A[res])
990                ASSIGN_A(Aarg, rpp_A[arg])
991
992                FOR_0_LE_l_LT_p
993                if  (0 == ARES) {
994                    HOV_INC(Ares, k1)
995                    HOV_INC(Aarg, k1)
996                } else {
997                    aTmp = ARES;
998                    ARES_INC = 0.0;
999                    MAXDEC(AARG,aTmp);
1000                    AARG_INC_O;
1001                    FOR_0_LE_i_LT_k
1002                    { aTmp = ARES;
1003                      ARES_INC = 0.0;
1004                      AARG_INC += aTmp;
1005                    }
1006                }
1007
1008                GET_TAYL(res,k,p)
1009                break;
1010
1011                /*--------------------------------------------------------------------------*/
1012            case plus_d_a:             /* Add an adouble and a double    plus_d_a */
1013                /* (+) */
1014                res   = get_locint_r();
1015                arg   = get_locint_r();
1016                coval = get_val_r();
1017
1018                ASSIGN_A(Ares, rpp_A[res])
1019                ASSIGN_A(Aarg, rpp_A[arg])
1020
1021                FOR_0_LE_l_LT_p
1022                if  (0 == ARES) {
1023                    HOV_INC(Ares, k1)
1024                    HOV_INC(Aarg, k1)
1025                } else {
1026                    aTmp = ARES;
1027                    ARES_INC = 0.0;
1028                    MAXDEC(AARG,aTmp);
1029                    AARG_INC_O;
1030                    FOR_0_LE_i_LT_k
1031                    { aTmp = ARES;
1032                      ARES_INC = 0.0;
1033                      AARG_INC += aTmp;
1034                    }
1035                }
1036
1037                GET_TAYL(res,k,p)
1038                break;
1039
1040                /*--------------------------------------------------------------------------*/
1041            case min_a_a:              /* Subtraction of two adoubles    min_a_a */
1042                /* (-) */
1043                res  = get_locint_r();
1044                arg2 = get_locint_r();
1045                arg1 = get_locint_r();
1046
1047                ASSIGN_A(Ares,  rpp_A[res])
1048                ASSIGN_A(Aarg1, rpp_A[arg1])
1049                ASSIGN_A(Aarg2, rpp_A[arg2])
1050
1051                FOR_0_LE_l_LT_p
1052                if  (0 == ARES) {
1053                    HOV_INC(Ares,  k1)
1054                    HOV_INC(Aarg1, k1)
1055                    HOV_INC(Aarg2, k1)
1056                } else {
1057                    aTmp = ARES;
1058                    ARES_INC = 0.0;
1059                    MAXDEC(AARG1,aTmp);
1060                    MAXDEC(AARG2,aTmp);
1061                    AARG2_INC_O;
1062                    AARG1_INC_O;
1063                    FOR_0_LE_i_LT_k
1064                    { aTmp = ARES;
1065                      ARES_INC = 0.0;
1066                      AARG1_INC += aTmp;
1067                      AARG2_INC -= aTmp;
1068                    }
1069                }
1070
1071                GET_TAYL(res,k,p)
1072                break;
1073
1074                /*--------------------------------------------------------------------------*/
1075            case min_d_a:                /* Subtract an adouble from a    min_d_a */
1076                /* double (-) */
1077                res   = get_locint_r();
1078                arg   = get_locint_r();
1079                coval = get_val_r();
1080
1081                ASSIGN_A(Ares, rpp_A[res])
1082                ASSIGN_A(Aarg, rpp_A[arg])
1083
1084                FOR_0_LE_l_LT_p
1085                if (0 == ARES) {
1086                    HOV_INC(Ares, k1)
1087                    HOV_INC(Aarg, k1)
1088                } else {
1089                    aTmp = ARES;
1090                    ARES_INC = 0.0;
1091                    MAXDEC(AARG,aTmp);
1092                    AARG_INC_O;
1093                    FOR_0_LE_i_LT_k
1094                    { aTmp = ARES;
1095                      ARES_INC = 0.0;
1096                      AARG_INC -= aTmp;
1097                    }
1098                }
1099
1100                GET_TAYL(res,k,p)
1101                break;
1102
1103                /*--------------------------------------------------------------------------*/
1104            case mult_a_a:               /* Multiply two adoubles (*)    mult_a_a */
1105                res  = get_locint_r();
1106                arg2 = get_locint_r();
1107                arg1 = get_locint_r();
1108
1109                GET_TAYL(res,k,p)
1110
1111                ASSIGN_A(Ares,  rpp_A[res])
1112                ASSIGN_A(Aarg2, rpp_A[arg2])
1113                ASSIGN_A(Aarg1, rpp_A[arg1])
1114                ASSIGN_T(Targ1, rpp_T[arg1])
1115                ASSIGN_T(Targ2, rpp_T[arg2])
1116
1117                FOR_0_LE_l_LT_p
1118                if (0 == ARES) {
1119                    HOV_INC(Aarg1, k1)
1120                    HOV_INC(Aarg2, k1)
1121                    HOV_INC(Ares,  k1)
1122                } else {
1123                    comp = (ARES > 2.0) ? ARES : 2.0 ;
1124                    ARES_INC = 0.0;
1125                    MAXDEC(AARG1,comp);
1126                    MAXDEC(AARG2,comp);
1127                    AARG1_INC_O;
1128                    AARG2_INC_O;
1129
1130                    copyAndZeroset(k,Ares,rp_Atemp);
1131                    inconv(k,rp_Atemp,Targ1,Aarg2);
1132                    inconv(k,rp_Atemp,Targ2,Aarg1);
1133
1134                    HOV_INC(Ares,  k)
1135                    HOV_INC(Aarg1, k)
1136                    HOV_INC(Aarg2, k)
1137                    HOS_OV_INC(Targ1, k)
1138                    HOS_OV_INC(Targ2, k)
1139                }
1140                break;
1141
1142                /*--------------------------------------------------------------------------*/
1143                /* olvo 991122: new op_code with recomputation */
1144            case eq_plus_prod:   /* increment a product of           eq_plus_prod */
1145                /* two adoubles (*) */
1146                res  = get_locint_r();
1147                arg2 = get_locint_r();
1148                arg1 = get_locint_r();
1149
1150
1151                ASSIGN_A(Ares,  rpp_A[res])
1152                ASSIGN_A(Aarg2, rpp_A[arg2])
1153                ASSIGN_A(Aarg1, rpp_A[arg1])
1154                ASSIGN_T(Targ1, rpp_T[arg1])
1155                ASSIGN_T(Targ2, rpp_T[arg2])
1156
1157                /* RECOMPUTATION */
1158                ASSIGN_T( Tres,  rpp_T[res])
1159#if !defined(_HOS_OV_)
1160                deconv1(k,Targ1,Targ2,Tres);
1161#endif
1162
1163                FOR_0_LE_l_LT_p {
1164#if defined(_HOS_OV_)
1165                deconv1(k,Targ1,Targ2,Tres);
1166#endif
1167                if (0 == ARES) {
1168                    HOV_INC(Aarg1, k1)
1169                    HOV_INC(Aarg2, k1)
1170                    HOV_INC(Ares,  k1)
1171                } else {
1172                    comp = (ARES > 2.0) ? ARES : 2.0 ;
1173                    ARES_INC = comp;
1174                    MAXDEC(AARG1,comp);
1175                    MAXDEC(AARG2,comp);
1176                    AARG1_INC_O;
1177                    AARG2_INC_O;
1178
1179                    inconv(k,Ares,Targ1,Aarg2);
1180                    inconv(k,Ares,Targ2,Aarg1);
1181
1182                    HOV_INC(Ares,  k)
1183                    HOV_INC(Aarg1, k)
1184                    HOV_INC(Aarg2, k)
1185                    HOS_OV_INC(Targ1, k)
1186                    HOS_OV_INC(Targ2, k)
1187                    HOS_OV_INC(Tres, k)
1188                }
1189                }
1190                break;
1191
1192                /*--------------------------------------------------------------------------*/
1193                /* olvo 991122: new op_code with recomputation */
1194            case eq_min_prod:   /* decrement a product of             eq_min_prod */
1195                /* two adoubles (*) */
1196                res  = get_locint_r();
1197                arg2 = get_locint_r();
1198                arg1 = get_locint_r();
1199
1200
1201                ASSIGN_A(Ares,  rpp_A[res])
1202                ASSIGN_A(Aarg2, rpp_A[arg2])
1203                ASSIGN_A(Aarg1, rpp_A[arg1])
1204                ASSIGN_T(Targ1, rpp_T[arg1])
1205                ASSIGN_T(Targ2, rpp_T[arg2])
1206
1207                /* RECOMPUTATION */
1208                ASSIGN_T( Tres,  rpp_T[res])
1209#if !defined(_HOS_OV_)
1210                inconv1(k,Targ1,Targ2,Tres);
1211#endif
1212
1213                FOR_0_LE_l_LT_p {
1214#if defined(_HOS_OV_)
1215                inconv1(k,Targ1,Targ2,Tres);
1216#endif
1217                if (0 == ARES) {
1218                    HOV_INC(Aarg1, k1)
1219                    HOV_INC(Aarg2, k1)
1220                    HOV_INC(Ares,  k1)
1221                } else {
1222                    comp = (ARES > 2.0) ? ARES : 2.0 ;
1223                    ARES_INC = comp;
1224                    MAXDEC(AARG1,comp);
1225                    MAXDEC(AARG2,comp);
1226                    AARG1_INC_O;
1227                    AARG2_INC_O;
1228
1229                    deconv1(k,Ares,Targ1,Aarg2);
1230                    deconv1(k,Ares,Targ2,Aarg1);
1231
1232                    HOV_INC(Ares,  k)
1233                    HOV_INC(Aarg1, k)
1234                    HOV_INC(Aarg2, k)
1235                    HOS_OV_INC(Targ1, k)
1236                    HOS_OV_INC(Targ2, k)
1237                    HOS_OV_INC(Tres, k)
1238                }
1239                }
1240                break;
1241
1242                /*--------------------------------------------------------------------------*/
1243            case mult_d_a:         /* Multiply an adouble by a double    mult_d_a */
1244                /* (*) */
1245                res   = get_locint_r();
1246                arg   = get_locint_r();
1247                coval = get_val_r();
1248
1249                ASSIGN_A(Ares, rpp_A[res])
1250                ASSIGN_A(Aarg, rpp_A[arg])
1251
1252                FOR_0_LE_l_LT_p
1253                if (0 == ARES) {
1254                    HOV_INC(Ares, k1)
1255                    HOV_INC(Aarg, k1)
1256                } else {
1257                    aTmp = ARES;
1258                    ARES_INC = 0.0;
1259                    MAXDEC(AARG,aTmp);
1260                    AARG_INC_O;
1261                    FOR_0_LE_i_LT_k
1262                    { aTmp = ARES;
1263                      ARES_INC = 0.0;
1264                      AARG_INC += coval * aTmp;
1265                    }
1266                }
1267
1268                GET_TAYL(res,k,p)
1269                break;
1270
1271                /*--------------------------------------------------------------------------*/
1272            case mult_a_p:         /* Multiply an adouble by a double    mult_a_p */
1273                /* (*) */
1274                res   = get_locint_r();
1275                arg1  = get_locint_r();
1276                arg   = get_locint_r();
1277                coval = ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.paramstore[arg1];
1278
1279                ASSIGN_A(Ares, rpp_A[res])
1280                ASSIGN_A(Aarg, rpp_A[arg])
1281
1282                FOR_0_LE_l_LT_p
1283                if (0 == ARES) {
1284                    HOV_INC(Ares, k1)
1285                    HOV_INC(Aarg, k1)
1286                } else {
1287                    aTmp = ARES;
1288                    ARES_INC = 0.0;
1289                    MAXDEC(AARG,aTmp);
1290                    AARG_INC_O;
1291                    FOR_0_LE_i_LT_k
1292                    { aTmp = ARES;
1293                      ARES_INC = 0.0;
1294                      AARG_INC += coval * aTmp;
1295                    }
1296                }
1297
1298                GET_TAYL(res,k,p)
1299                break;
1300
1301                /*--------------------------------------------------------------------------*/
1302            case div_a_a:           /* Divide an adouble by an adouble    div_a_a */
1303                /* (/) */
1304                res  = get_locint_r();
1305                arg2 = get_locint_r();
1306                arg1 = get_locint_r();
1307
1308                ASSIGN_A(Ares,  rpp_A[res])
1309                ASSIGN_A(Aarg2, rpp_A[arg2])
1310                ASSIGN_A(Aarg1, rpp_A[arg1])
1311                ASSIGN_T(Tres,  rpp_T[res])
1312                ASSIGN_T(Targ2, rpp_T[arg2])
1313
1314                /* olvo 980922 allows reflexive operation */
1315                if (arg2 == res) {
1316                    FOR_0_LE_l_LT_pk
1317                    rp_Ttemp2[l] = Tres[l];
1318                    Tres = rp_Ttemp2;
1319                    GET_TAYL(res,k,p)
1320                }
1321
1322                VEC_COMPUTED_INIT
1323                FOR_0_LE_l_LT_p
1324                { if (0 == ARES) {
1325                  HOV_INC(Ares,  k1)
1326                      HOV_INC(Aarg1, k1)
1327                      HOV_INC(Aarg2, k1)
1328                  } else {
1329                      aTmp = ARES;
1330                      ARES_INC = 0.0;
1331                      MAXDEC(AARG1,3.0);
1332                      MAXDEC(AARG1,aTmp);
1333                      MAXDEC(AARG2,3.0);
1334                      MAXDEC(AARG2,aTmp);
1335                      AARG1_INC_O;
1336                      AARG2_INC_O;
1337
1338                      VEC_COMPUTED_CHECK
1339                      recipr(k,1.0,Targ2,rp_Ttemp);
1340                      conv0(k ,rp_Ttemp,
1341                           Tres, rp_Atemp2);
1342                      VEC_COMPUTED_END
1343                      copyAndZeroset(k,Ares,rp_Atemp);
1344                      inconv(k, rp_Atemp,
1345                             rp_Ttemp, Aarg1);
1346                      deconv(k, rp_Atemp,
1347                             rp_Atemp2, Aarg2);
1348
1349                      HOV_INC(Ares,  k)
1350                      HOV_INC(Aarg1, k)
1351                      HOV_INC(Aarg2, k)
1352                      HOS_OV_INC(Tres, k)
1353                      HOS_OV_INC(Targ2, k)
1354                  }
1355            }
1356
1357                if (res != arg2)
1358                    GET_TAYL(res,k,p)
1359                    break;
1360
1361                /*--------------------------------------------------------------------------*/
1362            case div_d_a:             /* Division double - adouble (/)    div_d_a */
1363                res   = get_locint_r();
1364                arg   = get_locint_r();
1365                coval = get_val_r();
1366
1367                ASSIGN_A(Ares, rpp_A[res])
1368                ASSIGN_A(Aarg, rpp_A[arg])
1369                ASSIGN_T(Tres, rpp_T[res])
1370                ASSIGN_T(Targ, rpp_T[arg])
1371
1372                /* olvo 980922 allows reflexive operation */
1373                if (arg == res) {
1374                    FOR_0_LE_l_LT_pk
1375                    rp_Ttemp2[l] = Tres[l];
1376                    Tres = rp_Ttemp2;
1377                    GET_TAYL(arg,k,p)
1378                }
1379
1380                VEC_COMPUTED_INIT
1381                FOR_0_LE_l_LT_p
1382                { if (0 == ARES) {
1383                  HOV_INC(Ares, k1)
1384                      HOV_INC(Aarg, k1)
1385                  } else {
1386                      aTmp = ARES;
1387                      ARES_INC = 0.0;
1388                      MAXDEC(AARG,aTmp);
1389                      MAXDEC(AARG,3.0);
1390                      AARG_INC_O;
1391
1392                      VEC_COMPUTED_CHECK
1393                      recipr(k,1.0,Targ,rp_Ttemp);
1394                      conv0(k, rp_Ttemp,
1395                           Tres, rp_Atemp);
1396                      VEC_COMPUTED_END
1397                      deconvZeroR(k,Ares,rp_Atemp,Aarg);
1398
1399                      HOV_INC(Ares, k)
1400                      HOV_INC(Aarg, k)
1401                      HOS_OV_INC(Tres, k)
1402                      HOS_OV_INC(Targ, k)
1403                  }
1404            }
1405
1406                if (arg != res)
1407                GET_TAYL(res,k,p)
1408                break;
1409
1410
1411                /****************************************************************************/
1412            case div_p_a:             /* Division double - adouble (/)    div_p_a */
1413                res   = get_locint_r();
1414                arg1  = get_locint_r();
1415                arg   = get_locint_r();
1416                coval = ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.paramstore[arg1];
1417
1418                ASSIGN_A(Ares, rpp_A[res])
1419                ASSIGN_A(Aarg, rpp_A[arg])
1420                ASSIGN_T(Tres, rpp_T[res])
1421                ASSIGN_T(Targ, rpp_T[arg])
1422
1423                /* olvo 980922 allows reflexive operation */
1424                if (arg == res) {
1425                    FOR_0_LE_l_LT_pk
1426                    rp_Ttemp2[l] = Tres[l];
1427                    Tres = rp_Ttemp2;
1428                    GET_TAYL(arg,k,p)
1429                }
1430
1431                VEC_COMPUTED_INIT
1432                FOR_0_LE_l_LT_p
1433                { if (0 == ARES) {
1434                  HOV_INC(Ares, k1)
1435                      HOV_INC(Aarg, k1)
1436                  } else {
1437                      aTmp = ARES;
1438                      ARES_INC = 0.0;
1439                      MAXDEC(AARG,aTmp);
1440                      MAXDEC(AARG,3.0);
1441                      AARG_INC_O;
1442
1443                      VEC_COMPUTED_CHECK
1444                      recipr(k,1.0,Targ,rp_Ttemp);
1445                      conv0(k, rp_Ttemp,
1446                           Tres, rp_Atemp);
1447                      VEC_COMPUTED_END
1448                      deconvZeroR(k,Ares,rp_Atemp,Aarg);
1449
1450                      HOV_INC(Ares, k)
1451                      HOV_INC(Aarg, k)
1452                      HOS_OV_INC(Tres, k)
1453                      HOS_OV_INC(Targ, k)
1454                  }
1455            }
1456
1457                if (arg != res)
1458                GET_TAYL(res,k,p)
1459                break;
1460
1461
1462                /****************************************************************************/
1463                /*                                                         SIGN  OPERATIONS */
1464
1465                /*--------------------------------------------------------------------------*/
1466            case pos_sign_a:                                        /* pos_sign_a */
1467                res   = get_locint_r();
1468                arg   = get_locint_r();
1469
1470                ASSIGN_A(Ares, rpp_A[res])
1471                ASSIGN_A(Aarg, rpp_A[arg])
1472
1473                FOR_0_LE_l_LT_p
1474                if  (0 == ARES) {
1475                    HOV_INC(Ares, k1)
1476                    HOV_INC(Aarg, k1)
1477                } else {
1478                    aTmp = ARES;
1479                    ARES_INC = 0.0;
1480                    MAXDEC(AARG,aTmp);
1481                    AARG_INC_O;
1482                    FOR_0_LE_i_LT_k
1483                    { aTmp = ARES;
1484                      ARES_INC = 0.0;
1485                      AARG_INC += aTmp;
1486                    }
1487                }
1488
1489                GET_TAYL(res,k,p)
1490                break;
1491
1492                /*--------------------------------------------------------------------------*/
1493            case neg_sign_a:                                        /* neg_sign_a */
1494                res   = get_locint_r();
1495                arg   = get_locint_r();
1496
1497                ASSIGN_A(Ares, rpp_A[res])
1498                ASSIGN_A(Aarg, rpp_A[arg])
1499
1500                FOR_0_LE_l_LT_p
1501                if  (0 == ARES) {
1502                    HOV_INC(Ares, k1)
1503                    HOV_INC(Aarg, k1)
1504                } else {
1505                    aTmp = ARES;
1506                    ARES_INC = 0.0;
1507                    MAXDEC(AARG,aTmp);
1508                    AARG_INC_O;
1509                    FOR_0_LE_i_LT_k
1510                    { aTmp = ARES;
1511                      ARES_INC = 0.0;
1512                      AARG_INC -= aTmp;
1513                    }
1514                }
1515
1516                GET_TAYL(res,k,p)
1517                break;
1518
1519
1520                /****************************************************************************/
1521                /*                                                         UNARY OPERATIONS */
1522
1523                /*--------------------------------------------------------------------------*/
1524            case exp_op:                          /* exponent operation    exp_op */
1525                res = get_locint_r();
1526                arg = get_locint_r();
1527
1528                ASSIGN_A(Ares, rpp_A[res])
1529                ASSIGN_A(Aarg, rpp_A[arg])
1530                ASSIGN_T(Tres, rpp_T[res])
1531                ASSIGN_T(Targ, rpp_T[arg])
1532
1533                FOR_0_LE_l_LT_p
1534                { if (0 == ARES) {
1535                  HOV_INC(Aarg, k1)
1536                      HOV_INC(Ares, k1)
1537                  } else {
1538                      aTmp = ARES;
1539                      ARES_INC = 0.0;
1540                      MAXDEC(AARG,aTmp);
1541                      MAXDEC(AARG,4.0);
1542                      AARG_INC_O;
1543
1544                      inconv0(k,Ares,Tres,Aarg);
1545
1546                      HOV_INC(Ares, k)
1547                      HOV_INC(Aarg, k)
1548                      HOS_OV_INC(Tres, k)
1549                  }
1550            }
1551
1552                GET_TAYL(res,k,p)
1553                break;
1554
1555                /*--------------------------------------------------------------------------*/
1556            case sin_op:                              /* sine operation    sin_op */
1557                res  = get_locint_r();
1558                arg2 = get_locint_r();
1559                arg1 = get_locint_r();
1560
1561                ASSIGN_A(Ares,  rpp_A[res])
1562                ASSIGN_A(Aarg1, rpp_A[arg1])
1563                ASSIGN_T(Targ2, rpp_T[arg2])
1564
1565                FOR_0_LE_l_LT_p
1566                { if (0 == ARES) {
1567                  HOV_INC(Aarg1, k1)
1568                      HOV_INC(Ares,  k1)
1569                  } else {
1570                      aTmp = ARES;
1571                      ARES_INC = 0.0;
1572                      MAXDEC(AARG1,aTmp);
1573                      MAXDEC(AARG1,4.0);
1574                      AARG1_INC_O;
1575
1576                      inconv0(k,Ares,Targ2,Aarg1);
1577
1578                      HOV_INC(Ares,  k)
1579                      HOV_INC(Aarg1, k)
1580                      HOS_OV_INC(Targ2, k)
1581                  }
1582            }
1583
1584                GET_TAYL(res,k,p)
1585                GET_TAYL(arg2,k,p) /* olvo 980710 covalue */
1586                /* NOTE: rpp_A[arg2] should be 0 already */
1587                break;
1588
1589                /*--------------------------------------------------------------------------*/
1590            case cos_op:                            /* cosine operation    cos_op */
1591                res  = get_locint_r();
1592                arg2 = get_locint_r();
1593                arg1 = get_locint_r();
1594
1595                ASSIGN_A(Ares,  rpp_A[res])
1596                ASSIGN_A(Aarg1, rpp_A[arg1])
1597                ASSIGN_T(Targ2, rpp_T[arg2])
1598
1599                FOR_0_LE_l_LT_p
1600                { if (0 == ARES) {
1601                  HOV_INC(Aarg1, k1)
1602                      HOV_INC(Ares,  k1)
1603                  } else {
1604                      aTmp = ARES;
1605                      ARES_INC = 0.0;
1606                      MAXDEC(AARG1,aTmp);
1607                      MAXDEC(AARG1,4.0);
1608                      AARG1_INC_O;
1609
1610                      deconv0(k,Ares,Targ2,Aarg1);
1611
1612                      HOV_INC(Ares,  k)
1613                      HOV_INC(Aarg1, k)
1614                      HOS_OV_INC(Targ2, k)
1615                  }
1616            }
1617
1618                GET_TAYL(res,k,p)
1619                GET_TAYL(arg2,k,p) /* olvo 980710 covalue */
1620                /* NOTE: rpp_A[arg2] should be 0 already */
1621                break;
1622                /*xxx*/
1623                /*--------------------------------------------------------------------------*/
1624            case atan_op:                                             /* atan_op  */
1625            case asin_op:                                             /* asin_op  */
1626            case acos_op:                                             /* acos_op  */
1627            case asinh_op:                                            /* asinh_op */
1628            case acosh_op:                                            /* acosh_op */
1629            case atanh_op:                                            /* atanh_op */
1630            case erf_op:                                              /* erf_op   */
1631                res  = get_locint_r();
1632                arg2 = get_locint_r();
1633                arg1 = get_locint_r();
1634
1635                GET_TAYL(res,k,p)
1636
1637                ASSIGN_A(Ares,  rpp_A[res])
1638                ASSIGN_A(Aarg1, rpp_A[arg1])
1639                ASSIGN_T(Targ2, rpp_T[arg2])
1640
1641                FOR_0_LE_l_LT_p
1642                { if (0 == ARES) {
1643                  HOV_INC(Aarg1, k1)
1644                      HOV_INC(Ares,  k1)
1645                  } else {
1646                      aTmp = ARES;
1647                      ARES_INC = 0.0;
1648                      MAXDEC(AARG1,aTmp);
1649                      MAXDEC(AARG1,4.0);
1650                      AARG1_INC_O;
1651
1652                      inconv0(k,Ares,Targ2,Aarg1);
1653
1654                      HOV_INC(Aarg1, k)
1655                      HOV_INC(Ares,  k)
1656                      HOS_OV_INC(Targ2, k)
1657                  }
1658            }
1659                break;
1660
1661                /*--------------------------------------------------------------------------*/
1662            case log_op:                                                /* log_op */
1663                res = get_locint_r();
1664                arg = get_locint_r();
1665
1666                GET_TAYL(res,k,p)
1667
1668                ASSIGN_A(Ares, rpp_A[res])
1669                ASSIGN_A(Aarg, rpp_A[arg])
1670                ASSIGN_T(Targ, rpp_T[arg])
1671
1672                VEC_COMPUTED_INIT
1673                FOR_0_LE_l_LT_p
1674                { if (0 == ARES) {
1675                  HOV_INC(Aarg, k1)
1676                      HOV_INC(Ares, k1)
1677                  } else {
1678                      aTmp = ARES;
1679                      ARES_INC = 0.0;
1680                      MAXDEC(AARG,aTmp);
1681                      MAXDEC(AARG,4.0);
1682                      AARG_INC_O;
1683
1684                      VEC_COMPUTED_CHECK
1685                      recipr(k,1.0,Targ,rp_Ttemp);
1686                      VEC_COMPUTED_END
1687                      inconv0(k,Ares,rp_Ttemp,Aarg);
1688
1689                      HOV_INC(Ares, k)
1690                      HOV_INC(Aarg, k)
1691                      HOS_OV_INC(Targ, k)
1692                  }
1693            }
1694                break;
1695
1696                /*--------------------------------------------------------------------------*/
1697            case pow_op:                                                /* pow_op */
1698                res   = get_locint_r();
1699                arg   = get_locint_r();
1700                coval = get_val_r();
1701
1702                ASSIGN_T(Targ, rpp_T[arg])
1703                ASSIGN_T(Tres, rpp_T[res])
1704                ASSIGN_A(Ares, rpp_A[res])
1705                ASSIGN_A(Aarg, rpp_A[arg])
1706
1707                /* olvo 980921 allows reflexive operation */
1708                if (arg == res) {
1709                    FOR_0_LE_l_LT_pk
1710                    rp_Ttemp2[l] = Tres[l];
1711                    Tres = rp_Ttemp2;
1712                    GET_TAYL(arg,k,p)
1713                }
1714
1715                VEC_COMPUTED_INIT
1716                FOR_0_LE_l_LT_p
1717                if (0 == ARES) {
1718                    HOV_INC(Aarg, k1)
1719                    HOV_INC(Ares, k1)
1720                } else {
1721                    aTmp = ARES;
1722                    ARES_INC = 0.0;
1723                    MAXDEC(AARG,aTmp);
1724                    MAXDEC(AARG,4.0);
1725                    AARG_INC_O;
1726
1727                    VEC_COMPUTED_CHECK
1728                    if (fabs(Targ[0]) > ADOLC_EPS) {
1729                        divide(k,Tres,Targ,rp_Ttemp);
1730                        for (i=0;i<k;i++) {
1731                            rp_Ttemp[i] *= coval;
1732                            /*                 printf(" EPS i %d %f\n",i,rp_Ttemp[i]); */
1733                        }
1734                        inconv0(k,Ares,rp_Ttemp,Aarg);
1735                    } else {
1736                        if (coval <= 0.0) {
1737                            FOR_0_LE_i_LT_k
1738                            {
1739                                Aarg[i] = make_nan();
1740                                Ares[i] = 0;
1741                            }
1742                        } else {
1743                            /* coval not a whole number */
1744                            if (coval - floor(coval) != 0) {
1745                                i = 0;
1746                                FOR_0_LE_i_LT_k
1747                                {
1748                                    if (coval - i > 1) {
1749                                    Aarg[i] = 0;
1750                                        Ares[i] = 0;
1751                                    }
1752                                    if ((coval - i < 1) && (coval - i > 0)) {
1753                                    Aarg[i] = make_inf();
1754                                        Ares[i] = 0;
1755                                    }
1756                                    if (coval - i < 0) {
1757                                    Aarg[i] = make_nan();
1758                                        Ares[i] = 0;
1759                                    }
1760                                }
1761                            } else {
1762                                if (coval == 1) {
1763                                    FOR_0_LE_i_LT_k
1764                                    { /* ! no tempory */
1765                                        Aarg[i] += Ares[i];
1766                                        Ares[i] = 0.0;
1767                                    }
1768                                } else {
1769                                    /* coval is an int > 1 */
1770                                    /* the following is not efficient but at least it works */
1771                                    /* it reformulates x^n into x* ... *x n times */
1772
1773                                    copyAndZeroset(k,Ares,rp_Atemp);
1774                                    inconv(k,rp_Atemp,Targ,Aarg);
1775                                    inconv(k,rp_Atemp,Targ,Aarg);
1776                                    if (coval == 3) {
1777                                        conv(k,Aarg,Targ,rp_Atemp);
1778                                        FOR_0_LE_i_LT_k
1779                                        Aarg[i] = 2.0 * rp_Atemp[i];
1780                                   }
1781                                }
1782                            }
1783                        }
1784                    }
1785                    VEC_COMPUTED_END
1786
1787                    HOV_INC(Ares, k)
1788                    HOV_INC(Aarg, k)
1789                    HOS_OV_INC(Tres, k)
1790                    HOS_OV_INC(Targ, k)
1791                }
1792
1793                GET_TAYL(res,k,p)
1794                break;
1795
1796                /*--------------------------------------------------------------------------*/
1797            case pow_op_p:                                                /* pow_op_p */
1798                res   = get_locint_r();
1799                arg1  = get_locint_r();
1800                arg   = get_locint_r();
1801                coval = ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.paramstore[arg1];
1802
1803                ASSIGN_T(Targ, rpp_T[arg])
1804                ASSIGN_T(Tres, rpp_T[res])
1805                ASSIGN_A(Ares, rpp_A[res])
1806                ASSIGN_A(Aarg, rpp_A[arg])
1807
1808                /* olvo 980921 allows reflexive operation */
1809                if (arg == res) {
1810                    FOR_0_LE_l_LT_pk
1811                    rp_Ttemp2[l] = Tres[l];
1812                    Tres = rp_Ttemp2;
1813                    GET_TAYL(arg,k,p)
1814                }
1815
1816                VEC_COMPUTED_INIT
1817                FOR_0_LE_l_LT_p
1818                if (0 == ARES) {
1819                    HOV_INC(Aarg, k1)
1820                    HOV_INC(Ares, k1)
1821                } else {
1822                    aTmp = ARES;
1823                    ARES_INC = 0.0;
1824                    MAXDEC(AARG,aTmp);
1825                    MAXDEC(AARG,4.0);
1826                    AARG_INC_O;
1827
1828                    VEC_COMPUTED_CHECK
1829                    if (fabs(Targ[0]) > ADOLC_EPS) {
1830                        divide(k,Tres,Targ,rp_Ttemp);
1831                        for (i=0;i<k;i++) {
1832                            rp_Ttemp[i] *= coval;
1833                            /*                 printf(" EPS i %d %f\n",i,rp_Ttemp[i]); */
1834                        }
1835                        inconv0(k,Ares,rp_Ttemp,Aarg);
1836                    } else {
1837                        if (coval <= 0.0) {
1838                            FOR_0_LE_i_LT_k
1839                            {
1840                                Aarg[i] = make_nan();
1841                                Ares[i] = 0;
1842                            }
1843                        } else {
1844                            /* coval not a whole number */
1845                            if (coval - floor(coval) != 0) {
1846                                i = 0;
1847                                FOR_0_LE_i_LT_k
1848                                {
1849                                    if (coval - i > 1) {
1850                                    Aarg[i] = 0;
1851                                        Ares[i] = 0;
1852                                    }
1853                                    if ((coval - i < 1) && (coval - i > 0)) {
1854                                    Aarg[i] = make_inf();
1855                                        Ares[i] = 0;
1856                                    }
1857                                    if (coval - i < 0) {
1858                                    Aarg[i] = make_nan();
1859                                        Ares[i] = 0;
1860                                    }
1861                                }
1862                            } else {
1863                                if (coval == 1) {
1864                                    FOR_0_LE_i_LT_k
1865                                    { /* ! no tempory */
1866                                        Aarg[i] += Ares[i];
1867                                        Ares[i] = 0.0;
1868                                    }
1869                                } else {
1870                                    /* coval is an int > 1 */
1871                                    /* the following is not efficient but at least it works */
1872                                    /* it reformulates x^n into x* ... *x n times */
1873
1874                                    copyAndZeroset(k,Ares,rp_Atemp);
1875                                    inconv(k,rp_Atemp,Targ,Aarg);
1876                                    inconv(k,rp_Atemp,Targ,Aarg);
1877                                    if (coval == 3) {
1878                                        conv(k,Aarg,Targ,rp_Atemp);
1879                                        FOR_0_LE_i_LT_k
1880                                        Aarg[i] = 2.0 * rp_Atemp[i];
1881                                   }
1882                                }
1883                            }
1884                        }
1885                    }
1886                    VEC_COMPUTED_END
1887
1888                    HOV_INC(Ares, k)
1889                    HOV_INC(Aarg, k)
1890                    HOS_OV_INC(Tres, k)
1891                    HOS_OV_INC(Targ, k)
1892                }
1893
1894                GET_TAYL(res,k,p)
1895                break;
1896
1897                /*--------------------------------------------------------------------------*/
1898            case sqrt_op:                                              /* sqrt_op */
1899                res = get_locint_r();
1900                arg = get_locint_r();
1901
1902                ASSIGN_A(Ares, rpp_A[res])
1903                ASSIGN_A(Aarg, rpp_A[arg])
1904                ASSIGN_T(Tres, rpp_T[res])
1905
1906                VEC_COMPUTED_INIT
1907                FOR_0_LE_l_LT_p
1908                if (0 == ARES) {
1909                    HOV_INC(Aarg, k1)
1910                    HOV_INC(Ares, k1)
1911                } else {
1912                    aTmp = ARES;
1913                    ARES_INC = 0.0;
1914                    MAXDEC(AARG,aTmp);
1915                    MAXDEC(AARG,4.0);
1916                    AARG_INC_O;
1917
1918                    VEC_COMPUTED_CHECK
1919                    recipr(k,0.5,Tres,rp_Ttemp);
1920                    VEC_COMPUTED_END
1921                    inconv0(k,Ares,rp_Ttemp,Aarg);
1922
1923                    HOV_INC(Ares, k)
1924                    HOV_INC(Aarg, k)
1925                    HOS_OV_INC(Tres,k)
1926                }
1927
1928                GET_TAYL(res,k,p)
1929                break;
1930
1931                /*--------------------------------------------------------------------------*/
1932            case gen_quad:                                            /* gen_quad */
1933                res   = get_locint_r();
1934                arg2  = get_locint_r();
1935                arg1  = get_locint_r();
1936                coval = get_val_r();
1937                coval = get_val_r();
1938
1939                ASSIGN_A(Ares,  rpp_A[res])
1940                ASSIGN_A(Aarg1, rpp_A[arg1])
1941                ASSIGN_T(Targ2, rpp_T[arg2])
1942
1943                FOR_0_LE_l_LT_p
1944                if (0 == ARES) {
1945                    HOV_INC(Aarg1, k1)
1946                    HOV_INC(Ares,  k1)
1947                } else {
1948                    aTmp = ARES;
1949                    ARES_INC = 0.0;
1950                    MAXDEC(AARG1,aTmp);
1951                    MAXDEC(AARG1,4.0);
1952                    AARG1_INC_O;
1953
1954                    inconv0(k,Ares,Targ2,Aarg1);
1955
1956                    HOV_INC(Aarg1, k)
1957                    HOV_INC(Ares,  k)
1958                    HOS_OV_INC(Targ2,  k)
1959                }
1960
1961                GET_TAYL(res,k,p)
1962                break;
1963
1964                /*--------------------------------------------------------------------------*/
1965            case min_op:                                                /* min_op */
1966
1967#ifdef _HOS_OV_
1968
1969                fprintf(DIAG_OUT," operation min_op not implemented for hos_ov");
1970                break;
1971#endif
1972                res   = get_locint_r();
1973                arg2  = get_locint_r();
1974                arg1  = get_locint_r();
1975                coval = get_val_r();
1976
1977                GET_TAYL(res,k,p)
1978
1979                ASSIGN_A(Aarg1, rpp_A[arg1])
1980                ASSIGN_A(Aarg2, rpp_A[arg2])
1981                ASSIGN_A(Ares,  rpp_A[res])
1982                ASSIGN_T(Targ1, rpp_T[arg1])
1983                ASSIGN_T(Targ2, rpp_T[arg2])
1984                ASSIGN_A(AP1,   NULL)
1985                ASSIGN_A(AP2,   Ares)
1986
1987                if (Targ1[0] > Targ2[0]) {
1988                    FOR_0_LE_l_LT_p
1989                    { if ((coval) && (*AP2))
1990                      MINDEC(ret_c,2);
1991                      HOV_INC(AP2,k1)
1992                    }
1993                    AP1 = Aarg2;
1994                arg = 0;
1995            } else
1996                if (Targ1[0] < Targ2[0]) {
1997                        FOR_0_LE_l_LT_p
1998                        { if ((!coval) && (*AP2))
1999                          MINDEC(ret_c,2);
2000                          HOV_INC(AP2,k1)
2001                        }
2002                        AP1 = Aarg1;
2003                    arg = 0;
2004                } else /* both are equal */ /* must be changed for hos_ov, but how? */
2005                    /* seems to influence the return value */
2006                    for (i=1;i<k;i++) {
2007                            if (Targ1[i] > Targ2[i]) {
2008                                FOR_0_LE_l_LT_p
2009                                { if (*AP2)
2010                                  MINDEC(ret_c,1);
2011                                  HOV_INC(AP2,k1)
2012                                }
2013                                AP1 = Aarg2;
2014                            arg = i+1;
2015                        } else
2016                            if (Targ1[i] < Targ2[i]) {
2017                                    FOR_0_LE_l_LT_p
2018                                    { if (*AP2)
2019                                      MINDEC(ret_c,1);
2020                                      HOV_INC(AP2,k1)
2021                                    }
2022                                    AP1 = Aarg1;
2023                                arg = i+1;
2024                            }
2025                        if (AP1 != NULL)
2026                                break;
2027                        }
2028
2029                if (AP1 != NULL)
2030                    FOR_0_LE_l_LT_p
2031                    { if (0 == ARES) {
2032                      HOV_INC(AP1, k1)
2033                          HOV_INC(Ares,k1);
2034                      } else {
2035                          aTmp = ARES;
2036                          ARES_INC = 0.0;
2037                          if (arg)  /* we are at the tie */
2038                              *AP1 = 5.0;
2039                          else
2040                              MAXDEC(*AP1,aTmp);
2041                          AP1++;
2042                          for (i=0;i<k;i++) {
2043                              aTmp = ARES;
2044                              ARES_INC = 0.0;
2045                              *AP1++ += aTmp;
2046                          }
2047                      }
2048            }
2049                else /* both are identical */
2050                {
2051                    FOR_0_LE_l_LT_p
2052                    { if (0 == ARES) {
2053                      HOV_INC(Aarg1,k1)
2054                          HOV_INC(Aarg2,k1)
2055                          HOV_INC(Ares, k1)
2056                      } else {
2057                          aTmp = ARES;
2058                          ARES_INC = 0.0;
2059                          MAXDEC(AARG1,aTmp);  /*assume sthg like fmin(x,x) */
2060                          MAXDEC(AARG2,aTmp);
2061                          AARG1_INC_O;
2062                          AARG2_INC_O;
2063                          for (i=0;i<k;i++) {
2064                              aTmp = ARES;
2065                              ARES_INC = 0.0;
2066                              AARG1_INC += aTmp/2;
2067                              AARG2_INC += aTmp/2;
2068                          }
2069                      }
2070                    }
2071                    if (arg1 != arg2)
2072                        MINDEC(ret_c,1);
2073                }
2074                break;
2075
2076
2077                /*--------------------------------------------------------------------------*/
2078            case abs_val:                                              /* abs_val */
2079                res   = get_locint_r();
2080                arg   = get_locint_r();
2081                coval = get_val_r();
2082                /* must be changed for hos_ov, but how? */
2083                /* seems to influence the return value  */
2084                GET_TAYL(res,k,p)
2085
2086                ASSIGN_A(Ares, rpp_A[res])
2087                ASSIGN_A(Aarg, rpp_A[arg])
2088                ASSIGN_T(Targ, rpp_T[arg])
2089
2090                FOR_0_LE_l_LT_q
2091                {
2092                    x[l] = 0.0;
2093                    jj[l] = 0;
2094                    for (i=0;i<k;i++)
2095                    if ( (x[l] == 0.0) && (Targ[i] != 0.0) ) {
2096                        jj[l] = i;
2097                            if (Targ[i] < 0.0)
2098                                x[l] = -1.0;
2099                            else
2100                                x[l] = 1.0;
2101                        }
2102                    HOS_OV_INC(Targ,k)
2103            }
2104                ASSIGN_T(Targ, rpp_T[arg])
2105                FOR_0_LE_l_LT_p
2106                { if (0 == ARES) {
2107                  HOV_INC(Aarg, k1)
2108                      HOV_INC(Ares, k1)
2109                  } else {
2110                      if (Targ[0] == 0.0) {
2111                          ARES_INC = 0.0;
2112                          AARG_INC = 5.0;
2113                      } else {
2114                          aTmp = ARES;
2115                          ARES_INC = 0.0;
2116                          MAXDEC(AARG,aTmp);
2117                          AARG_INC_O;
2118                      }
2119                      if(Targ[0] == 0.0)
2120                          MINDEC(ret_c,1);
2121                      for (i=0;i<jj[l];i++)
2122                          ARES_INC = 0.0;
2123                      Aarg += jj[l];
2124                      for (i=jj[l];i<k;i++) {
2125                          aTmp = ARES;
2126                          ARES_INC = 0.0;
2127                          if ( (coval) && (x[l]<0) && (aTmp) )
2128                              MINDEC(ret_c,2);
2129                          if ( (!coval) && (x[l]>0) && (aTmp))
2130                              MINDEC(ret_c,2);
2131                          AARG_INC += x[l] * aTmp;
2132                      }
2133                  }
2134                  HOS_OV_INC(Targ,k)
2135            }
2136                break;
2137
2138                /*--------------------------------------------------------------------------*/
2139            case ceil_op:                                              /* ceil_op */
2140                res   = get_locint_r();
2141                arg   = get_locint_r();
2142                coval = get_val_r();
2143
2144                GET_TAYL(res,k,p)
2145
2146                coval = (coval != ceil(*rpp_T[arg]) );
2147
2148                ASSIGN_A(Ares, rpp_A[res])
2149
2150                FOR_0_LE_l_LT_p
2151                if (0 == ARES) {
2152                    HOV_INC(Aarg,  k1)
2153                    HOV_INC(Ares,  k1)
2154                } else {
2155                    ARES_INC = 0.0;
2156                    AARG_INC = 5.0;
2157                    FOR_0_LE_i_LT_k
2158                    { if ((coval) && (ARES))
2159                      MINDEC(ret_c,2);
2160                      ARES_INC = 0.0;
2161                    }
2162                    HOV_INC(Aarg, k)
2163                    }
2164                break;
2165
2166                /*--------------------------------------------------------------------------*/
2167            case floor_op:                                            /* floor_op */
2168                res   = get_locint_r();
2169                arg   = get_locint_r();
2170                coval = get_val_r();
2171
2172                GET_TAYL(res,k,p)
2173
2174                coval = ( coval != floor(*rpp_T[arg]) );
2175
2176                ASSIGN_A(Ares, rpp_A[res])
2177                ASSIGN_A(Aarg, rpp_A[arg])
2178
2179                FOR_0_LE_l_LT_p
2180                if (0 == ARES) {
2181                    HOV_INC(Aarg, k1)
2182                    HOV_INC(Ares, k1)
2183                } else {
2184                    ARES = 0.0;
2185                    AARG_INC = 5.0;
2186                    FOR_0_LE_i_LT_k
2187                    { if ( (coval) && (ARES) )
2188                      MINDEC(ret_c,2);
2189                      ARES_INC = 0.0;
2190                    }
2191                    HOV_INC(Aarg, k)
2192                    }
2193                break;
2194
2195
2196                /****************************************************************************/
2197                /*                                                             CONDITIONALS */
2198
2199                /*--------------------------------------------------------------------------*/
2200            case cond_assign:                                      /* cond_assign */
2201                res   = get_locint_r();
2202                arg2  = get_locint_r();
2203                arg1  = get_locint_r();
2204                arg   = get_locint_r();
2205                coval = get_val_r();
2206
2207                GET_TAYL(res,k,p)
2208
2209                ASSIGN_A(Aarg1, rpp_A[arg1])
2210                ASSIGN_A(Ares,  rpp_A[res])
2211                ASSIGN_A(Aarg2, rpp_A[arg2])
2212                ASSIGN_T(Targ,  rpp_T[arg])
2213
2214                /* olvo 980925 changed code a little bit */
2215                if (TARG > 0.0) {
2216                    if (res != arg1)
2217                        FOR_0_LE_l_LT_p
2218                        { if (0 == ARES) {
2219                          HOV_INC(Ares,  k1)
2220                              HOV_INC(Aarg1, k1)
2221                          } else {
2222                              if (coval <= 0.0)
2223                                  MINDEC(ret_c,2);
2224                              MAXDEC(AARG1,ARES);
2225                              ARES_INC = 0.0;
2226                              AARG1_INC_O;
2227                              FOR_0_LE_i_LT_k
2228                              { AARG1_INC += ARES;
2229                                ARES_INC = 0;
2230                              }
2231                          }
2232                    }
2233                    else
2234                        FOR_0_LE_l_LT_p {
2235                            if ((coval <= 0.0) && (ARES))
2236                            MINDEC(ret_c,2);
2237                            HOV_INC(Ares,  k1)
2238                        }
2239                    } else /* TARG <= 0.0 */
2240            {
2241                if (res != arg2)
2242                        FOR_0_LE_l_LT_p
2243                        { if (0 == ARES) {
2244                          HOV_INC(Ares,  k1)
2245                              HOV_INC(Aarg2, k1)
2246                          } else {
2247                              if (TARG == 0.0) /* we are at the tie */
2248                              { MINDEC(ret_c,0);
2249                                  AARG1 = 5.0;
2250                                  AARG2_INC = 5.0;
2251                              } else {
2252                                  if (coval <= 0.0)
2253                                      MINDEC(ret_c,2);
2254                                  MAXDEC(AARG2,ARES);
2255                                  AARG2_INC_O;
2256                              }
2257                              ARES_INC = 0.0;
2258
2259                              FOR_0_LE_i_LT_k
2260                              { AARG2_INC += ARES;
2261                                ARES_INC = 0;
2262                              }
2263                          }
2264                      HOV_INC(Aarg1, k1)
2265                    } else
2266                        FOR_0_LE_l_LT_p {
2267                            if (ARES) {
2268                            if (TARG == 0.0) /* we are at the tie */
2269                                { MINDEC(ret_c,0);
2270                                    AARG1 = 5.0;
2271                                    AARG2 = 5.0;
2272                                } else
2273                                    if (coval <= 0.0)
2274                                        MINDEC(ret_c,2);
2275                            }
2276                        HOV_INC(Ares,  k1)
2277                        HOV_INC(Aarg1, k1)
2278                        HOV_INC(Aarg2, k1)
2279                    }
2280                }
2281                break;
2282
2283            case cond_eq_assign:                                      /* cond_eq_assign */
2284                res   = get_locint_r();
2285                arg2  = get_locint_r();
2286                arg1  = get_locint_r();
2287                arg   = get_locint_r();
2288                coval = get_val_r();
2289
2290                GET_TAYL(res,k,p)
2291
2292                ASSIGN_A(Aarg1, rpp_A[arg1])
2293                ASSIGN_A(Ares,  rpp_A[res])
2294                ASSIGN_A(Aarg2, rpp_A[arg2])
2295                ASSIGN_T(Targ,  rpp_T[arg])
2296
2297                /* olvo 980925 changed code a little bit */
2298                if (TARG >= 0.0) {
2299                    if (res != arg1)
2300                        FOR_0_LE_l_LT_p
2301                        { if (0 == ARES) {
2302                          HOV_INC(Ares,  k1)
2303                              HOV_INC(Aarg1, k1)
2304                          } else {
2305                              if (coval < 0.0)
2306                                  MINDEC(ret_c,2);
2307                              MAXDEC(AARG1,ARES);
2308                              ARES_INC = 0.0;
2309                              AARG1_INC_O;
2310                              FOR_0_LE_i_LT_k
2311                              { AARG1_INC += ARES;
2312                                ARES_INC = 0;
2313                              }
2314                          }
2315                    }
2316                    else
2317                        FOR_0_LE_l_LT_p {
2318                            if ((coval < 0.0) && (ARES))
2319                            MINDEC(ret_c,2);
2320                            HOV_INC(Ares,  k1)
2321                        }
2322                    } else /* TARG < 0.0 */
2323            {
2324                if (res != arg2)
2325                        FOR_0_LE_l_LT_p
2326                        { if (0 == ARES) {
2327                          HOV_INC(Ares,  k1)
2328                              HOV_INC(Aarg2, k1)
2329                          } else {
2330                                  if (coval < 0.0)
2331                                      MINDEC(ret_c,2);
2332                                  MAXDEC(AARG2,ARES);
2333                                  AARG2_INC_O;
2334                              ARES_INC = 0.0;
2335
2336                              FOR_0_LE_i_LT_k
2337                              { AARG2_INC += ARES;
2338                                ARES_INC = 0;
2339                              }
2340                          }
2341                      HOV_INC(Aarg1, k1)
2342                    } else
2343                        FOR_0_LE_l_LT_p {
2344                            if (ARES) {
2345                                    if (coval < 0.0)
2346                                        MINDEC(ret_c,2);
2347                            }
2348                        HOV_INC(Ares,  k1)
2349                        HOV_INC(Aarg1, k1)
2350                        HOV_INC(Aarg2, k1)
2351                    }
2352                }
2353                break;
2354
2355                /*--------------------------------------------------------------------------*/
2356            case cond_assign_s:                                  /* cond_assign_s */
2357                res   = get_locint_r();
2358                arg1  = get_locint_r();
2359                arg   = get_locint_r();
2360                coval = get_val_r();
2361
2362                GET_TAYL(res,k,p)
2363
2364                ASSIGN_A(Aarg1, rpp_A[arg1])
2365                ASSIGN_A(Ares,  rpp_A[res])
2366                ASSIGN_T(Targ,  rpp_T[arg])
2367
2368                /* olvo 980925 changed code a little bit */
2369                if (TARG == 0.0) /* we are at the tie */
2370                { FOR_0_LE_l_LT_p
2371                    { if  (ARES)
2372                      AARG1 = 5.0;
2373                      HOV_INC(Aarg1, k1)
2374                      HOV_INC(Ares,  k1)
2375                    }
2376                    MINDEC(ret_c,0);
2377                } else
2378                    if (TARG > 0.0) {
2379                        if (res != arg1)
2380                            FOR_0_LE_l_LT_p
2381                            { if  (0 == ARES) {
2382                              HOV_INC(Ares,  k1)
2383                                  HOV_INC(Aarg1, k1)
2384                              } else {
2385                                  if (coval <= 0.0)
2386                                      MINDEC(ret_c,2);
2387                                  MAXDEC(AARG1,ARES);
2388                                  ARES_INC = 0.0;
2389                                  AARG1_INC_O;
2390                                  FOR_0_LE_i_LT_k
2391                                  { (AARG1_INC) += ARES;
2392                                    ARES_INC = 0;
2393                                  }
2394                              }
2395                        }
2396                        else
2397                            FOR_0_LE_l_LT_p {
2398                                if ((coval <= 0.0) && (ARES))
2399                                MINDEC(ret_c,2);
2400                                HOV_INC(Ares,  k1)
2401                            }
2402                        }
2403            break;
2404            case cond_eq_assign_s:                                  /* cond_eq_assign_s */
2405                res   = get_locint_r();
2406                arg1  = get_locint_r();
2407                arg   = get_locint_r();
2408                coval = get_val_r();
2409
2410                GET_TAYL(res,k,p)
2411
2412                ASSIGN_A(Aarg1, rpp_A[arg1])
2413                ASSIGN_A(Ares,  rpp_A[res])
2414                ASSIGN_T(Targ,  rpp_T[arg])
2415
2416                /* olvo 980925 changed code a little bit */
2417                    if (TARG >= 0.0) {
2418                        if (res != arg1)
2419                            FOR_0_LE_l_LT_p
2420                            { if  (0 == ARES) {
2421                              HOV_INC(Ares,  k1)
2422                                  HOV_INC(Aarg1, k1)
2423                              } else {
2424                                  if (coval < 0.0)
2425                                      MINDEC(ret_c,2);
2426                                  MAXDEC(AARG1,ARES);
2427                                  ARES_INC = 0.0;
2428                                  AARG1_INC_O;
2429                                  FOR_0_LE_i_LT_k
2430                                  { (AARG1_INC) += ARES;
2431                                    ARES_INC = 0;
2432                                  }
2433                              }
2434                        }
2435                        else
2436                            FOR_0_LE_l_LT_p {
2437                                if ((coval < 0.0) && (ARES))
2438                                MINDEC(ret_c,2);
2439                                HOV_INC(Ares,  k1)
2440                            }
2441                        }
2442            break;
2443                /*--------------------------------------------------------------------------*/
2444                /* NEW CONDITIONALS */
2445                /*--------------------------------------------------------------------------*/
2446#if defined(ADOLC_ADVANCED_BRANCHING)
2447            case neq_a_a:
2448            case eq_a_a:
2449            case le_a_a:
2450            case ge_a_a:
2451            case lt_a_a:
2452            case gt_a_a:
2453            case neq_a_p:
2454            case eq_a_p:
2455            case le_a_p:
2456            case ge_a_p:
2457            case lt_a_p:
2458            case gt_a_p:
2459                res = get_locint_r();
2460                arg1 = get_locint_r();
2461                arg = get_locint_r();
2462                coval = get_val_r();
2463                ASSIGN_A(Ares, rpp_A[res])
2464
2465                FOR_0_LE_l_LT_pk1
2466                ARES_INC = 0.0;
2467
2468                GET_TAYL(res,k,p)
2469                break;
2470#endif
2471
2472                /*--------------------------------------------------------------------------*/
2473            case subscript:
2474                coval = get_val_r();
2475                {
2476                    size_t idx, numval = (size_t)trunc(fabs(coval));
2477                    locint vectorloc;
2478                    res = get_locint_r();
2479                    vectorloc = get_locint_r();
2480                    arg = get_locint_r();
2481                    ASSIGN_T(Targ, rpp_T[arg])
2482                    idx = (size_t)trunc(fabs(TARG));
2483                    if (idx >= numval)
2484                        fprintf(DIAG_OUT, "ADOL-C warning: index out of bounds while subscripting n=%zu, idx=%zu\n", numval, idx);
2485                    arg1 = vectorloc+idx;
2486                    ASSIGN_A(Aarg1, rpp_A[arg1])
2487                    ASSIGN_A(Ares, rpp_A[res])
2488
2489                    FOR_0_LE_l_LT_p
2490                    if (0 == ARES) {
2491                        HOV_INC(Aarg1, k1)
2492                        HOV_INC(Ares, k1)
2493                    } else {
2494                        MAXDEC(AARG1,ARES);
2495                        AARG1_INC_O;
2496                        ARES_INC = 0.0;
2497                        FOR_0_LE_i_LT_k
2498                        {
2499                            AARG1_INC += ARES;
2500                            ARES_INC = 0.0;
2501                        }
2502                    }
2503                    GET_TAYL(res,k,p)
2504                }
2505                break;
2506
2507            case subscript_ref:
2508                coval = get_val_r();
2509                {
2510                    size_t idx, numval = (size_t)trunc(fabs(coval));
2511                    locint vectorloc;
2512                    res = get_locint_r();
2513                    vectorloc = get_locint_r();
2514                    arg = get_locint_r();
2515                    ASSIGN_T(Targ, rpp_T[arg])
2516                    ASSIGN_T(Tres, rpp_T[res])
2517                    idx = (size_t)trunc(fabs(TARG));
2518                    if (idx >= numval)
2519                        fprintf(DIAG_OUT, "ADOL-C warning: index out of bounds while subscripting (ref) n=%zu, idx=%zu\n", numval, idx);
2520                    arg1 = (size_t)trunc(fabs(TRES));
2521                    /*
2522                     * This is actually NOP
2523                     * basically all we need is that arg1 == vectorloc[idx]
2524                     * so doing a check here is probably good
2525                     */
2526                    if (arg1 != vectorloc+idx) {
2527                        fprintf(DIAG_OUT, "ADOL-C error: indexed active position does not match referenced position\nindexed = %zu, referenced = %d\n", vectorloc+idx, arg1);
2528                        adolc_exit(-2,"",__func__,__FILE__,__LINE__);
2529                    }
2530                    GET_TAYL(res,k,p)
2531                }
2532                break;
2533
2534            case ref_copyout:
2535                res = get_locint_r();
2536                arg1 = get_locint_r();
2537
2538                ASSIGN_T(Targ1, rpp_T[arg1])
2539                arg = (size_t)trunc(fabs(TARG1));
2540
2541                ASSIGN_A(Ares, rpp_A[res])
2542                ASSIGN_A(Aarg, rpp_A[arg])
2543
2544                FOR_0_LE_l_LT_p
2545                if (0 == ARES) {
2546                    HOV_INC(Aarg, k1)
2547                    HOV_INC(Ares, k1)
2548                } else {
2549                    MAXDEC(AARG,ARES);
2550                    AARG_INC_O;
2551                    ARES_INC = 0.0;
2552                    FOR_0_LE_i_LT_k
2553                    {
2554                        AARG_INC += ARES;
2555                        ARES_INC = 0.0;
2556                    }
2557                }
2558                GET_TAYL(res,k,p)
2559                break;
2560
2561            case ref_incr_a:                        /* Increment an adouble    incr_a */
2562            case ref_decr_a:                        /* Increment an adouble    decr_a */
2563                arg1   = get_locint_r();
2564
2565                ASSIGN_T(Targ1, rpp_T[arg1])
2566                res = (size_t)trunc(fabs(TARG1));
2567
2568                GET_TAYL(res,k,p)
2569                break;
2570
2571            case ref_assign_d:      /* assign an adouble variable a    assign_d */
2572                /* double value. (=) */
2573                coval = get_val_r();
2574                /* fallthrough */
2575            case ref_assign_d_zero: /* assign an adouble a        assign_d_zero */
2576            case ref_assign_d_one:  /* double value. (=)           assign_d_one */
2577                arg1   = get_locint_r();
2578
2579                ASSIGN_T(Targ1, rpp_T[arg1])
2580                res = (size_t)trunc(fabs(TARG1));
2581
2582                ASSIGN_A(Ares, rpp_A[res])
2583
2584                FOR_0_LE_l_LT_pk1
2585                ARES_INC = 0.0;
2586
2587                GET_TAYL(res,k,p)
2588                break;
2589
2590            case ref_assign_p:
2591                arg    = get_locint_r();
2592                arg1   = get_locint_r();
2593                coval = ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.paramstore[arg];
2594
2595                ASSIGN_T(Targ1, rpp_T[arg1])
2596                res = (size_t)trunc(fabs(TARG1));
2597
2598                ASSIGN_A(Ares, rpp_A[res])
2599
2600                FOR_0_LE_l_LT_pk1
2601                ARES_INC = 0.0;
2602
2603                GET_TAYL(res,k,p)
2604                break;
2605
2606            case ref_assign_a:     /* assign an adouble variable an    assign_a */
2607                /* adouble value. (=) */
2608                arg1 = get_locint_r();
2609                arg = get_locint_r();
2610
2611                ASSIGN_T(Targ1, rpp_T[arg1])
2612                res = (size_t)trunc(fabs(TARG1));
2613
2614                ASSIGN_A(Aarg, rpp_A[arg])
2615                ASSIGN_A(Ares, rpp_A[res])
2616
2617                FOR_0_LE_l_LT_p
2618                if  (0 == ARES) {
2619                    HOV_INC(Aarg, k1)
2620                    HOV_INC(Ares, k1)
2621                } else {
2622                    MAXDEC(AARG,ARES);
2623                    AARG_INC_O;
2624                    ARES_INC = 0.0;
2625                    FOR_0_LE_i_LT_k
2626                    { /* ! no tempory */
2627                        AARG_INC += ARES;
2628                        ARES_INC = 0.0;
2629                    }
2630                }
2631
2632                GET_TAYL(res,k,p)
2633                break;
2634
2635            case ref_assign_ind:       /* assign an adouble variable an    assign_ind */
2636                /* independent double value (<<=) */
2637                arg1 = get_locint_r();
2638                ASSIGN_T(Targ1, rpp_T[arg1])
2639                res = (size_t)trunc(fabs(TARG1));
2640                ASSIGN_A(Ares, rpp_A[res])
2641
2642                FOR_0_LE_l_LT_p
2643                {
2644#ifdef _HOV_
2645                    if (nonzero) /* ??? question: why here? */
2646                    nonzero[l][indexi] = (int)ARES;
2647#endif /* _HOV_ */
2648                    ARES_INC_O;
2649                    FOR_0_LE_i_LT_k
2650                        RESULTS(l,indexi,i) = ARES_INC;
2651                }
2652
2653                GET_TAYL(res,k,p)
2654                    indexi--;
2655                break;
2656
2657        case ref_eq_plus_d:            /* Add a floating point to an    eq_plus_d */
2658            /* adouble. (+=) */
2659                arg1   = get_locint_r();
2660                ASSIGN_T(Targ1, rpp_T[arg1])
2661                res = (size_t)trunc(fabs(TARG1));
2662                coval = get_val_r();
2663
2664                GET_TAYL(res,k,p)
2665                break;
2666
2667        case ref_eq_plus_p:            /* Add a floating point to an    eq_plus_d */
2668            /* adouble. (+=) */
2669                arg1   = get_locint_r();
2670                arg    = get_locint_r();
2671                ASSIGN_T(Targ1, rpp_T[arg1])
2672                res = (size_t)trunc(fabs(TARG1));
2673                coval = ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.paramstore[arg];
2674
2675                GET_TAYL(res,k,p)
2676                break;
2677
2678            case ref_eq_plus_a:             /* Add an adouble to another    eq_plus_a */
2679                /* adouble. (+=) */
2680                arg1 = get_locint_r();
2681                arg = get_locint_r();
2682
2683                ASSIGN_T(Targ1, rpp_T[arg1])
2684                res = (size_t)trunc(fabs(TARG1));
2685                ASSIGN_A(Ares, rpp_A[res])
2686                ASSIGN_A(Aarg, rpp_A[arg])
2687
2688                FOR_0_LE_l_LT_p
2689                if  (0 == ARES) {
2690                    HOV_INC(Ares, k1)
2691                    HOV_INC(Aarg, k1)
2692                } else {
2693                    MAXDEC(AARG,ARES);
2694                    AARG_INC_O;
2695                    ARES_INC_O;
2696                    FOR_0_LE_i_LT_k
2697                    AARG_INC += ARES_INC;
2698                }
2699
2700                GET_TAYL(res,k,p)
2701                break;
2702
2703            case ref_eq_min_d:       /* Subtract a floating point from an    eq_min_d */
2704                /* adouble. (-=) */
2705                arg1   = get_locint_r();
2706                ASSIGN_T(Targ1, rpp_T[arg1])
2707                res = (size_t)trunc(fabs(TARG1));
2708                coval = get_val_r();
2709
2710                GET_TAYL(res,k,p)
2711                break;
2712
2713           case ref_eq_min_p:            /* Add a floating point to an    eq_min_d */
2714            /* adouble. (-=) */
2715                arg1   = get_locint_r();
2716                arg    = get_locint_r();
2717                ASSIGN_T(Targ1, rpp_T[arg1])
2718                res = (size_t)trunc(fabs(TARG1));
2719                coval = ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.paramstore[arg];
2720
2721                GET_TAYL(res,k,p)
2722                break;
2723
2724            case ref_eq_min_a:        /* Subtract an adouble from another    eq_min_a */
2725                /* adouble. (-=) */
2726                arg1 = get_locint_r();
2727                arg = get_locint_r();
2728               
2729                ASSIGN_T(Targ1, rpp_T[arg1])
2730                res = (size_t)trunc(fabs(TARG1));
2731                ASSIGN_A(Ares, rpp_A[res])
2732                ASSIGN_A(Aarg, rpp_A[arg])
2733
2734                FOR_0_LE_l_LT_p
2735                if  (0==ARES) {
2736                    HOV_INC(Ares, k1)
2737                    HOV_INC(Aarg, k1)
2738                } else {
2739                    MAXDEC(AARG,ARES);
2740                    AARG_INC_O;
2741                    ARES_INC_O;
2742                    FOR_0_LE_i_LT_k
2743                    AARG_INC -= ARES_INC;
2744                }
2745
2746                GET_TAYL(res,k,p)
2747                break;
2748
2749            case ref_eq_mult_d:              /* Multiply an adouble by a    eq_mult_d */
2750                /* flaoting point. (*=) */
2751                arg1   = get_locint_r();
2752                ASSIGN_T(Targ1, rpp_T[arg1])
2753                res = (size_t)trunc(fabs(TARG1));
2754                coval = get_val_r();
2755
2756                ASSIGN_A(Ares, rpp_A[res])
2757
2758                FOR_0_LE_l_LT_p
2759                if ( 0 == ARES_INC )
2760                    HOV_INC(Ares, k)
2761                    else
2762                        FOR_0_LE_i_LT_k
2763                        ARES_INC *= coval;
2764
2765                GET_TAYL(res,k,p)
2766                break;
2767
2768            case ref_eq_mult_p:              /* Multiply an adouble by a    eq_mult_p */
2769                /* flaoting point. (*=) */
2770                arg1   = get_locint_r();
2771                arg    = get_locint_r();
2772                ASSIGN_T(Targ1, rpp_T[arg1])
2773                res = (size_t)trunc(fabs(TARG1));
2774                coval = ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.paramstore[arg];
2775
2776                ASSIGN_A(Ares, rpp_A[res])
2777
2778                FOR_0_LE_l_LT_p
2779                if ( 0 == ARES_INC )
2780                    HOV_INC(Ares, k)
2781                    else
2782                        FOR_0_LE_i_LT_k
2783                        ARES_INC *= coval;
2784
2785                GET_TAYL(res,k,p)
2786                break;
2787
2788        case ref_eq_mult_a:       /* Multiply one adouble by another    eq_mult_a */
2789                /* (*=) */
2790                arg1 = get_locint_r();
2791                arg = get_locint_r();
2792                ASSIGN_T(Targ1, rpp_T[arg1])
2793                res = (size_t)trunc(fabs(TARG1));
2794
2795                GET_TAYL(res,k,p)
2796
2797                ASSIGN_A(Ares, rpp_A[res])
2798                ASSIGN_A(Aarg, rpp_A[arg])
2799                ASSIGN_A(Aqo,  rp_Atemp)
2800                ASSIGN_T(Tres, rpp_T[res])
2801                ASSIGN_T(Targ, rpp_T[arg])
2802
2803                FOR_0_LE_l_LT_p {
2804                    if (0 == ARES) {
2805                    HOV_INC(Aarg, k1)
2806                        HOV_INC(Ares, k1)
2807                    } else {
2808                        MAXDEC(ARES,2.0);
2809                        MAXDEC(AARG,ARES);
2810                        AARG_INC_O;
2811                        ARES_INC_O;
2812                        conv(k,Ares,Targ,rp_Atemp);
2813                        if(arg != res) {
2814                            inconv(k,Ares,Tres,Aarg);
2815                            FOR_0_LE_i_LT_k
2816                            ARES_INC = AQO_INC;
2817                        } else
2818                            FOR_0_LE_i_LT_k
2819                            ARES_INC = 2.0 * AQO_INC;
2820                        HOV_INC(Aarg,k)
2821                        HOS_OV_INC(Tres,k)
2822                        HOS_OV_INC(Targ,k)
2823                        HOS_OV_ASSIGN_A(Aqo,  rp_Atemp)
2824                    }
2825            }
2826                break;
2827
2828            case vec_copy:
2829
2830                res = get_locint_r();
2831                size = get_locint_r();
2832                arg = get_locint_r();
2833
2834                for(qq=0;qq<size;qq++) {
2835
2836                ASSIGN_A(Aarg, rpp_A[arg+qq])
2837                ASSIGN_A(Ares, rpp_A[res+qq])
2838
2839                FOR_0_LE_l_LT_p
2840                if  (0 == ARES) {
2841                    HOV_INC(Aarg, k1)
2842                    HOV_INC(Ares, k1)
2843                } else {
2844                    MAXDEC(AARG,ARES);
2845                    AARG_INC_O;
2846                    ARES_INC = 0.0;
2847                    FOR_0_LE_i_LT_k
2848                    { /* ! no tempory */
2849                        AARG_INC += ARES;
2850                        ARES_INC = 0.0;
2851                    }
2852                }
2853
2854                GET_TAYL(res+qq,k,p)
2855
2856                }
2857
2858                break;
2859
2860        case vec_dot:
2861                res = get_locint_r();
2862                size = get_locint_r();
2863                arg2 = get_locint_r();
2864                arg1 = get_locint_r();
2865                for (qq=0;qq<size;qq++) {
2866                ASSIGN_A(Ares,  rpp_A[res])
2867                ASSIGN_A(Aarg2, rpp_A[arg2+qq])
2868                ASSIGN_A(Aarg1, rpp_A[arg1+qq])
2869                ASSIGN_T(Targ1, rpp_T[arg1+qq])
2870                ASSIGN_T(Targ2, rpp_T[arg2+qq])
2871                FOR_0_LE_l_LT_p {
2872                if (0 == ARES) {
2873                    HOV_INC(Aarg1, k1)
2874                    HOV_INC(Aarg2, k1)
2875                    HOV_INC(Ares,  k1)
2876                } else {
2877                    comp = (ARES > 2.0) ? ARES : 2.0 ;
2878                    ARES_INC = comp;
2879                    MAXDEC(AARG1,comp);
2880                    MAXDEC(AARG2,comp);
2881                    AARG1_INC_O;
2882                    AARG2_INC_O;
2883
2884                    inconv(k,Ares,Targ1,Aarg2);
2885                    inconv(k,Ares,Targ2,Aarg1);
2886
2887                    HOV_INC(Ares,  k)
2888                    HOV_INC(Aarg1, k)
2889                    HOV_INC(Aarg2, k)
2890                    HOS_OV_INC(Targ1, k)
2891                    HOS_OV_INC(Targ2, k)
2892                    HOS_OV_INC(Tres, k)
2893                }
2894                }
2895                }
2896                GET_TAYL(res,k,p)
2897                break;
2898
2899        case vec_axpy:
2900                res = get_locint_r();
2901                size = get_locint_r();
2902                arg2 = get_locint_r();
2903                arg1 = get_locint_r();
2904                arg = get_locint_r();
2905                for (qq=0;qq<size;qq++) {
2906                ASSIGN_A(Ares,  rpp_A[res+qq])
2907                ASSIGN_A(Aarg,  rpp_A[arg])
2908                ASSIGN_A(Aarg2, rpp_A[arg2+qq])
2909                ASSIGN_A(Aarg1, rpp_A[arg1+qq])
2910                ASSIGN_T(Targ,  rpp_T[arg])
2911                ASSIGN_T(Targ1, rpp_T[arg1+qq])
2912                if (0 == ARES) {
2913                    HOV_INC(Aarg, k1)
2914                    HOV_INC(Aarg1, k1)
2915                    HOV_INC(Aarg2, k1)
2916                    HOV_INC(Ares,  k1)
2917                } else {
2918                    comp = (ARES > 2.0) ? ARES : 2.0 ;
2919                    MAXDEC(AARG2,ARES);
2920                    ARES_INC = 0.0;
2921                    MAXDEC(AARG,comp);
2922                    MAXDEC(AARG1,comp);
2923                    AARG_INC_O;
2924                    AARG1_INC_O;
2925                    AARG2_INC_O;
2926                    copyAndZeroset(k,Ares,rp_Atemp);
2927                    inconv(k,rp_Atemp,Targ1,Aarg);
2928                    inconv(k,rp_Atemp,Targ,Aarg1);
2929                    FOR_0_LE_i_LT_k
2930                        AARG2_INC += *rp_Atemp++;
2931
2932                    HOV_INC(Ares,k)
2933                    HOV_INC(Aarg, k)
2934                    HOV_INC(Aarg1,k)
2935                    HOS_OV_INC(Targ,k)
2936                    HOS_OV_INC(Targ1,k)
2937                }
2938                GET_TAYL(res+qq,k,p)
2939                }
2940                break;               
2941
2942            case ref_cond_assign:                                      /* cond_assign */
2943            {   
2944                revreal *Tref;
2945                locint ref   = get_locint_r();
2946                arg2  = get_locint_r();
2947                arg1  = get_locint_r();
2948                arg   = get_locint_r();
2949                coval = get_val_r();
2950               
2951                ASSIGN_T(Tref, rpp_T[ref])
2952
2953#ifdef _HIGHER_ORDER_
2954#define TREF  *Tref
2955#else
2956#define TREF   rpp_T[ref]
2957#endif   
2958
2959                res = (size_t)trunc(fabs(TREF));
2960#undef TREF
2961                GET_TAYL(res,k,p)
2962
2963                ASSIGN_A(Aarg1, rpp_A[arg1])
2964                ASSIGN_A(Ares,  rpp_A[res])
2965                ASSIGN_A(Aarg2, rpp_A[arg2])
2966                ASSIGN_T(Targ,  rpp_T[arg])
2967
2968                /* olvo 980925 changed code a little bit */
2969                if (TARG > 0.0) {
2970                    if (res != arg1)
2971                        FOR_0_LE_l_LT_p
2972                        { if (0 == ARES) {
2973                          HOV_INC(Ares,  k1)
2974                              HOV_INC(Aarg1, k1)
2975                          } else {
2976                              if (coval <= 0.0)
2977                                  MINDEC(ret_c,2);
2978                              MAXDEC(AARG1,ARES);
2979                              ARES_INC = 0.0;
2980                              AARG1_INC_O;
2981                              FOR_0_LE_i_LT_k
2982                              { AARG1_INC += ARES;
2983                                ARES_INC = 0;
2984                              }
2985                          }
2986                    }
2987                    else
2988                        FOR_0_LE_l_LT_p {
2989                            if ((coval <= 0.0) && (ARES))
2990                            MINDEC(ret_c,2);
2991                            HOV_INC(Ares,  k1)
2992                        }
2993                    } else /* TARG <= 0.0 */
2994            {
2995                if (res != arg2)
2996                        FOR_0_LE_l_LT_p
2997                        { if (0 == ARES) {
2998                          HOV_INC(Ares,  k1)
2999                              HOV_INC(Aarg2, k1)
3000                          } else {
3001                              if (TARG == 0.0) /* we are at the tie */
3002                              { MINDEC(ret_c,0);
3003                                  AARG1 = 5.0;
3004                                  AARG2_INC = 5.0;
3005                              } else {
3006                                  if (coval <= 0.0)
3007                                      MINDEC(ret_c,2);
3008                                  MAXDEC(AARG2,ARES);
3009                                  AARG2_INC_O;
3010                              }
3011                              ARES_INC = 0.0;
3012
3013                              FOR_0_LE_i_LT_k
3014                              { AARG2_INC += ARES;
3015                                ARES_INC = 0;
3016                              }
3017                          }
3018                      HOV_INC(Aarg1, k1)
3019                    } else
3020                        FOR_0_LE_l_LT_p {
3021                            if (ARES) {
3022                            if (TARG == 0.0) /* we are at the tie */
3023                                { MINDEC(ret_c,0);
3024                                    AARG1 = 5.0;
3025                                    AARG2 = 5.0;
3026                                } else
3027                                    if (coval <= 0.0)
3028                                        MINDEC(ret_c,2);
3029                            }
3030                        HOV_INC(Ares,  k1)
3031                        HOV_INC(Aarg1, k1)
3032                        HOV_INC(Aarg2, k1)
3033                    }
3034                }
3035            }
3036                break;
3037            case ref_cond_eq_assign:                                      /* cond_eq_assign */
3038            {   
3039                revreal *Tref;
3040                locint ref   = get_locint_r();
3041                arg2  = get_locint_r();
3042                arg1  = get_locint_r();
3043                arg   = get_locint_r();
3044                coval = get_val_r();
3045               
3046                ASSIGN_T(Tref, rpp_T[ref])
3047
3048#ifdef _HIGHER_ORDER_
3049#define TREF  *Tref
3050#else
3051#define TREF   rpp_T[ref]
3052#endif   
3053
3054                res = (size_t)trunc(fabs(TREF));
3055#undef TREF
3056                GET_TAYL(res,k,p)
3057
3058                ASSIGN_A(Aarg1, rpp_A[arg1])
3059                ASSIGN_A(Ares,  rpp_A[res])
3060                ASSIGN_A(Aarg2, rpp_A[arg2])
3061                ASSIGN_T(Targ,  rpp_T[arg])
3062
3063                /* olvo 980925 changed code a little bit */
3064                if (TARG >= 0.0) {
3065                    if (res != arg1)
3066                        FOR_0_LE_l_LT_p
3067                        { if (0 == ARES) {
3068                          HOV_INC(Ares,  k1)
3069                              HOV_INC(Aarg1, k1)
3070                          } else {
3071                              if (coval < 0.0)
3072                                  MINDEC(ret_c,2);
3073                              MAXDEC(AARG1,ARES);
3074                              ARES_INC = 0.0;
3075                              AARG1_INC_O;
3076                              FOR_0_LE_i_LT_k
3077                              { AARG1_INC += ARES;
3078                                ARES_INC = 0;
3079                              }
3080                          }
3081                    }
3082                    else
3083                        FOR_0_LE_l_LT_p {
3084                            if ((coval < 0.0) && (ARES))
3085                            MINDEC(ret_c,2);
3086                            HOV_INC(Ares,  k1)
3087                        }
3088                    } else /* TARG < 0.0 */
3089            {
3090                if (res != arg2)
3091                        FOR_0_LE_l_LT_p
3092                        { if (0 == ARES) {
3093                          HOV_INC(Ares,  k1)
3094                              HOV_INC(Aarg2, k1)
3095                          } else {
3096                                  if (coval < 0.0)
3097                                      MINDEC(ret_c,2);
3098                                  MAXDEC(AARG2,ARES);
3099                                  AARG2_INC_O;
3100                              ARES_INC = 0.0;
3101
3102                              FOR_0_LE_i_LT_k
3103                              { AARG2_INC += ARES;
3104                                ARES_INC = 0;
3105                              }
3106                          }
3107                      HOV_INC(Aarg1, k1)
3108                    } else
3109                        FOR_0_LE_l_LT_p {
3110                            if (ARES) {
3111                                    if (coval < 0.0)
3112                                        MINDEC(ret_c,2);
3113                            }
3114                        HOV_INC(Ares,  k1)
3115                        HOV_INC(Aarg1, k1)
3116                        HOV_INC(Aarg2, k1)
3117                    }
3118                }
3119            }
3120                break;
3121
3122            case ref_cond_assign_s:                                  /* cond_assign_s */
3123                arg2   = get_locint_r();
3124                arg1  = get_locint_r();
3125                arg   = get_locint_r();
3126                coval = get_val_r();
3127               
3128                ASSIGN_T(Targ2, rpp_T[arg2])
3129                res = (size_t)trunc(fabs(TARG2));
3130
3131                GET_TAYL(res,k,p)
3132
3133                ASSIGN_A(Aarg1, rpp_A[arg1])
3134                ASSIGN_A(Ares,  rpp_A[res])
3135                ASSIGN_T(Targ,  rpp_T[arg])
3136
3137                /* olvo 980925 changed code a little bit */
3138                if (TARG == 0.0) /* we are at the tie */
3139                { FOR_0_LE_l_LT_p
3140                    { if  (ARES)
3141                      AARG1 = 5.0;
3142                      HOV_INC(Aarg1, k1)
3143                      HOV_INC(Ares,  k1)
3144                    }
3145                    MINDEC(ret_c,0);
3146                } else
3147                    if (TARG > 0.0) {
3148                        if (res != arg1)
3149                            FOR_0_LE_l_LT_p
3150                            { if  (0 == ARES) {
3151                              HOV_INC(Ares,  k1)
3152                                  HOV_INC(Aarg1, k1)
3153                              } else {
3154                                  if (coval <= 0.0)
3155                                      MINDEC(ret_c,2);
3156                                  MAXDEC(AARG1,ARES);
3157                                  ARES_INC = 0.0;
3158                                  AARG1_INC_O;
3159                                  FOR_0_LE_i_LT_k
3160                                  { (AARG1_INC) += ARES;
3161                                    ARES_INC = 0;
3162                                  }
3163                              }
3164                        }
3165                        else
3166                            FOR_0_LE_l_LT_p {
3167                                if ((coval <= 0.0) && (ARES))
3168                                MINDEC(ret_c,2);
3169                                HOV_INC(Ares,  k1)
3170                            }
3171                        }
3172            break;
3173
3174            case ref_cond_eq_assign_s:                                  /* cond_eq_assign_s */
3175                arg2   = get_locint_r();
3176                arg1  = get_locint_r();
3177                arg   = get_locint_r();
3178                coval = get_val_r();
3179               
3180                ASSIGN_T(Targ2, rpp_T[arg2])
3181                res = (size_t)trunc(fabs(TARG2));
3182
3183                GET_TAYL(res,k,p)
3184
3185                ASSIGN_A(Aarg1, rpp_A[arg1])
3186                ASSIGN_A(Ares,  rpp_A[res])
3187                ASSIGN_T(Targ,  rpp_T[arg])
3188
3189                /* olvo 980925 changed code a little bit */
3190                    if (TARG >= 0.0) {
3191                        if (res != arg1)
3192                            FOR_0_LE_l_LT_p
3193                            { if  (0 == ARES) {
3194                              HOV_INC(Ares,  k1)
3195                                  HOV_INC(Aarg1, k1)
3196                              } else {
3197                                  if (coval < 0.0)
3198                                      MINDEC(ret_c,2);
3199                                  MAXDEC(AARG1,ARES);
3200                                  ARES_INC = 0.0;
3201                                  AARG1_INC_O;
3202                                  FOR_0_LE_i_LT_k
3203                                  { (AARG1_INC) += ARES;
3204                                    ARES_INC = 0;
3205                                  }
3206                              }
3207                        }
3208                        else
3209                            FOR_0_LE_l_LT_p {
3210                                if ((coval < 0.0) && (ARES))
3211                                MINDEC(ret_c,2);
3212                                HOV_INC(Ares,  k1)
3213                            }
3214                        }
3215            break;
3216
3217                /****************************************************************************/
3218                /*                                                          REMAINING STUFF */
3219
3220                /*--------------------------------------------------------------------------*/
3221            case take_stock_op:                                  /* take_stock_op */
3222                res = get_locint_r();
3223                size = get_locint_r();
3224                get_val_v_r(size);
3225
3226                res += size;
3227                for (ls=size;ls>0;ls--) {
3228                    res--;
3229
3230                    ASSIGN_A( Ares, rpp_A[res])
3231
3232                    FOR_0_LE_l_LT_pk1
3233                    ARES_INC = 0.0;
3234                }
3235                break;
3236
3237                /*--------------------------------------------------------------------------*/
3238            case death_not:                                          /* death_not */
3239                arg2 = get_locint_r();
3240                arg1 = get_locint_r();
3241
3242                for (j=arg1;j<=arg2;j++) {
3243                    ASSIGN_A(Aarg1, rpp_A[j])
3244
3245                    FOR_0_LE_l_LT_p
3246                    for (i=0; i<k1; i++)
3247                        AARG1_INC = 0.0;
3248                }
3249               
3250                for (j=arg1;j<=arg2;j++)
3251                    GET_TAYL(j,k,p)
3252
3253                break;
3254
3255                /*--------------------------------------------------------------------------*/
3256            default:                                                   /* default */
3257                /*             Die here, we screwed up     */
3258
3259                fprintf(DIAG_OUT,"ADOL-C fatal error in " GENERATED_FILENAME " ("
3260                        __FILE__
3261                        ") : no such operation %d\n", operation);
3262                adolc_exit(-1,"",__func__,__FILE__,__LINE__);
3263                break;
3264        } /* endswitch */
3265
3266        /* Get the next operation */   
3267        operation=get_op_r();
3268#if defined(ADOLC_DEBUG)
3269        ++countPerOperation[operation];
3270#endif /* ADOLC_DEBUG */
3271    }
3272
3273#if defined(ADOLC_DEBUG)
3274    printf("\nTape contains:\n");
3275    for (v = 0; v < 256; ++v)
3276        if (countPerOperation[v] > 0)
3277            printf("operation %3d: %6d time(s) - %6d taylors read (%10.2f per operation)\n", v, countPerOperation[v], taylorPerOperation[v], (double)taylorPerOperation[v] / (double)countPerOperation[v]);
3278    printf("\n");
3279#endif /* ADOLC_DEBUG */
3280
3281    /* clean up */
3282    free((char*)*rpp_T);
3283    free((char*) rpp_T);
3284    free(*rpp_A);
3285    free(rpp_A);
3286    free(rp_Ttemp);
3287    free(rp_Ttemp2);
3288    free(rp_Atemp);
3289    free(rp_Atemp2);
3290
3291    free((char*) jj);
3292    free((char*) x);
3293
3294    ADOLC_CURRENT_TAPE_INFOS.workMode = ADOLC_NO_MODE;
3295    end_sweep();
3296
3297    return ret_c;
3298}
3299
3300/****************************************************************************/
3301/*                                                               THAT'S ALL */
3302
3303END_C_DECLS
Note: See TracBrowser for help on using the repository browser.