Back to home page

MITgcm

 
 

    


File indexing completed on 2023-05-28 05:10:58 UTC

view on githubraw file Latest commit b4daa243 on 2023-05-28 03:53:22 UTC
b4daa24319 Shre*0001 /*
                0002  * TAPENADE Automatic Differentiation Engine
                0003  * Copyright (C) 1999-2021 Inria
                0004  * See the LICENSE.md file in the project root for more information.
                0005  *
                0006  */
                0007 
                0008 #include <stdlib.h>
                0009 #include <stdio.h>
                0010 #include <string.h>
                0011 
                0012 /* The size of a BLOCK in characters. Suggested 16384. Should try 2^16=65536 */
                0013 #define ONE_BLOCK_SIZE 65536
                0014 
                0015 /* The main stack is a double-chain of DoubleChainedBlock objects.
                0016  * Each DoubleChainedBlock holds an array[ONE_BLOCK_SIZE] of char. */
                0017 typedef struct _DoubleChainedBlock{
                0018   unsigned int rank ;
                0019   struct _DoubleChainedBlock *prev ;
                0020   char                       *contents ;
                0021   struct _DoubleChainedBlock *next ;
                0022 } DoubleChainedBlock ;
                0023 
                0024 char initContents[ONE_BLOCK_SIZE] = {'\0'} ;
                0025 DoubleChainedBlock initBlock = {0,NULL,initContents,NULL} ;
                0026 static DoubleChainedBlock *curStack = &initBlock ;
                0027 static char               *curStackTop = initContents ;
                0028 
                0029 static unsigned long int maintraffic = 0 ;
                0030 
                0031 void setCurLocation(unsigned long int location) {
                0032   unsigned int targetRank = (unsigned int)location/ONE_BLOCK_SIZE ;
                0033   unsigned int targetOffset = (unsigned int)location%ONE_BLOCK_SIZE ;
                0034   if (targetRank>curStack->rank)
                0035     while (targetRank>curStack->rank) curStack = curStack->next ;
                0036   else if (targetRank<curStack->rank)
                0037     while (targetRank<curStack->rank) curStack = curStack->prev ;
                0038   curStackTop = curStack->contents + targetOffset ;
                0039 }
                0040 
                0041 unsigned long int getCurLocation() {
                0042   return (curStackTop-curStack->contents)+curStack->rank*ONE_BLOCK_SIZE ;
                0043 }
                0044 
                0045 void showLocation(unsigned long int location) {
                0046   printf("%1i.%05i", (unsigned int)location/ONE_BLOCK_SIZE, (unsigned int)location%ONE_BLOCK_SIZE) ;
                0047 }
                0048 
                0049 /*************** REPEATED ACCESS MECHANISM *********************/
                0050 
                0051 typedef struct _StackRepeatCell {
                0052   int hasBackPop ;
                0053   unsigned long int backPop ;
                0054   unsigned long int resume ;
                0055   unsigned long int freePush ;
                0056   struct _StackRepeatCell *previous ;
                0057 } StackRepeatCell ;
                0058 
                0059 StackRepeatCell *stackRepeatTop = NULL ;
                0060 
                0061 void showStackRepeatsRec(StackRepeatCell *inRepeatStack) {
                0062   if (inRepeatStack->previous) {showStackRepeatsRec(inRepeatStack->previous) ; printf(" ; ") ;}
                0063   printf("<") ;
                0064   if (inRepeatStack->hasBackPop) showLocation(inRepeatStack->backPop) ;
                0065   printf("|") ;
                0066   showLocation(inRepeatStack->resume) ;
                0067   printf("|") ;
                0068   showLocation(inRepeatStack->freePush) ;
                0069   printf(">") ;
                0070 }
                0071 
                0072 void showStackRepeats() {
                0073   showStackRepeatsRec(stackRepeatTop) ;
                0074 }
                0075 
                0076 void showStack() {
                0077   DoubleChainedBlock *inStack = &initBlock ;
                0078   int i ;
                0079   while (inStack) {
                0080     printf("[%1i] ",inStack->rank) ;
                0081     for (i=0 ; i<ONE_BLOCK_SIZE ; ++i) {
                0082       if (i!=0 && i%4==0) printf(".") ;
                0083       if (inStack==curStack && &(inStack->contents[i])==curStackTop) printf(" | ") ;
                0084       printf("%02x",(unsigned char)inStack->contents[i]) ;
                0085     }
                0086     inStack = inStack->next ;
                0087     if (inStack) printf("\n        ") ;
                0088   }
                0089   printf("\n        REPEATS:") ;
                0090   if (stackRepeatTop)
                0091     showStackRepeats() ;
                0092   else
                0093     printf(" none!") ;
                0094   printf("\n") ;
                0095 }
                0096 
                0097 void showStackSize(int i4i, int i8i, int r4i, int r8i, int c8i, int c16i, int s1i, int biti, int ptri, int pos) {
                0098   printf(" --%5i--> <",pos) ;
                0099   showLocation(getCurLocation()) ;
                0100   printf(">%1i.%1i.%1i.%1i.%1i.%1i.%1i.%1i.%1i\n",i4i, i8i, r4i, r8i, c8i, c16i, s1i, biti, ptri) ;
                0101 }
                0102 
                0103 void adStack_showPeakSize() {
                0104   DoubleChainedBlock *inStack = &initBlock ;
                0105   int i = 0 ;
                0106   while (inStack) {
                0107     inStack = inStack->next ;
                0108     ++i ;
                0109   }
                0110   printf("Peak stack size (%1i blocks): %1llu bytes\n",
                0111          i, ((long long int)i)*((long long int)ONE_BLOCK_SIZE)) ;
                0112 }
                0113 
                0114 void showTotalTraffic(unsigned long long int localtraffic) {
                0115   printf("Total pushed traffic %1llu bytes\n", maintraffic+localtraffic) ;
                0116 }
                0117 
                0118 /** If we are in a protected, read-only section, memorize location as "backPop"
                0119  * and go to the "freePush" location */
                0120 void checkPushInReadOnly() {
                0121   if (stackRepeatTop) {
                0122     unsigned long int current = getCurLocation() ;
                0123     if (current<stackRepeatTop->freePush) {
                0124       stackRepeatTop->hasBackPop = 1 ;
                0125       stackRepeatTop->backPop = current ;
                0126       setCurLocation(stackRepeatTop->freePush) ;
                0127 /*       printf(" FREEPUSH(") ;                   //Trace */
                0128 /*       showLocation(stackRepeatTop->backPop) ;  //Trace */
                0129 /*       printf("=>") ;                           //Trace */
                0130 /*       showLocation(stackRepeatTop->freePush) ; //Trace */
                0131 /*       printf(")") ;                            //Trace */
                0132     }
                0133   }
                0134 }
                0135 
                0136 /** If current location is the "freePush" location,
                0137  * go back to its "backPop" location, which is in a protected, read-only section */
                0138 void checkPopToReadOnly() {
                0139   if (stackRepeatTop && stackRepeatTop->hasBackPop) {
                0140     unsigned long int current = getCurLocation() ;
                0141     if (current==stackRepeatTop->freePush) {
                0142       setCurLocation(stackRepeatTop->backPop) ;
                0143       stackRepeatTop->hasBackPop = 0 ;
                0144 /*       printf(" BACKPOP(") ;                    //Trace */
                0145 /*       showLocation(stackRepeatTop->freePush) ; //Trace */
                0146 /*       printf("=>") ;                           //Trace */
                0147 /*       showLocation(stackRepeatTop->backPop) ;  //Trace */
                0148 /*       printf(")") ;                            //Trace */
                0149     }
                0150   }
                0151 }
                0152 
                0153 // A global for communication from startStackRepeat1() to startStackRepeat2():
                0154 StackRepeatCell *newRepeatCell = NULL ;
                0155 
                0156 void startStackRepeat1() {
                0157   // Create (push) a new "stack" repeat level:
                0158   newRepeatCell = (StackRepeatCell *)malloc(sizeof(StackRepeatCell)) ;
                0159   newRepeatCell->previous = stackRepeatTop ;
                0160   newRepeatCell->hasBackPop = 0 ;
                0161   // Store current location as the "resume" location:
                0162   unsigned long int current = getCurLocation() ;
                0163   newRepeatCell->resume = current ;
                0164   // Move to the "freePush" location if there is one:
                0165   if (stackRepeatTop && current<stackRepeatTop->freePush)
                0166     setCurLocation(stackRepeatTop->freePush) ;
                0167 }
                0168 
                0169 void startStackRepeat2() {
                0170   // Store current stack location as the "freePush" location:
                0171   newRepeatCell->freePush = getCurLocation() ;
                0172   // Reset current location to stored "resume" location:
                0173   setCurLocation(newRepeatCell->resume) ;
                0174   // Make this new repeat level the current repeat level:
                0175   stackRepeatTop = newRepeatCell ;
                0176 /*   printf("\n+Rep ") ; showStackRepeats() ; printf("\n") ; //Trace */
                0177 }
                0178 
                0179 void resetStackRepeat1() {
                0180 /*   printf("\n>Rep ") ; showStackRepeats() ; printf("\n") ; //Trace */
                0181   // If we are in a nested checkpoint, force exit from it:
                0182   if (stackRepeatTop->hasBackPop) {
                0183     //setCurLocation(stackRepeatTop->backPop) ; //correct but useless code
                0184     stackRepeatTop->hasBackPop = 0 ;
                0185   }
                0186   // Go to repeat location of current repeat level
                0187   setCurLocation(stackRepeatTop->freePush) ;
                0188 }
                0189 
                0190 void resetStackRepeat2() {
                0191   // Reset current location to "ResumeLocation":
                0192   setCurLocation(stackRepeatTop->resume) ;
                0193 }
                0194 
                0195 void endStackRepeat() {
                0196 /*   printf("\n-Rep ") ; showStackRepeats() ; printf("\n") ; //Trace */
                0197   // If we are in a nested checkpoint, go back to its "backPop" (read-only) location:
                0198   if (stackRepeatTop->hasBackPop) {
                0199     setCurLocation(stackRepeatTop->backPop) ;
                0200     //stackRepeatTop->hasBackPop = 0 ; //correct but useless code
                0201   }
                0202   // Remove (pop) top "stack" repeat level:
                0203   StackRepeatCell *oldRepeatCell = stackRepeatTop ;
                0204   stackRepeatTop = stackRepeatTop->previous ;
                0205   free(oldRepeatCell) ;
                0206   // current location may have moved back ; check if we must move further back:
                0207   checkPopToReadOnly() ;
                0208 }
                0209 
                0210 /******************* PUSH/POP MECHANISM *******************/
                0211 
                0212 /* PUSHes "nbChars" consecutive chars from a location starting at address "x".
                0213  * Checks that there is enough space left to hold "nbChars" chars.
                0214  * Otherwise, allocates the necessary space. */
                0215 void pushNArray(char *x, unsigned int nbChars, int checkReadOnly) {
                0216   if (checkReadOnly) checkPushInReadOnly() ;
                0217   if (checkReadOnly) maintraffic += nbChars ;
                0218 /* unsigned long int lfrom = getCurLocation() ; //Trace */
                0219   unsigned int nbmax = ONE_BLOCK_SIZE-(curStackTop-(curStack->contents)) ;
                0220   if (nbChars <= nbmax) {
                0221     memcpy(curStackTop,x,nbChars) ;
                0222     curStackTop+=nbChars ;
                0223   } else {
                0224     char *inx = x+(nbChars-nbmax) ;
                0225     if (nbmax>0) memcpy(curStackTop,inx,nbmax) ;
                0226     while (inx>x) {
                0227       if (curStack->next)
                0228         curStack = curStack->next ;
                0229       else {
                0230         /* Create new block: */
                0231     DoubleChainedBlock *newStack ;
                0232     char *contents = (char *)malloc(ONE_BLOCK_SIZE*sizeof(char)) ;
                0233     newStack = (DoubleChainedBlock*)malloc(sizeof(DoubleChainedBlock)) ;
                0234     if ((contents == NULL) || (newStack == NULL)) {
                0235       DoubleChainedBlock *stack = curStack ;
                0236       int nbBlocks = (stack?-1:0) ;
                0237       while(stack) {
                0238           stack = stack->prev ;
                0239           nbBlocks++ ;
                0240       }
                0241       printf("Out of memory (allocated %i blocks of %i bytes)\n",
                0242          nbBlocks, ONE_BLOCK_SIZE) ;
                0243           exit(0);
                0244     }
                0245         curStack->next = newStack ;
                0246     newStack->prev = curStack ;
                0247         newStack->rank = curStack->rank + 1 ;
                0248     newStack->next = NULL ;
                0249     newStack->contents = contents ;
                0250     curStack = newStack ;
                0251         /* new block created! */
                0252       }
                0253       inx -= ONE_BLOCK_SIZE ;
                0254       if(inx>x)
                0255         memcpy(curStack->contents,inx,ONE_BLOCK_SIZE) ;
                0256       else {
                0257         unsigned int nbhead = (inx-x)+ONE_BLOCK_SIZE ;
                0258         curStackTop = curStack->contents ;
                0259         memcpy(curStackTop,x,nbhead) ;
                0260         curStackTop += nbhead ;
                0261       }
                0262     }
                0263   }
                0264 /* unsigned long int lto = getCurLocation() ; //Trace */
                0265 /* printf("pushNArray(") ;                    //Trace */
                0266 /* showLocation(lfrom) ;                      //Trace */
                0267 /* printf("=>") ;                             //Trace */
                0268 /* showLocation(lto) ;                        //Trace */
                0269 /* printf(")") ;                              //Trace */
                0270 }
                0271 
                0272 /* POPs "nbChars" consecutive chars to a location starting at address "x".
                0273  * Checks that there is enough data to fill "nbChars" chars.
                0274  * Otherwise, pops as many blocks as necessary. */
                0275 void popNArray(char *x, unsigned int nbChars, int checkReadOnly) {
                0276 /* unsigned long int lfrom = getCurLocation() ; //Trace */
                0277   unsigned int nbmax = curStackTop-(curStack->contents) ;
                0278   if (nbChars <= nbmax) {
                0279     curStackTop-=nbChars ;
                0280     memcpy(x,curStackTop,nbChars);
                0281   } else {
                0282     char *tlx = x+nbChars ;
                0283     if (nbmax>0) memcpy(x,curStack->contents,nbmax) ;
                0284     x+=nbmax ;
                0285     while (x<tlx) {
                0286       curStack = curStack->prev ;
                0287       if (curStack==NULL) printf("Popping from an empty stack!!!\n") ;
                0288       if (x+ONE_BLOCK_SIZE<tlx) {
                0289     memcpy(x,curStack->contents,ONE_BLOCK_SIZE) ;
                0290     x += ONE_BLOCK_SIZE ;
                0291       } else {
                0292     unsigned int nbtail = tlx-x ;
                0293     curStackTop = (curStack->contents)+ONE_BLOCK_SIZE-nbtail ;
                0294     memcpy(x,curStackTop,nbtail) ;
                0295     x = tlx ;
                0296       }
                0297     }
                0298   }
                0299 /* unsigned long int lto = getCurLocation() ; //Trace */
                0300 /* printf("popNArray(") ;                     //Trace */
                0301 /* showLocation(lfrom) ;                      //Trace */
                0302 /* printf("=>") ;                             //Trace */
                0303 /* showLocation(lto) ;                        //Trace */
                0304 /* printf(")") ;                              //Trace */
                0305   if (checkReadOnly) checkPopToReadOnly() ;
                0306 }
                0307 
                0308 typedef struct {float r,i;} ccmplx ;
                0309 typedef struct {double dr, di;} cdcmplx ;
                0310 
                0311 void pushInteger4Array(int *x, int n) {
                0312   pushNArray((char *)x,(unsigned int)(n*4), 1) ;
                0313 }
                0314 
                0315 void popInteger4Array(int *x, int n) {
                0316   popNArray((char *)x,(unsigned int)(n*4), 1) ;
                0317 }
                0318 
                0319 void pushInteger8Array(long *x, int n) {
                0320   pushNArray((char *)x,(unsigned int)(n*8), 1) ;
                0321 }
                0322 
                0323 void popInteger8Array(long *x, int n) {
                0324   popNArray((char *)x,(unsigned int)(n*8), 1) ;
                0325 }
                0326 
                0327 void pushReal4Array(float *x, int n) {
                0328   pushNArray((char *)x,(unsigned int)(n*4), 1) ;
                0329 }
                0330 
                0331 void popReal4Array(float *x, int n) {
                0332   popNArray((char *)x,(unsigned int)(n*4), 1) ;
                0333 }
                0334 
                0335 void pushReal8Array(double *x, int n) {
                0336   pushNArray((char *)x,(unsigned int)(n*8), 1) ;
                0337 }
                0338 
                0339 void popReal8Array(double *x, int n) {
                0340   popNArray((char *)x,(unsigned int)(n*8), 1) ;
                0341 }
                0342 
                0343 void pushComplex8Array(ccmplx *x, int n) {
                0344   pushNArray((char *)x,(unsigned int)(n*8), 1) ;
                0345 }
                0346 
                0347 void popComplex8Array(ccmplx *x, int n) {
                0348   popNArray((char *)x,(unsigned int)(n*8), 1) ;
                0349 }
                0350 
                0351 void pushComplex16Array(cdcmplx *x, int n) {
                0352   pushNArray((char *)x,(unsigned int)(n*16), 1) ;
                0353 }
                0354 
                0355 void popComplex16Array(cdcmplx *x, int n) {
                0356   popNArray((char *)x,(unsigned int)(n*16), 1) ;
                0357 }
                0358 
                0359 void pushCharacterArray(char *x, int n) {
                0360   pushNArray(x,(unsigned int)n, 1) ;
                0361 }
                0362 
                0363 void popCharacterArray(char *x, int n) {
                0364   popNArray(x,(unsigned int)n, 1) ;
                0365 }
                0366 
                0367 /* ********* Useful only for testpushpop.f90. Should go away! ********* */
                0368 
                0369 void showpushpopsequence_(int *op, int *index, int* nbobjects, int* sorts, int* sizes) {
                0370   char *prefix = "" ;
                0371   if (*op==1) prefix = "+" ;
                0372   else if (*op==-1) prefix = "-" ;
                0373   else if (*op==2) prefix = "+s" ;
                0374   else if (*op==-2) prefix = "-s" ;
                0375   else if (*op==-3) prefix = "Ls" ;
                0376   printf("%s%02i", prefix, *index) ;
                0377   // Comment the rest for compact display:
                0378   printf(":") ;
                0379   int i ;
                0380   for (i=0 ; i<*nbobjects ; ++i) {
                0381     switch (sorts[i]) {
                0382     case 1:
                0383       printf(" I4") ;
                0384       break ;
                0385     case 2:
                0386       printf(" I8") ;
                0387       break ;
                0388     case 3:
                0389       printf(" R4") ;
                0390       break ;
                0391     case 4:
                0392       printf(" R8") ;
                0393       break ;
                0394     case 5:
                0395       printf(" C8") ;
                0396       break ;
                0397     case 6:
                0398       printf(" C16") ;
                0399       break ;
                0400     case 7:
                0401       printf(" char") ;
                0402       break ;
                0403     case 8:
                0404       printf(" bit") ;
                0405       break ;
                0406     case 9:
                0407       printf(" PTR") ;
                0408       break ;
                0409     }
                0410     if (sizes[i]!=0) printf("[%1i]",sizes[i]) ;
                0411   }
                0412 }
                0413 
                0414 /****************** INTERFACE CALLED FROM FORTRAN *******************/
                0415 
                0416 void showstack_() {
                0417   showStack() ;
                0418 }
                0419 
                0420 void showstacksize_(int *i4i, int *i8i, int *r4i, int *r8i, int *c8i, int *c16i, int *s1i, int *biti, int *ptri, int *pos) {
                0421   showStackSize(*i4i,*i8i,*r4i,*r8i,*c8i,*c16i,*s1i,*biti,*ptri, *pos) ;
                0422 }
                0423 
                0424 void adstack_showpeaksize_() {
                0425   adStack_showPeakSize() ;
                0426 }
                0427 
                0428 void adstack_showpeaksize__() {
                0429   adStack_showPeakSize() ;
                0430 }
                0431 
                0432 void showtotaltraffic_(unsigned long long int *traffic) {
                0433   showTotalTraffic(*traffic) ;
                0434 }
                0435 
                0436 void startstackrepeat1_() {
                0437   startStackRepeat1() ;
                0438 }
                0439 
                0440 void startstackrepeat2_() {
                0441   startStackRepeat2() ;
                0442 }
                0443 
                0444 void resetstackrepeat1_() {
                0445   resetStackRepeat1() ;
                0446 }
                0447 
                0448 void resetstackrepeat2_() {
                0449   resetStackRepeat2() ;
                0450 }
                0451 
                0452 void endstackrepeat_() {
                0453   endStackRepeat() ;
                0454 }
                0455 
                0456 void pushnarray_(char *x, unsigned int *nbChars, int *checkReadOnly) {
                0457   pushNArray(x, *nbChars, *checkReadOnly) ;
                0458 }
                0459 
                0460 void popnarray_(char *x, unsigned int *nbChars, int *checkReadOnly) {
                0461   popNArray(x, *nbChars, *checkReadOnly) ;
                0462 }
                0463 
                0464 void pushinteger4array_(int *ii, int *ll) {
                0465   pushInteger4Array(ii, *ll) ;
                0466 }
                0467 
                0468 void popinteger4array_(int *ii, int *ll) {
                0469   popInteger4Array(ii, *ll) ;
                0470 }
                0471 
                0472 void pushinteger8array_(long *ii, int *ll) {
                0473   pushInteger8Array(ii, *ll) ;
                0474 }
                0475 
                0476 void popinteger8array_(long *ii, int *ll) {
                0477   popInteger8Array(ii, *ll) ;
                0478 }
                0479 
                0480 void pushreal4array_(float *ii, int *ll) {
                0481   pushReal4Array(ii, *ll) ;
                0482 }
                0483 
                0484 void popreal4array_(float *ii, int *ll) {
                0485   popReal4Array(ii, *ll) ;
                0486 }
                0487 
                0488 void pushreal8array_(double *ii, int *ll) {
                0489   pushReal8Array(ii, *ll) ;
                0490 }
                0491 
                0492 void popreal8array_(double *ii, int *ll) {
                0493   popReal8Array(ii, *ll) ;
                0494 }
                0495 
                0496 void pushcomplex8array_(ccmplx *ii, int *ll) {
                0497   pushComplex8Array(ii, *ll) ;
                0498 }
                0499 
                0500 void popcomplex8array_(ccmplx *ii, int *ll) {
                0501   popComplex8Array(ii, *ll) ;
                0502 }
                0503 
                0504 void pushcomplex16array_(cdcmplx *ii, int *ll) {
                0505   pushComplex16Array(ii, *ll) ;
                0506 }
                0507 
                0508 void popcomplex16array_(cdcmplx *ii, int *ll) {
                0509   popComplex16Array(ii, *ll) ;
                0510 }
                0511 
                0512 void pushcharacterarray_(char *ii, int *ll) {
                0513   pushCharacterArray(ii, *ll) ;
                0514 }
                0515 
                0516 void popcharacterarray_(char *ii, int *ll) {
                0517   popCharacterArray(ii, *ll) ;
                0518 }
                0519 
                0520 void pushbooleanarray_(char *x, unsigned int *n) {
                0521   pushNArray(x,(*n*4), 1) ;
                0522 }
                0523 
                0524 void popbooleanarray_(char *x, unsigned int *n) {
                0525   popNArray(x,(*n*4), 1) ;
                0526 }