source: branches/MPI/ADOL-C/src/fo_rev.c @ 424

Last change on this file since 424 was 424, checked in by kulshres, 9 years ago

update with 'mpi' branch from 'gitclone'

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

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