source: stable/2.4/ADOL-C/src/tape_handling.cpp @ 397

Last change on this file since 397 was 353, checked in by kulshres, 7 years ago

check for zero blocks separately

iter=indexFree.erase(iter) advances iter to the next element, then the
for loop advances it again. So we potentially miss a block, which could
either be a tail block or a zero block. Just stop as soon as we find tail
and we remove the zero blocks after allocating the new tail block, if any
are remaining.

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

  • Property svn:keywords set to Author Date Id Revision
File size: 49.2 KB
Line 
1/*----------------------------------------------------------------------------
2 ADOL-C -- Automatic Differentiation by Overloading in C++
3 File:     tape_handling.cpp
4 Revision: $Id: tape_handling.cpp 353 2012-09-24 15:11:31Z kulshres $
5 Contents: management of tape infos
6
7 Copyright (c) Andreas Kowarz, Andrea Walther, Kshitij Kulshreshtha,
8               Benjamin Letschert
9
10 This file is part of ADOL-C. This software is provided as open source.
11 Any use, reproduction, or distribution of the software constitutes
12 recipient's acceptance of the terms of the accompanying license file.
13 
14---------------------------------------------------------------------------*/
15#include "taping_p.h"
16#include "checkpointing_p.h"
17#include <adolc/revolve.h>
18
19#include <cassert>
20#include <limits>
21#include <iostream>
22#include <string.h>
23#ifdef HAVE_UNISTD_H
24#include <unistd.h>
25#endif
26#include <vector>
27#include <stack>
28#include <errno.h>
29
30using namespace std;
31
32#ifdef SPARSE
33BEGIN_C_DECLS
34extern void freeSparseJacInfos(double *y, double **B, unsigned int **JP, void *g, 
35                               void *jr1d, int seed_rows, int seed_clms, int depen);
36extern void freeSparseHessInfos(double **Hcomp, double ***Xppp, double ***Yppp, double ***Zppp, 
37                                double **Upp, unsigned int **HP,
38                                void *g, void *hr, int p, int indep);
39END_C_DECLS
40#endif
41
42GlobalTapeVarsCL::GlobalTapeVarsCL() {
43  store = NULL;
44  storeSize = 0;
45  numLives = 0;
46  storeManagerPtr = new StoreManagerLocintBlock(store, storeSize, numLives);
47}
48
49GlobalTapeVarsCL::~GlobalTapeVarsCL() {
50  if (storeManagerPtr) {
51    delete storeManagerPtr;
52    storeManagerPtr = NULL;
53  }
54}
55
56const GlobalTapeVarsCL& GlobalTapeVarsCL::operator=(const GlobalTapeVarsCL& gtv) {
57    storeSize = gtv.storeSize;
58    numLives = gtv.numLives;
59    maxLoc = gtv.maxLoc;
60    operationBufferSize = gtv.operationBufferSize;
61    locationBufferSize = gtv.locationBufferSize;
62    valueBufferSize = gtv.valueBufferSize;
63    taylorBufferSize = gtv.taylorBufferSize;
64    maxNumberTaylorBuffers = gtv.maxNumberTaylorBuffers;
65    inParallelRegion = gtv.inParallelRegion;
66    newTape = gtv.newTape;
67    branchSwitchWarning = gtv.branchSwitchWarning;
68    currentTapeInfosPtr = gtv.currentTapeInfosPtr;
69    store = new double[storeSize];
70    memcpy(store, gtv.store, storeSize*sizeof(double));
71    storeManagerPtr = new
72        StoreManagerLocintBlock(
73            dynamic_cast<StoreManagerLocintBlock*>(gtv.storeManagerPtr),
74            store, storeSize, numLives);
75    return *this;
76}
77
78StoreManagerLocint::StoreManagerLocint(double * &storePtr, size_t &size, size_t &numlives) : 
79    storePtr(storePtr),
80    indexFree(0),
81    head(0),
82    maxsize(size), currentfill(numlives)
83{
84#ifdef ADOLC_DEBUG
85    std::cerr << "StoreManagerInteger::StoreManagerInteger()\n";
86#endif
87}
88
89StoreManagerLocint::~StoreManagerLocint() 
90{
91#ifdef ADOLC_DEBUG
92    std::cerr << "StoreManagerInteger::~StoreManagerInteger()\n";
93#endif
94    if (storePtr) {
95        delete[] storePtr;
96        storePtr = 0;
97    }
98    if (indexFree) {
99        delete[] indexFree;
100        indexFree = 0;
101    }
102    maxsize = 0;
103    currentfill = 0;
104    head = 0;
105}
106
107StoreManagerLocint::StoreManagerLocint(const StoreManagerLocint *const stm,
108                                       double * &storePtr, size_t &size, size_t &numlives) : 
109    storePtr(storePtr),
110    maxsize(size), currentfill(numlives)
111{
112#ifdef ADOLC_DEBUG
113    std::cerr << "StoreManagerInteger::StoreManagerInteger()\n";
114#endif
115    head = stm->head;
116    indexFree = new locint[maxsize];
117    for (int i = 0; i < maxsize; i++)
118        indexFree[i] = stm->indexFree[i];
119}
120
121locint StoreManagerLocint::next_loc() {
122    if (head == 0) {
123      grow();
124    }
125    assert(head);
126    locint const result = head;
127    head = indexFree[head];
128    ++currentfill;
129#ifdef ADOLC_DEBUG
130    std::cerr << "next_loc: " << result << " fill: " << size() << "max: " << maxSize() << endl;
131#endif
132    return result;
133}
134
135void StoreManagerLocint::free_loc(locint loc) {
136    assert(0 < loc && loc < maxsize);
137    indexFree[loc] = head;
138    head = loc;
139    --currentfill;
140#ifdef ADOLC_DEBUG
141    std::cerr << "free_loc: " << loc << " fill: " << size() << "max: " << maxSize() << endl;
142#endif
143}
144
145void StoreManagerLocint::grow() {
146    if (maxsize == 0) maxsize += initialSize;
147    size_t const oldMaxsize = maxsize;
148    maxsize *= 2;
149
150    if (maxsize > std::numeric_limits<locint>::max()) {
151      // encapsulate this error message
152      fprintf(DIAG_OUT,"\nADOL-C error:\n");
153      fprintf(DIAG_OUT,"maximal number (%d) of live active variables exceeded\n\n", 
154              std::numeric_limits<locint>::max());
155      exit(-3);
156    }
157
158#ifdef ADOLC_DEBUG
159    std::cerr << "StoreManagerInteger::grow(): increase size from " << oldMaxsize
160         << " to " << maxsize << " entries (currently " << size() << " entries used)\n";
161    assert(oldMaxsize == initialSize or size() == oldMaxsize);
162#endif
163
164    double *const oldStore = storePtr;
165    locint *const oldIndex = indexFree;
166
167#if defined(ADOLC_DEBUG)
168    std::cerr << "StoreManagerInteger::grow(): allocate " << maxsize * sizeof(double) << " B doubles " 
169         << "and " << maxsize * sizeof(locint) << " B locints\n";
170#endif
171    storePtr = new double[maxsize];
172    indexFree = new locint[maxsize];
173    // we use index 0 as end-of-list marker
174    size_t i = 1;
175    storePtr[0] =  std::numeric_limits<double>::quiet_NaN();
176
177    if (oldMaxsize != initialSize) { // not the first time
178#if defined(ADOLC_DEBUG)
179      std::cerr << "StoreManagerInteger::grow(): copy values\n";
180#endif
181      for (size_t j = i; j < oldMaxsize; ++j) {
182        indexFree[j] = oldIndex[j];
183      }
184      for (size_t j = i; j < oldMaxsize; ++j) {
185        storePtr[j] = oldStore[j];
186      }
187
188      // reset i to start of new slots (upper half)
189      i = oldMaxsize;
190
191#if defined(ADOLC_DEBUG)
192      std::cerr << "StoreManagerInteger::grow(): free " << oldMaxsize * sizeof(double)
193                << " + " << oldMaxsize * sizeof(locint) << " B\n";
194#endif
195      delete [] oldStore;
196      delete [] oldIndex;
197    }
198
199    head = i;
200    // create initial linked list for new slots
201    for ( ; i < maxsize-1; ++i) {
202      indexFree[i] = i + 1;
203    }
204    indexFree[i] = 0; // end marker
205    assert(i == maxsize-1);
206}
207
208
209/****************************************************************************/
210/* Returns the next free location in "adouble" memory.                      */
211/****************************************************************************/
212locint next_loc() {
213  ADOLC_OPENMP_THREAD_NUMBER;
214  ADOLC_OPENMP_GET_THREAD_NUMBER;
215  return ADOLC_GLOBAL_TAPE_VARS.storeManagerPtr->next_loc();
216}
217
218/****************************************************************************/
219/* frees the specified location in "adouble" memory                         */
220/****************************************************************************/
221void free_loc(locint loc) {
222  ADOLC_OPENMP_THREAD_NUMBER;
223  ADOLC_OPENMP_GET_THREAD_NUMBER;
224  ADOLC_GLOBAL_TAPE_VARS.storeManagerPtr->free_loc(loc);
225}
226
227/* vector of tape infos for all tapes in use */
228vector<TapeInfos *> ADOLC_TAPE_INFOS_BUFFER_DECL;
229
230/* stack of pointers to tape infos
231 * represents the order of tape usage when doing nested taping */
232stack<TapeInfos *> ADOLC_TAPE_STACK_DECL;
233
234/* the main tape info buffer and its fallback */
235TapeInfos ADOLC_CURRENT_TAPE_INFOS_DECL;
236TapeInfos ADOLC_CURRENT_TAPE_INFOS_FALLBACK_DECL;
237
238/* global tapeing variables */
239GlobalTapeVars ADOLC_GLOBAL_TAPE_VARS_DECL;
240
241#if defined(_OPENMP)
242static vector<TapeInfos *> *tapeInfosBuffer_s;
243static stack<TapeInfos *>  *tapeStack_s;
244static TapeInfos           *currentTapeInfos_s;
245static TapeInfos           *currentTapeInfos_fallBack_s;
246static GlobalTapeVars      *globalTapeVars_s;
247static ADOLC_BUFFER_TYPE   *ADOLC_extDiffFctsBuffer_s;
248static stack<StackElement> *ADOLC_checkpointsStack_s;
249static revolve_nums        *revolve_numbers_s;
250
251static vector<TapeInfos *> *tapeInfosBuffer_p;
252static stack<TapeInfos *>  *tapeStack_p;
253static TapeInfos           *currentTapeInfos_p;
254static TapeInfos           *currentTapeInfos_fallBack_p;
255static GlobalTapeVars      *globalTapeVars_p;
256static ADOLC_BUFFER_TYPE   *ADOLC_extDiffFctsBuffer_p;
257static stack<StackElement> *ADOLC_checkpointsStack_p;
258static revolve_nums        *revolve_numbers_p;
259#endif
260
261/*--------------------------------------------------------------------------*/
262/* This function sets the flag "newTape" if either a taylor buffer has been */
263/* created or a taping process has been performed. Calling the function is  */
264/* also useful to "convince" the linker of including the cleaner part into  */
265/* the binary when linking statically!                                      */
266/*--------------------------------------------------------------------------*/
267void markNewTape() {
268    ADOLC_OPENMP_THREAD_NUMBER;
269    ADOLC_OPENMP_GET_THREAD_NUMBER;
270    ADOLC_GLOBAL_TAPE_VARS.newTape = 1;
271}
272
273/* inits the struct for the new tape */
274void initTapeInfos(TapeInfos *newTapeInfos) {
275    char *ptr;
276
277    ptr = (char *)newTapeInfos;
278    for (unsigned int i = 0; i < sizeof(TapeInfos) -
279            sizeof(PersistantTapeInfos); ++i) ptr[i] = 0;
280}
281
282/* as above but keep allocated buffers if possible */
283void initTapeInfos_keep(TapeInfos *newTapeInfos) {
284    unsigned char *opBuffer = newTapeInfos->opBuffer;
285    locint *locBuffer = newTapeInfos->locBuffer;
286    double *valBuffer = newTapeInfos->valBuffer;
287    revreal *tayBuffer = newTapeInfos->tayBuffer;
288    FILE *tay_file = newTapeInfos->tay_file;
289
290    initTapeInfos(newTapeInfos);
291
292    newTapeInfos->opBuffer = opBuffer;
293    newTapeInfos->locBuffer = locBuffer;
294    newTapeInfos->valBuffer = valBuffer;
295    newTapeInfos->tayBuffer = tayBuffer;
296    newTapeInfos->tay_file = tay_file;
297}
298
299/* inits a new tape and updates the tape stack (called from start_trace)
300 * - returns 0 without error
301 * - returns 1 if tapeID was already/still in use */
302int initNewTape(short tapeID) {
303    TapeInfos *newTapeInfos = NULL;
304    bool newTI = false;
305    int retval = 0;
306
307    ADOLC_OPENMP_THREAD_NUMBER;
308    ADOLC_OPENMP_GET_THREAD_NUMBER;
309
310    /* check if tape is in use */
311    vector<TapeInfos *>::iterator tiIter;
312    if (!ADOLC_TAPE_INFOS_BUFFER.empty()) {
313        for (tiIter=ADOLC_TAPE_INFOS_BUFFER.begin();
314                tiIter!=ADOLC_TAPE_INFOS_BUFFER.end();
315                ++tiIter) {
316            if ((*tiIter)->tapeID==tapeID) {
317                newTapeInfos=*tiIter;
318                if ((*tiIter)->inUse != 0) {
319                    if ((*tiIter)->tapingComplete == 0)
320                        fail(ADOLC_TAPING_TAPE_STILL_IN_USE);
321                    if ( (*tiIter)->stats[OP_FILE_ACCESS]  == 0 &&
322                            (*tiIter)->stats[LOC_FILE_ACCESS] == 0 &&
323                            (*tiIter)->stats[VAL_FILE_ACCESS] == 0  ) {
324#              if defined(ADOLC_DEBUG)
325                        fprintf(DIAG_OUT, "\nADOL-C warning: Tape %d existed in main memory"
326                                " only and gets overwritten!\n\n", tapeID);
327#              endif
328                        /* free associated resources */
329                        retval = 1;
330                    }
331                }
332                if ((*tiIter)->tay_file != NULL)
333                    rewind((*tiIter)->tay_file);
334                initTapeInfos_keep(*tiIter);
335                (*tiIter)->tapeID = tapeID;
336#ifdef SPARSE
337                freeSparseJacInfos(newTapeInfos->pTapeInfos.sJinfos.y,
338                                   newTapeInfos->pTapeInfos.sJinfos.B,
339                                   newTapeInfos->pTapeInfos.sJinfos.JP,
340                                   newTapeInfos->pTapeInfos.sJinfos.g,
341                                   newTapeInfos->pTapeInfos.sJinfos.jr1d,
342                                   newTapeInfos->pTapeInfos.sJinfos.seed_rows,
343                                   newTapeInfos->pTapeInfos.sJinfos.seed_clms,
344                                   newTapeInfos->pTapeInfos.sJinfos.depen);
345                freeSparseHessInfos(newTapeInfos->pTapeInfos.sHinfos.Hcomp, 
346                                    newTapeInfos->pTapeInfos.sHinfos.Xppp, 
347                                    newTapeInfos->pTapeInfos.sHinfos.Yppp, 
348                                    newTapeInfos->pTapeInfos.sHinfos.Zppp, 
349                                    newTapeInfos->pTapeInfos.sHinfos.Upp, 
350                                    newTapeInfos->pTapeInfos.sHinfos.HP,
351                                    newTapeInfos->pTapeInfos.sHinfos.g, 
352                                    newTapeInfos->pTapeInfos.sHinfos.hr, 
353                                    newTapeInfos->pTapeInfos.sHinfos.p, 
354                                    newTapeInfos->pTapeInfos.sHinfos.indep);   
355                newTapeInfos->pTapeInfos.inJacSparseUse=0;
356                newTapeInfos->pTapeInfos.inHessSparseUse=0;
357                newTapeInfos->pTapeInfos.sJinfos.B=NULL;
358                newTapeInfos->pTapeInfos.sJinfos.y=NULL;
359                newTapeInfos->pTapeInfos.sJinfos.g=NULL;
360                newTapeInfos->pTapeInfos.sJinfos.jr1d=NULL;
361                newTapeInfos->pTapeInfos.sJinfos.Seed=NULL;
362                newTapeInfos->pTapeInfos.sJinfos.JP=NULL;
363                newTapeInfos->pTapeInfos.sJinfos.depen=0;
364                newTapeInfos->pTapeInfos.sJinfos.nnz_in=0;
365                newTapeInfos->pTapeInfos.sJinfos.seed_rows=0;
366                newTapeInfos->pTapeInfos.sJinfos.seed_clms=0;
367                newTapeInfos->pTapeInfos.sHinfos.Zppp=NULL;
368                newTapeInfos->pTapeInfos.sHinfos.Yppp=NULL;
369                newTapeInfos->pTapeInfos.sHinfos.Xppp=NULL;
370                newTapeInfos->pTapeInfos.sHinfos.Upp=NULL;
371                newTapeInfos->pTapeInfos.sHinfos.Hcomp=NULL;
372                newTapeInfos->pTapeInfos.sHinfos.HP=NULL;
373                newTapeInfos->pTapeInfos.sHinfos.g=NULL;
374                newTapeInfos->pTapeInfos.sHinfos.hr=NULL;
375                newTapeInfos->pTapeInfos.sHinfos.nnz_in=0;
376                newTapeInfos->pTapeInfos.sHinfos.indep=0;
377                newTapeInfos->pTapeInfos.sHinfos.p=0;
378#endif
379                break;
380            }
381        }
382    }
383
384    /* create new info struct and initialize it */
385    if (newTapeInfos == NULL) {
386        newTapeInfos = new TapeInfos(tapeID);
387        newTI = true;
388    }
389    newTapeInfos->traceFlag=1;
390    newTapeInfos->inUse=1;
391#ifdef SPARSE
392    newTapeInfos->pTapeInfos.inJacSparseUse=0;
393    newTapeInfos->pTapeInfos.inHessSparseUse=0;
394    newTapeInfos->pTapeInfos.sJinfos.B=NULL;
395    newTapeInfos->pTapeInfos.sJinfos.y=NULL;
396    newTapeInfos->pTapeInfos.sJinfos.g=NULL;
397    newTapeInfos->pTapeInfos.sJinfos.jr1d=NULL;
398    newTapeInfos->pTapeInfos.sJinfos.Seed=NULL;
399    newTapeInfos->pTapeInfos.sJinfos.JP=NULL;
400    newTapeInfos->pTapeInfos.sJinfos.depen=0;
401    newTapeInfos->pTapeInfos.sJinfos.nnz_in=0;
402    newTapeInfos->pTapeInfos.sJinfos.seed_rows=0;
403    newTapeInfos->pTapeInfos.sJinfos.seed_clms=0;
404    newTapeInfos->pTapeInfos.sHinfos.Zppp=NULL;
405    newTapeInfos->pTapeInfos.sHinfos.Yppp=NULL;
406    newTapeInfos->pTapeInfos.sHinfos.Xppp=NULL;
407    newTapeInfos->pTapeInfos.sHinfos.Upp=NULL;
408    newTapeInfos->pTapeInfos.sHinfos.Hcomp=NULL;
409    newTapeInfos->pTapeInfos.sHinfos.HP=NULL;
410    newTapeInfos->pTapeInfos.sHinfos.g=NULL;
411    newTapeInfos->pTapeInfos.sHinfos.hr=NULL;
412    newTapeInfos->pTapeInfos.sHinfos.nnz_in=0;
413    newTapeInfos->pTapeInfos.sHinfos.indep=0;
414    newTapeInfos->pTapeInfos.sHinfos.p=0;
415#endif
416
417    newTapeInfos->stats[OP_BUFFER_SIZE] =
418        ADOLC_GLOBAL_TAPE_VARS.operationBufferSize;
419    newTapeInfos->stats[LOC_BUFFER_SIZE] =
420        ADOLC_GLOBAL_TAPE_VARS.locationBufferSize;
421    newTapeInfos->stats[VAL_BUFFER_SIZE] =
422        ADOLC_GLOBAL_TAPE_VARS.valueBufferSize;
423    newTapeInfos->stats[TAY_BUFFER_SIZE] =
424        ADOLC_GLOBAL_TAPE_VARS.taylorBufferSize;
425
426    /* update tapeStack and save tapeInfos */
427    if (ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr != NULL) {
428        memcpy(ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr,
429                &ADOLC_CURRENT_TAPE_INFOS, sizeof(TapeInfos));
430        ADOLC_TAPE_STACK.push(ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr);
431    } else {
432        memcpy(&ADOLC_CURRENT_TAPE_INFOS_FALLBACK,
433                &ADOLC_CURRENT_TAPE_INFOS, sizeof(TapeInfos));
434        ADOLC_TAPE_STACK.push(&ADOLC_CURRENT_TAPE_INFOS_FALLBACK);
435    }
436    if (newTI) ADOLC_TAPE_INFOS_BUFFER.push_back(newTapeInfos);
437
438    /* set the new tape infos as current */
439    memcpy(&ADOLC_CURRENT_TAPE_INFOS, newTapeInfos, sizeof(TapeInfos));
440    ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr = newTapeInfos;
441
442    return retval;
443}
444
445/* opens an existing tape or creates a new handle for a tape on hard disk
446 * - called from init_for_sweep and init_rev_sweep */
447void openTape(short tapeID, char mode) {
448    TapeInfos *tempTapeInfos=NULL;
449
450    ADOLC_OPENMP_THREAD_NUMBER;
451    ADOLC_OPENMP_GET_THREAD_NUMBER;
452
453    /* check if tape information exist in memory */
454    vector<TapeInfos *>::iterator tiIter;
455    if (!ADOLC_TAPE_INFOS_BUFFER.empty()) {
456        for (tiIter=ADOLC_TAPE_INFOS_BUFFER.begin();
457                tiIter!=ADOLC_TAPE_INFOS_BUFFER.end();
458                ++tiIter) {
459            if ((*tiIter)->tapeID == tapeID) {
460                /* tape has been used before (in the current program) */
461                if ((*tiIter)->inUse == 0) {
462                    /* forward sweep */
463                    if ((*tiIter)->tay_file != NULL)
464                        rewind((*tiIter)->tay_file);
465                    initTapeInfos_keep(*tiIter);
466                    (*tiIter)->traceFlag=1;
467                    (*tiIter)->tapeID = tapeID;
468                    (*tiIter)->tapingComplete = 1;
469                    (*tiIter)->inUse = 1;
470                    read_tape_stats(*tiIter);
471               }
472                if (ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr != NULL) {
473                    memcpy(ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr,
474                            &ADOLC_CURRENT_TAPE_INFOS, sizeof(TapeInfos));
475                    ADOLC_TAPE_STACK.push(
476                            ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr);
477                } else {
478                    memcpy(&ADOLC_CURRENT_TAPE_INFOS_FALLBACK,
479                            &ADOLC_CURRENT_TAPE_INFOS, sizeof(TapeInfos));
480                    ADOLC_TAPE_STACK.push(&ADOLC_CURRENT_TAPE_INFOS_FALLBACK);
481                }
482                memcpy(&ADOLC_CURRENT_TAPE_INFOS, *tiIter, sizeof(TapeInfos));
483                ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr = *tiIter;
484                return;
485            }
486        }
487    }
488
489    /* tapeID not used so far */
490    if (mode == ADOLC_REVERSE) {
491        failAdditionalInfo1 = tapeID;
492        fail(ADOLC_REVERSE_NO_TAYLOR_STACK);
493    }
494
495    /* create new info struct and initialize it */
496    tempTapeInfos = new TapeInfos(tapeID);
497    tempTapeInfos->traceFlag=1;
498    tempTapeInfos->inUse = 1;
499    tempTapeInfos->tapingComplete = 1;
500    ADOLC_TAPE_INFOS_BUFFER.push_back(tempTapeInfos);
501
502    read_tape_stats(tempTapeInfos);
503    /* update tapeStack and save tapeInfos */
504    if (ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr != NULL) {
505        memcpy(ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr,
506                &ADOLC_CURRENT_TAPE_INFOS, sizeof(TapeInfos));
507        ADOLC_TAPE_STACK.push(ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr);
508    } else {
509        memcpy(&ADOLC_CURRENT_TAPE_INFOS_FALLBACK,
510                &ADOLC_CURRENT_TAPE_INFOS, sizeof(TapeInfos));
511        ADOLC_TAPE_STACK.push(&ADOLC_CURRENT_TAPE_INFOS_FALLBACK);
512    }
513
514    /* set the new tape infos as current */
515    memcpy(&ADOLC_CURRENT_TAPE_INFOS, tempTapeInfos, sizeof(TapeInfos));
516    ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr = tempTapeInfos;
517}
518
519/* release the current tape and give control to the previous one */
520void releaseTape() {
521    ADOLC_OPENMP_THREAD_NUMBER;
522    ADOLC_OPENMP_GET_THREAD_NUMBER;
523
524    /* if operations, locations and constants tapes have been written and value
525     * stack information have not been created tapeInfos are no longer needed*/
526    if (ADOLC_CURRENT_TAPE_INFOS.keepTaylors            == 0 &&
527            ADOLC_CURRENT_TAPE_INFOS.stats[OP_FILE_ACCESS]  == 1 &&
528            ADOLC_CURRENT_TAPE_INFOS.stats[LOC_FILE_ACCESS] == 1 &&
529            ADOLC_CURRENT_TAPE_INFOS.stats[VAL_FILE_ACCESS] == 1 ) {
530        ADOLC_CURRENT_TAPE_INFOS.inUse = 0;
531    }
532
533    memcpy(ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr,
534            &ADOLC_CURRENT_TAPE_INFOS, sizeof(TapeInfos));
535    ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr = ADOLC_TAPE_STACK.top();
536    memcpy(&ADOLC_CURRENT_TAPE_INFOS,
537            ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr, sizeof(TapeInfos));
538    ADOLC_TAPE_STACK.pop();
539    if (ADOLC_TAPE_STACK.empty())
540        ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr = NULL;
541}
542
543/* updates the tape infos for the given ID - a tapeInfos struct is created
544 * and registered if non is found but its state will remain "not in use" */
545TapeInfos *getTapeInfos(short tapeID) {
546    TapeInfos *tapeInfos;
547    vector<TapeInfos *>::iterator tiIter;
548
549    ADOLC_OPENMP_THREAD_NUMBER;
550    ADOLC_OPENMP_GET_THREAD_NUMBER;
551
552    /* check if TapeInfos for tapeID exist */
553    if (!ADOLC_TAPE_INFOS_BUFFER.empty()) {
554        for (tiIter=ADOLC_TAPE_INFOS_BUFFER.begin();
555                tiIter!=ADOLC_TAPE_INFOS_BUFFER.end();
556                ++tiIter) {
557            if ((*tiIter)->tapeID==tapeID) {
558                tapeInfos=*tiIter;
559                if (tapeInfos->inUse==0) read_tape_stats(tapeInfos);
560                return tapeInfos;
561            }
562        }
563    }
564    /* create new TapeInfos, initialize and update tapeInfosBuffer */
565    tapeInfos = new TapeInfos(tapeID);
566    ADOLC_TAPE_INFOS_BUFFER.push_back(tapeInfos);
567    tapeInfos->traceFlag=1;
568    tapeInfos->inUse=0;
569#ifdef SPARSE
570    tapeInfos->pTapeInfos.inJacSparseUse=0;
571    tapeInfos->pTapeInfos.inHessSparseUse=0;
572#endif
573    tapeInfos->tapingComplete = 1;
574    read_tape_stats(tapeInfos);
575    return tapeInfos;
576}
577
578#ifdef SPARSE
579/* updates the tape infos on sparse Jac for the given ID  */
580void setTapeInfoJacSparse(short tapeID, SparseJacInfos sJinfos) {
581    TapeInfos *tapeInfos;
582    vector<TapeInfos *>::iterator tiIter;
583
584    ADOLC_OPENMP_THREAD_NUMBER;
585    ADOLC_OPENMP_GET_THREAD_NUMBER;
586
587    /* check if TapeInfos for tapeID exist */
588    if (!ADOLC_TAPE_INFOS_BUFFER.empty()) {
589        for (tiIter=ADOLC_TAPE_INFOS_BUFFER.begin();
590                tiIter!=ADOLC_TAPE_INFOS_BUFFER.end();
591                ++tiIter) {
592            if ((*tiIter)->tapeID==tapeID) {
593                tapeInfos=*tiIter;
594                // free memory of tape entry that had been used previously
595                freeSparseJacInfos(tapeInfos->pTapeInfos.sJinfos.y,
596                        tapeInfos->pTapeInfos.sJinfos.B,
597                        tapeInfos->pTapeInfos.sJinfos.JP,
598                        tapeInfos->pTapeInfos.sJinfos.g,
599                        tapeInfos->pTapeInfos.sJinfos.jr1d,
600                        tapeInfos->pTapeInfos.sJinfos.seed_rows,
601                        tapeInfos->pTapeInfos.sJinfos.seed_clms,
602                        tapeInfos->pTapeInfos.sJinfos.depen);
603                tapeInfos->pTapeInfos.sJinfos.y=sJinfos.y;
604                tapeInfos->pTapeInfos.sJinfos.Seed=sJinfos.Seed;
605                tapeInfos->pTapeInfos.sJinfos.B=sJinfos.B;
606                tapeInfos->pTapeInfos.sJinfos.JP=sJinfos.JP;
607                tapeInfos->pTapeInfos.sJinfos.depen=sJinfos.depen;
608                tapeInfos->pTapeInfos.sJinfos.nnz_in=sJinfos.nnz_in;
609                tapeInfos->pTapeInfos.sJinfos.seed_clms=sJinfos.seed_clms;
610                tapeInfos->pTapeInfos.sJinfos.seed_rows=sJinfos.seed_rows;
611                tapeInfos->pTapeInfos.sJinfos.g=sJinfos.g;
612                tapeInfos->pTapeInfos.sJinfos.jr1d=sJinfos.jr1d;
613            }
614        }
615    }
616}
617#endif
618
619#ifdef SPARSE
620/* updates the tape infos on sparse Hess for the given ID  */
621void setTapeInfoHessSparse(short tapeID, SparseHessInfos sHinfos) {
622    TapeInfos *tapeInfos;
623    vector<TapeInfos *>::iterator tiIter;
624
625    ADOLC_OPENMP_THREAD_NUMBER;
626    ADOLC_OPENMP_GET_THREAD_NUMBER;
627
628    /* check if TapeInfos for tapeID exist */
629    if (!ADOLC_TAPE_INFOS_BUFFER.empty()) {
630        for (tiIter=ADOLC_TAPE_INFOS_BUFFER.begin();
631                tiIter!=ADOLC_TAPE_INFOS_BUFFER.end();
632                ++tiIter) {
633            if ((*tiIter)->tapeID==tapeID) {
634                tapeInfos=*tiIter;
635                // free memory of tape entry that had been used previously
636                    freeSparseHessInfos(tapeInfos->pTapeInfos.sHinfos.Hcomp, 
637                                        tapeInfos->pTapeInfos.sHinfos.Xppp, 
638                                        tapeInfos->pTapeInfos.sHinfos.Yppp, 
639                                        tapeInfos->pTapeInfos.sHinfos.Zppp, 
640                                        tapeInfos->pTapeInfos.sHinfos.Upp, 
641                                        tapeInfos->pTapeInfos.sHinfos.HP,
642                                        tapeInfos->pTapeInfos.sHinfos.g, 
643                                        tapeInfos->pTapeInfos.sHinfos.hr, 
644                                        tapeInfos->pTapeInfos.sHinfos.p, 
645                                        tapeInfos->pTapeInfos.sHinfos.indep);   
646                    tapeInfos->pTapeInfos.sHinfos.Hcomp=sHinfos.Hcomp;
647                    tapeInfos->pTapeInfos.sHinfos.Xppp=sHinfos.Xppp;
648                    tapeInfos->pTapeInfos.sHinfos.Yppp=sHinfos.Yppp;
649                    tapeInfos->pTapeInfos.sHinfos.Zppp=sHinfos.Zppp;
650                    tapeInfos->pTapeInfos.sHinfos.Upp=sHinfos.Upp;
651                    tapeInfos->pTapeInfos.sHinfos.HP=sHinfos.HP;
652                    tapeInfos->pTapeInfos.sHinfos.indep=sHinfos.indep;
653                    tapeInfos->pTapeInfos.sHinfos.nnz_in=sHinfos.nnz_in;
654                    tapeInfos->pTapeInfos.sHinfos.p=sHinfos.p;
655                    tapeInfos->pTapeInfos.sHinfos.g=sHinfos.g;
656                    tapeInfos->pTapeInfos.sHinfos.hr=sHinfos.hr;
657            }
658        }
659    }
660}
661#endif
662
663void init() {
664    ADOLC_OPENMP_THREAD_NUMBER;
665    errno = 0;
666    ADOLC_OPENMP_GET_THREAD_NUMBER;
667
668#if defined(_OPENMP)
669    tapeInfosBuffer = new vector<TapeInfos *>;
670    tapeStack = new stack<TapeInfos *>;
671    currentTapeInfos = new TapeInfos;
672    currentTapeInfos->tapingComplete = 1;
673    currentTapeInfos_fallBack = new TapeInfos;
674    globalTapeVars = new GlobalTapeVars;
675    ADOLC_extDiffFctsBuffer = new ADOLC_BUFFER_TYPE;
676    ADOLC_checkpointsStack = new stack<StackElement>;
677    revolve_numbers = new revolve_nums;
678#endif /* _OPENMP */
679
680    ADOLC_CURRENT_TAPE_INFOS.traceFlag = 0;
681    ADOLC_CURRENT_TAPE_INFOS.keepTaylors = 0;
682
683    ADOLC_GLOBAL_TAPE_VARS.maxLoc=1;
684    for (uint i=0; i<sizeof(locint)*8-1; ++i) {
685        ADOLC_GLOBAL_TAPE_VARS.maxLoc<<=1;
686        ++ADOLC_GLOBAL_TAPE_VARS.maxLoc;
687    }
688    ADOLC_GLOBAL_TAPE_VARS.inParallelRegion = 0;
689    ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr = NULL;
690    ADOLC_GLOBAL_TAPE_VARS.branchSwitchWarning = 1;
691
692    adolc_id.adolc_ver    = ADOLC_VERSION;
693    adolc_id.adolc_sub    = ADOLC_SUBVERSION;
694    adolc_id.adolc_lvl    = ADOLC_PATCHLEVEL;
695    adolc_id.locint_size  = sizeof(locint);
696    adolc_id.revreal_size = sizeof(revreal);
697    adolc_id.address_size = sizeof(size_t);
698
699    ADOLC_EXT_DIFF_FCTS_BUFFER.init(init_CpInfos);
700}
701
702/* does things like closing/removing temporary files, ... */
703void cleanUp() {
704    ADOLC_OPENMP_THREAD_NUMBER;
705    ADOLC_OPENMP_GET_THREAD_NUMBER;
706
707    vector<TapeInfos *>::iterator tiIter;
708    if (!ADOLC_TAPE_INFOS_BUFFER.empty()) {
709        for (tiIter=ADOLC_TAPE_INFOS_BUFFER.begin();
710                tiIter!=ADOLC_TAPE_INFOS_BUFFER.end();
711                ++tiIter)
712        {
713            /* close open files though they may be incomplete */
714            if ((*tiIter)->op_file!=NULL)
715            {
716                fclose((*tiIter)->op_file);
717                (*tiIter)->op_file = NULL;
718            }
719            if ((*tiIter)->val_file!=NULL)
720            {
721                fclose((*tiIter)->val_file);
722                (*tiIter)->val_file = NULL;
723            }
724            if ((*tiIter)->loc_file!=NULL)
725            {
726                fclose((*tiIter)->loc_file);
727                (*tiIter)->loc_file = NULL;
728            }
729            if ((*tiIter)->tay_file!=NULL) {
730                fclose((*tiIter)->tay_file);
731                (*tiIter)->tay_file = NULL;
732                remove((*tiIter)->pTapeInfos.tay_fileName);
733            }
734            if ((*tiIter)->opBuffer != NULL)
735            {
736                free((*tiIter)->opBuffer);
737                (*tiIter)->opBuffer = NULL;
738            }
739            if ((*tiIter)->valBuffer != NULL)
740            {
741                free((*tiIter)->valBuffer);
742                (*tiIter)->valBuffer = NULL;
743            }
744            if ((*tiIter)->locBuffer != NULL)
745            {
746                free((*tiIter)->locBuffer);
747                (*tiIter)->locBuffer = NULL;
748            }
749            if ((*tiIter)->tayBuffer != NULL)
750            {
751                free((*tiIter)->tayBuffer);
752                (*tiIter)->tayBuffer = NULL;
753            }
754
755#ifdef SPARSE
756            freeSparseJacInfos((*tiIter)->pTapeInfos.sJinfos.y,
757                               (*tiIter)->pTapeInfos.sJinfos.B,
758                               (*tiIter)->pTapeInfos.sJinfos.JP,
759                               (*tiIter)->pTapeInfos.sJinfos.g,
760                               (*tiIter)->pTapeInfos.sJinfos.jr1d,
761                               (*tiIter)->pTapeInfos.sJinfos.seed_rows,
762                               (*tiIter)->pTapeInfos.sJinfos.seed_clms,
763                               (*tiIter)->pTapeInfos.sJinfos.depen);
764            freeSparseHessInfos((*tiIter)->pTapeInfos.sHinfos.Hcomp, 
765                                (*tiIter)->pTapeInfos.sHinfos.Xppp, 
766                                (*tiIter)->pTapeInfos.sHinfos.Yppp, 
767                                (*tiIter)->pTapeInfos.sHinfos.Zppp, 
768                                (*tiIter)->pTapeInfos.sHinfos.Upp, 
769                                (*tiIter)->pTapeInfos.sHinfos.HP,
770                                (*tiIter)->pTapeInfos.sHinfos.g, 
771                                (*tiIter)->pTapeInfos.sHinfos.hr, 
772                                (*tiIter)->pTapeInfos.sHinfos.p, 
773                                (*tiIter)->pTapeInfos.sHinfos.indep);   
774#endif
775
776            /* remove "main" tape files if not all three have been written */
777            int filesWritten = (*tiIter)->stats[OP_FILE_ACCESS] +
778                (*tiIter)->stats[LOC_FILE_ACCESS] +
779                (*tiIter)->stats[VAL_FILE_ACCESS];
780            if ( (filesWritten > 0) && ((*tiIter)->pTapeInfos.keepTape == 0) )
781            {
782                /* try to remove all tapes (even those not written by this
783                 * run) => this ensures that there is no mixture of tapes from
784                 * different ADOLC runs */
785                if ( (*tiIter)->stats[OP_FILE_ACCESS] == 1 )
786                    remove((*tiIter)->pTapeInfos.op_fileName);
787                if ( (*tiIter)->stats[LOC_FILE_ACCESS] == 1 )
788                    remove((*tiIter)->pTapeInfos.loc_fileName);
789                if ( (*tiIter)->stats[VAL_FILE_ACCESS] == 1 )
790                    remove((*tiIter)->pTapeInfos.val_fileName);
791            }
792            if ((*tiIter)->pTapeInfos.op_fileName != NULL)
793            {
794                free((*tiIter)->pTapeInfos.op_fileName);
795                (*tiIter)->pTapeInfos.op_fileName = NULL;
796            }
797            if ((*tiIter)->pTapeInfos.val_fileName != NULL)
798            {
799                free((*tiIter)->pTapeInfos.val_fileName);
800                (*tiIter)->pTapeInfos.val_fileName = NULL;
801            }
802            if ((*tiIter)->pTapeInfos.loc_fileName != NULL)
803            {
804                free((*tiIter)->pTapeInfos.loc_fileName);
805                (*tiIter)->pTapeInfos.loc_fileName = NULL;
806            }
807            if ((*tiIter)->pTapeInfos.tay_fileName != NULL)
808            {
809                free((*tiIter)->pTapeInfos.tay_fileName);
810                (*tiIter)->pTapeInfos.tay_fileName = NULL;
811            }
812
813            delete *tiIter;
814        }
815    }
816
817    cp_clearStack();
818
819    if (ADOLC_GLOBAL_TAPE_VARS.store != NULL) {
820        delete[] ADOLC_GLOBAL_TAPE_VARS.store;
821        ADOLC_GLOBAL_TAPE_VARS.store = NULL;
822    }
823
824#if defined(_OPENMP)
825    if (ADOLC_GLOBAL_TAPE_VARS.inParallelRegion == 0) {
826        /* cleanup on program exit */
827        delete revolve_numbers;
828        delete ADOLC_checkpointsStack;
829        delete ADOLC_extDiffFctsBuffer;
830        delete globalTapeVars;
831        delete currentTapeInfos;
832        delete currentTapeInfos_fallBack;
833        delete tapeStack;
834        delete tapeInfosBuffer;
835    }
836#endif
837
838    ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
839    clearTapeBaseNames();
840}
841
842int removeTape(short tapeID, short type) {
843    TapeInfos *tapeInfos = NULL;
844    vector<TapeInfos *>::iterator tiIter;
845    ADOLC_OPENMP_THREAD_NUMBER;
846    ADOLC_OPENMP_GET_THREAD_NUMBER;
847
848    /* check if TapeInfos for tapeID exist */
849    if (!ADOLC_TAPE_INFOS_BUFFER.empty()) {
850        for (tiIter = ADOLC_TAPE_INFOS_BUFFER.begin();
851                tiIter != ADOLC_TAPE_INFOS_BUFFER.end();
852                ++tiIter)
853        {
854            if ((*tiIter)->tapeID == tapeID) {
855                tapeInfos = *tiIter;
856                if (tapeInfos->tapingComplete == 0) return -1;
857                ADOLC_TAPE_INFOS_BUFFER.erase(tiIter);
858                break;
859            }
860        }
861    }
862
863    if (tapeInfos == NULL) { // might be on disk only
864        tapeInfos = new TapeInfos(tapeID);
865        tapeInfos->tapingComplete = 1;
866    }
867
868    freeTapeResources(tapeInfos);
869#ifdef SPARSE
870    freeSparseJacInfos(tapeInfos->pTapeInfos.sJinfos.y,
871                       tapeInfos->pTapeInfos.sJinfos.B,
872                       tapeInfos->pTapeInfos.sJinfos.JP,
873                       tapeInfos->pTapeInfos.sJinfos.g,
874                       tapeInfos->pTapeInfos.sJinfos.jr1d,
875                       tapeInfos->pTapeInfos.sJinfos.seed_rows,
876                       tapeInfos->pTapeInfos.sJinfos.seed_clms,
877                       tapeInfos->pTapeInfos.sJinfos.depen);
878    freeSparseHessInfos(tapeInfos->pTapeInfos.sHinfos.Hcomp, 
879                        tapeInfos->pTapeInfos.sHinfos.Xppp, 
880                        tapeInfos->pTapeInfos.sHinfos.Yppp, 
881                        tapeInfos->pTapeInfos.sHinfos.Zppp, 
882                        tapeInfos->pTapeInfos.sHinfos.Upp, 
883                        tapeInfos->pTapeInfos.sHinfos.HP,
884                        tapeInfos->pTapeInfos.sHinfos.g, 
885                        tapeInfos->pTapeInfos.sHinfos.hr, 
886                        tapeInfos->pTapeInfos.sHinfos.p, 
887                        tapeInfos->pTapeInfos.sHinfos.indep);   
888#endif
889    ADOLC_OPENMP_RESTORE_THREAD_NUMBER;
890
891    if (type == ADOLC_REMOVE_COMPLETELY) {
892        remove(tapeInfos->pTapeInfos.op_fileName);
893        remove(tapeInfos->pTapeInfos.loc_fileName);
894        remove(tapeInfos->pTapeInfos.val_fileName);
895    }
896
897    free(tapeInfos->pTapeInfos.op_fileName);
898    free(tapeInfos->pTapeInfos.val_fileName);
899    free(tapeInfos->pTapeInfos.loc_fileName);
900    if (tapeInfos->pTapeInfos.tay_fileName != NULL)
901        free(tapeInfos->pTapeInfos.tay_fileName);
902
903    delete tapeInfos;
904
905    return 0;
906}
907
908/****************************************************************************/
909/* Initialization for the taping process. Creates buffers for this tape,    */
910/* sets files names, and calls appropriate setup routines.                  */
911/****************************************************************************/
912int trace_on(short tnum, int keepTaylors) {
913    int retval = 0;
914    ADOLC_OPENMP_THREAD_NUMBER;
915    ADOLC_OPENMP_GET_THREAD_NUMBER;
916
917    /* allocate memory for TapeInfos and update tapeStack */
918    retval = initNewTape(tnum);
919    ADOLC_CURRENT_TAPE_INFOS.keepTaylors=keepTaylors;
920    if (keepTaylors!=0) ADOLC_CURRENT_TAPE_INFOS.deg_save=1;
921    start_trace();
922    take_stock();               /* record all existing adoubles on the tape */
923    return retval;
924}
925
926int trace_on(short tnum, int keepTaylors,
927        uint obs, uint lbs, uint vbs, uint tbs)
928{
929    int retval = 0;
930    ADOLC_OPENMP_THREAD_NUMBER;
931    ADOLC_OPENMP_GET_THREAD_NUMBER;
932
933    /* allocate memory for TapeInfos and update tapeStack */
934    retval = initNewTape(tnum);
935    ADOLC_CURRENT_TAPE_INFOS.stats[OP_BUFFER_SIZE] = obs;
936    ADOLC_CURRENT_TAPE_INFOS.stats[LOC_BUFFER_SIZE] = lbs;
937    ADOLC_CURRENT_TAPE_INFOS.stats[VAL_BUFFER_SIZE] = vbs;
938    ADOLC_CURRENT_TAPE_INFOS.stats[TAY_BUFFER_SIZE] = tbs;
939    ADOLC_CURRENT_TAPE_INFOS.keepTaylors=keepTaylors;
940    if (keepTaylors!=0) ADOLC_CURRENT_TAPE_INFOS.deg_save=1;
941    start_trace();
942    take_stock();               /* record all existing adoubles on the tape */
943    return retval;
944}
945
946/****************************************************************************/
947/* Stop Tracing. Cleans up, and turns off trace_flag. Flag not equal zero   */
948/* enforces writing of the three main tape files (op+loc+val).              */
949/****************************************************************************/
950void trace_off(int flag) {
951    ADOLC_OPENMP_THREAD_NUMBER;
952    ADOLC_OPENMP_GET_THREAD_NUMBER;
953    ADOLC_CURRENT_TAPE_INFOS.pTapeInfos.keepTape = flag;
954    keep_stock();         /* copy remaining live variables + trace_flag = 0 */
955    stop_trace(flag);
956    cout.flush();
957    ADOLC_CURRENT_TAPE_INFOS.tapingComplete = 1;
958    releaseTape();
959}
960
961bool isTaping() {
962    ADOLC_OPENMP_THREAD_NUMBER;
963    ADOLC_OPENMP_GET_THREAD_NUMBER;
964    return ADOLC_CURRENT_TAPE_INFOS.traceFlag != 0;
965}
966
967
968
969/****************************************************************************/
970/* A class for initialization/finalization and OpenMP handling              */
971/****************************************************************************/
972class Keeper {
973    public:
974        inline Keeper() {
975            dummy = 0;
976            init();
977            readConfigFile();
978        }
979        inline ~Keeper() {
980            cleanUp();
981        }
982
983        inline void touch() {
984            dummy = 1;
985        }
986
987    private:
988        int dummy;
989};
990
991/* a static instance that does all work */
992static Keeper theKeeper;
993
994/**
995 * Hope to convince the linker to link the keeper code into the executable. */
996void initADOLC() {
997    theKeeper.touch();
998}
999
1000/****************************************************************************/
1001/****************************************************************************/
1002/* The following is necessary to provide a separate ADOL-C environment for  */
1003/* each OpenMP worker.                                                      */
1004/****************************************************************************/
1005/****************************************************************************/
1006#if defined(_OPENMP)
1007#include "adolc_openmp.h"
1008
1009ADOLC_OpenMP ADOLC_OpenMP_Handler;
1010ADOLC_OpenMP_NC ADOLC_OpenMP_Handler_NC;
1011int ADOLC_parallel_doCopy;
1012
1013static bool waitForMaster_begin = true;
1014static bool waitForMaster_end   = true;
1015static bool firstParallel       = true;
1016
1017/****************************************************************************/
1018/* Used by OpenMP to create a separate environment for every worker thread. */
1019/****************************************************************************/
1020void beginParallel() {
1021    ADOLC_OPENMP_THREAD_NUMBER;
1022#if defined(ADOLC_THREADSAVE_ERRNO)
1023    errno = omp_get_thread_num();
1024#endif
1025    ADOLC_OPENMP_GET_THREAD_NUMBER;
1026
1027    if (ADOLC_threadNumber == 0) { /* master only */
1028        int numThreads = omp_get_num_threads();
1029
1030        tapeInfosBuffer_s           = tapeInfosBuffer;
1031        tapeStack_s                 = tapeStack;
1032        currentTapeInfos_s          = currentTapeInfos;
1033        currentTapeInfos_fallBack_s = currentTapeInfos_fallBack;
1034        globalTapeVars_s            = globalTapeVars;
1035        ADOLC_extDiffFctsBuffer_s   = ADOLC_extDiffFctsBuffer;
1036        ADOLC_checkpointsStack_s    = ADOLC_checkpointsStack;
1037        revolve_numbers_s           = revolve_numbers;
1038
1039        if (firstParallel) {
1040            tapeInfosBuffer           = new vector<TapeInfos *>[numThreads];
1041            tapeStack                 = new stack<TapeInfos *>[numThreads];
1042            currentTapeInfos          = new TapeInfos[numThreads];
1043            currentTapeInfos_fallBack = new TapeInfos[numThreads];
1044            globalTapeVars            = new GlobalTapeVars[numThreads];
1045            ADOLC_extDiffFctsBuffer   = new ADOLC_BUFFER_TYPE[numThreads];
1046            ADOLC_checkpointsStack    = new stack<StackElement>[numThreads];
1047            revolve_numbers           = new revolve_nums[numThreads];
1048        } else {
1049            tapeInfosBuffer           = tapeInfosBuffer_p;
1050            tapeStack                 = tapeStack_p;
1051            currentTapeInfos          = currentTapeInfos_p;
1052            currentTapeInfos_fallBack = currentTapeInfos_fallBack_p;
1053            globalTapeVars            = globalTapeVars_p;
1054            ADOLC_extDiffFctsBuffer   = ADOLC_extDiffFctsBuffer_p;
1055            ADOLC_checkpointsStack    = ADOLC_checkpointsStack_p;
1056            revolve_numbers         = revolve_numbers_p;
1057        }
1058
1059        /* - set inParallelRegion for tmpGlobalTapeVars because it is source
1060         *   for initializing the parallel globalTapeVars structs
1061         * - inParallelRegion has to be set to one for all workers by master.
1062         *   This is necessary, to deter a speedy master from assuming all
1063         *   workers are done, in endParallel, before they even leaved
1064         *   beginParallel. */
1065        globalTapeVars_s[0].inParallelRegion = 1;
1066        for (int i = 0; i < numThreads; ++i)
1067            globalTapeVars[i].inParallelRegion = 1;
1068
1069        waitForMaster_end = true;
1070        waitForMaster_begin = false;
1071    } else 
1072        while (waitForMaster_begin) {
1073            usleep(1000); /* if anyone knows a better value, ... :-) */
1074        }
1075
1076    if (firstParallel) {
1077        ADOLC_EXT_DIFF_FCTS_BUFFER.init(init_CpInfos);
1078
1079        /* Use assignment operator instead of open coding
1080         * this copies the store and the storemanager too
1081         */
1082        ADOLC_GLOBAL_TAPE_VARS = *globalTapeVars_s;
1083
1084        ADOLC_GLOBAL_TAPE_VARS.newTape = 0;
1085        ADOLC_CURRENT_TAPE_INFOS.tapingComplete = 1;
1086        ADOLC_GLOBAL_TAPE_VARS.currentTapeInfosPtr = NULL;
1087    } else {
1088        if (ADOLC_parallel_doCopy) {
1089            ADOLC_GLOBAL_TAPE_VARS.storeSize = globalTapeVars_s->storeSize;
1090            ADOLC_GLOBAL_TAPE_VARS.numLives = globalTapeVars_s->numLives;
1091           
1092            ADOLC_GLOBAL_TAPE_VARS.branchSwitchWarning = globalTapeVars_s->branchSwitchWarning;
1093
1094            /* deleting the storemanager deletes the store too */
1095            delete ADOLC_GLOBAL_TAPE_VARS.storeManagerPtr;
1096
1097            ADOLC_GLOBAL_TAPE_VARS.store = new
1098                double[ADOLC_GLOBAL_TAPE_VARS.storeSize];
1099            memcpy(ADOLC_GLOBAL_TAPE_VARS.store, globalTapeVars_s->store,
1100                    ADOLC_GLOBAL_TAPE_VARS.storeSize * sizeof(double));
1101            ADOLC_GLOBAL_TAPE_VARS.storeManagerPtr = new
1102                StoreManagerLocintBlock(
1103                    dynamic_cast<StoreManagerLocintBlock*>(globalTapeVars_s->storeManagerPtr),
1104                    ADOLC_GLOBAL_TAPE_VARS.store,
1105                    ADOLC_GLOBAL_TAPE_VARS.storeSize,
1106                    ADOLC_GLOBAL_TAPE_VARS.numLives);
1107        }
1108    }
1109}
1110
1111/****************************************************************************/
1112/* Used by OpenMP to destroy the separate environment of every worker.      */
1113/****************************************************************************/
1114/* There are n+1 instances of ADOLC_OpenMP => n within the parallel region
1115 * and one in the serial part! */
1116void endParallel() {
1117    ADOLC_OPENMP_THREAD_NUMBER;
1118    ADOLC_OPENMP_GET_THREAD_NUMBER;
1119
1120    /* do nothing if called at program exit (serial part) */
1121    if (ADOLC_threadNumber == 0 &&
1122            ADOLC_GLOBAL_TAPE_VARS.inParallelRegion == 0) return;
1123
1124    ADOLC_GLOBAL_TAPE_VARS.inParallelRegion = 0;
1125
1126    if (ADOLC_threadNumber == 0) { /* master only */
1127        int num;
1128        int numThreads = omp_get_num_threads();
1129        bool firstIt = true;
1130        do { /* wait until all slaves have left the parallel part */
1131            if (firstIt) firstIt = false;
1132            else usleep(1000); /* no busy waiting */
1133            num = 1;
1134            for (int i = 1; i < numThreads; ++i)
1135                if (globalTapeVars[i].inParallelRegion == 0) ++num;
1136        } while (num != numThreads);
1137
1138        firstParallel = false;
1139
1140        revolve_numbers_p           = revolve_numbers;
1141        ADOLC_checkpointsStack_p    = ADOLC_checkpointsStack;
1142        ADOLC_extDiffFctsBuffer_p   = ADOLC_extDiffFctsBuffer;
1143        globalTapeVars_p            = globalTapeVars;
1144        currentTapeInfos_p          = currentTapeInfos;
1145        currentTapeInfos_fallBack_p = currentTapeInfos_fallBack;
1146        tapeStack_p                 = tapeStack;
1147        tapeInfosBuffer_p           = tapeInfosBuffer;
1148
1149        revolve_numbers           = revolve_numbers_s;
1150        ADOLC_checkpointsStack    = ADOLC_checkpointsStack_s;
1151        ADOLC_extDiffFctsBuffer   = ADOLC_extDiffFctsBuffer_s;
1152        globalTapeVars            = globalTapeVars_s;
1153        currentTapeInfos          = currentTapeInfos_s;
1154        currentTapeInfos_fallBack = currentTapeInfos_fallBack_s;
1155        tapeStack                 = tapeStack_s;
1156        tapeInfosBuffer           = tapeInfosBuffer_s;
1157
1158        ADOLC_GLOBAL_TAPE_VARS.inParallelRegion = 0;
1159        waitForMaster_begin = true;
1160        waitForMaster_end = false;
1161    } else
1162        while (waitForMaster_end) {
1163            usleep(1000); // no busy waiting
1164        }
1165}
1166
1167#endif /* _OPENMP */
1168
1169TapeInfos::TapeInfos() {
1170    initTapeInfos(this);
1171}
1172
1173TapeInfos::TapeInfos(short _tapeID) {
1174    initTapeInfos(this);
1175    tapeID = _tapeID;
1176    pTapeInfos.op_fileName = createFileName(tapeID, OPERATIONS_TAPE);
1177    pTapeInfos.loc_fileName = createFileName(tapeID, LOCATIONS_TAPE);
1178    pTapeInfos.val_fileName = createFileName(tapeID, VALUES_TAPE);
1179    pTapeInfos.tay_fileName = NULL;
1180}
1181
1182StoreManagerLocintBlock::StoreManagerLocintBlock(double * &storePtr, size_t &size, size_t &numlives) :
1183    storePtr(storePtr),
1184    maxsize(size), currentfill(numlives)
1185{
1186    indexFree.clear();
1187#ifdef ADOLC_DEBUG
1188    std::cerr << "StoreManagerIntegerBlock::StoreManagerIntegerBlock()\n";
1189#endif
1190}
1191
1192StoreManagerLocintBlock::~StoreManagerLocintBlock()
1193{
1194#ifdef ADOLC_DEBUG
1195    std::cerr << "StoreManagerIntegerBlock::~StoreManagerIntegerBlock()\n";
1196#endif
1197    if (storePtr) {
1198     delete[] storePtr;
1199     storePtr = 0;
1200    }
1201    if (!indexFree.empty() ) {
1202        indexFree.clear();
1203    }
1204    maxsize = 0;
1205    currentfill = 0;
1206}
1207
1208StoreManagerLocintBlock::StoreManagerLocintBlock(
1209    const StoreManagerLocintBlock *const stm,
1210    double * &storePtr, size_t &size, size_t &numlives) :
1211    storePtr(storePtr),
1212    maxsize(size), currentfill(numlives)
1213{
1214#ifdef ADOLC_DEBUG
1215    std::cerr << "StoreManagerInteger::StoreManagerInteger()\n";
1216#endif
1217    indexFree.clear();
1218    list<struct FreeBlock>::const_iterator iter = stm->indexFree.begin();
1219    for (; iter != stm->indexFree.end(); iter++)
1220        indexFree.push_back( *iter );
1221}
1222
1223
1224locint StoreManagerLocintBlock::next_loc() {
1225    if ( indexFree.empty() )
1226        grow();
1227
1228    locint const result = indexFree.front().next;
1229    indexFree.front().next++;
1230    indexFree.front().size--;
1231
1232    ++currentfill;
1233
1234#ifdef ADOLC_DEBUG
1235    std::cerr << "StoreManagerLocintBlock::next_loc: result: " << result << " fill: " << size() << "max: " << maxSize() << endl;
1236    list<struct FreeBlock>::iterator iter = indexFree.begin();
1237    for( ; iter != indexFree.end(); iter++ )
1238       std::cerr << "INDEXFELD ( " << iter->next << " , " << iter->size << ")" << endl;
1239#endif
1240
1241    if (indexFree.front().size == 0) {
1242        if (indexFree.size() <= 1)
1243            grow();
1244        else
1245          indexFree.pop_front();
1246    }
1247
1248    return result;
1249}
1250
1251void StoreManagerLocintBlock::ensure_block(size_t n) {
1252    bool found = false;
1253
1254#ifdef ADOLC_DEBUG
1255    std::cerr << "StoreManagerLocintBlock::ensure_Block: required " << n << " ... ";
1256    std::cerr << "searching for big enough block " << endl;
1257#endif
1258    list<struct FreeBlock>::iterator iter = indexFree.begin();
1259    for (; iter != indexFree.end() ; iter++ ) {
1260        if ( iter->size >= n) {
1261            if (iter != indexFree.begin() ) {
1262                struct FreeBlock tmp(*iter);
1263                iter = indexFree.erase(iter);
1264                indexFree.push_front(tmp);
1265            }
1266            found = true;
1267            break;
1268        }
1269    }
1270    if (!found) {
1271#ifdef ADOLC_DEBUG
1272        std::cerr << "no big enough block...growing " << endl;
1273#endif
1274        grow(n);
1275    }
1276
1277#ifdef ADOLC_DEBUG
1278    std::cerr << "StoreManagerLocintBlock::ensure_Block: " << " fill: " << size() << "max: " << maxSize() <<  " ensure_Block (" << n << ")" << endl;
1279    iter = indexFree.begin();
1280    for( ; iter != indexFree.end(); iter++ )
1281        std::cerr << "INDEXFELD ( " << iter->next << " , " << iter->size << ")" << endl;
1282#endif
1283}
1284
1285void StoreManagerLocintBlock::grow(size_t minGrow) {
1286    // first figure out what eventual size we want
1287    size_t const oldMaxsize = maxsize;
1288
1289    if (maxsize == 0){
1290        maxsize = initialSize;
1291    } else {
1292        maxsize *= 2;
1293    }
1294
1295    if (minGrow > 0) {
1296        while (maxsize - oldMaxsize < minGrow) {
1297            maxsize *= 2;
1298        }
1299    }
1300
1301    if (maxsize > std::numeric_limits<locint>::max()) {
1302      // encapsulate this error message
1303      fprintf(DIAG_OUT,"\nADOL-C error:\n");
1304      fprintf(DIAG_OUT,"maximal number (%u) of live active variables exceeded\n\n",
1305           std::numeric_limits<locint>::max());
1306      exit(-3);
1307    }
1308
1309#ifdef ADOLC_DEBUG
1310    // index 0 is not used, means one slot less
1311    std::cerr << "StoreManagerIntegerBlock::grow(): increase size from " << oldMaxsize
1312      << " to " << maxsize << " entries (currently " << size() << " entries used)\n";
1313#endif
1314
1315    double *const oldStore = storePtr;
1316
1317#if defined(ADOLC_DEBUG)
1318    std::cerr << "StoreManagerInteger::grow(): allocate " << maxsize * sizeof(double) << " B doubles\n";
1319#endif
1320    storePtr = new double[maxsize];
1321    assert(storePtr);
1322    memset(storePtr, 0, maxsize*sizeof(double));
1323
1324    if (oldStore != NULL) { // not the first time
1325#if defined(ADOLC_DEBUG)
1326      std::cerr << "StoreManagerInteger::grow(): copy values\n";
1327#endif
1328
1329      memcpy(storePtr, oldStore, oldMaxsize*sizeof(double));
1330
1331#if defined(ADOLC_DEBUG)
1332      std::cerr << "StoreManagerInteger::grow(): free " << oldMaxsize * sizeof(double) << "\n";
1333#endif
1334      delete [] oldStore;
1335    }
1336
1337    bool foundTail = false;
1338    list<struct FreeBlock>::iterator iter = indexFree.begin();
1339    for (; iter != indexFree.end() ; iter++ ) {
1340         if (iter->next + iter->size == oldMaxsize ) {
1341             iter->size += (maxsize - oldMaxsize);
1342              struct FreeBlock tmp(*iter);
1343              iter = indexFree.erase(iter);
1344              indexFree.push_front(tmp);
1345              foundTail = true;
1346              break;
1347         }
1348    }
1349
1350    if (! foundTail) {
1351        struct FreeBlock tmp;
1352        tmp.next = oldMaxsize;
1353        tmp.size = (maxsize - oldMaxsize);
1354        indexFree.push_front(tmp);
1355    }
1356
1357    iter = indexFree.begin();
1358    while (iter != indexFree.end()) {
1359         if (iter->size == 0)
1360             iter=indexFree.erase(iter); // don't leave 0 blocks around
1361         else
1362             iter++;
1363    }
1364#ifdef ADOLC_DEBUG
1365    std::cerr << "Growing:" << endl;
1366    iter = indexFree.begin();
1367    for( ; iter != indexFree.end(); iter++ )
1368       std::cerr << "INDEXFELD ( " << iter->next << " , " << iter->size << ")" << endl;
1369#endif
1370}
1371
1372void StoreManagerLocintBlock::free_loc(locint loc) {
1373    assert( loc < maxsize);
1374
1375    list<struct FreeBlock>::iterator iter = indexFree.begin();
1376    for (; iter != indexFree.end() ; iter++ ) {
1377         if (loc+1 == iter->next || iter->next + iter->size == loc) {
1378              iter->size++;
1379              if (loc + 1 == iter->next)
1380                   iter->next = loc;
1381    // bringing the matched element to the front maybe a good idea
1382    // in case several contiguous adouble are deallcated right after
1383    // one another, e.g. advector
1384              struct FreeBlock tmp(*iter);
1385              iter = indexFree.erase(iter);
1386              indexFree.push_front(tmp);
1387              iter = indexFree.begin();
1388              break;
1389         }
1390    }
1391    if (iter == indexFree.end()) {
1392         struct FreeBlock tmp;
1393         tmp.next = loc;
1394         tmp.size = 1;
1395         indexFree.push_front(tmp);
1396    }
1397
1398    --currentfill;
1399#ifdef ADOLC_DEBUG
1400    std::cerr << "free_loc: " << loc << " fill: " << size() << "max: " << maxSize() << endl;
1401
1402    iter = indexFree.begin();
1403    for( ; iter != indexFree.end(); iter++ )
1404       std::cerr << "INDEXFELD ( " << iter->next << " , " << iter->size << ")" << endl;
1405#endif
1406}
1407
1408void ensureContiguousLocations(size_t n) {
1409    ADOLC_OPENMP_THREAD_NUMBER;
1410    ADOLC_OPENMP_GET_THREAD_NUMBER;
1411    ADOLC_GLOBAL_TAPE_VARS.storeManagerPtr->ensure_block(n);
1412}
Note: See TracBrowser for help on using the repository browser.