TAP_support/ADFirstAidKit/adContextCPX.c

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 <string.h>
                0009 #include <stdio.h>
                0010 #include <stdlib.h>
                0011 #include "adContextCPX.h"
                0012
                0013 static int dbad_mode, dbad_phase ;
                0014 static double dbad_ddeps = 1.e-6 ;
                0015 static double dbad_seed = 0.137 ;
                0016 static double dbad_currentSeed = 0.0 ;
                0017 static double dbad_condensed_val, dbad_condensed_tgt ;
                0018
                0019 double dbad_nextRandom() {
                0020   dbad_currentSeed += dbad_seed ;
                0021   if (dbad_currentSeed>=1.0) dbad_currentSeed-=1.0 ;
                0022   /* Return a value in range [1.0 2.0[ */
                0023   return dbad_currentSeed+1.0 ;
                0024 }
                0025
                0026 void adContextCpx_init(double epsilon, double seed) {
                0027   dbad_mode = 1 ;
                0028   dbad_ddeps = epsilon ;
                0029   dbad_seed = seed ;
                0030   char* phase = getenv("DBAD_PHASE") ;
                0031   if (phase==NULL) {
                0032     printf("Please set DBAD_PHASE environment variable to 1 (perturbed) or 2 (tangent)\n") ;
                0033     exit(0) ;
                0034   } else if (strcmp(phase,"2")==0) {
                0035     printf("Tangent code,  seed=%7.1e\n", seed) ;
                0036     printf("=============================================\n") ;
                0037     dbad_phase = 2 ;
                0038     dbad_currentSeed = 0.0 ;
                0039   } else if (strcmp(phase,"1")==0) {
                0040     printf("Perturbed run, seed=%7.1e, epsilon=%7.1e\n", seed, epsilon) ;
                0041     printf("=============================================\n") ;
                0042     dbad_phase = 1 ;
                0043     dbad_currentSeed = 0.0 ;
                0044   } else if (strcmp(phase,"99")==0) {
                0045     printf("INTERNAL INTERFACE TESTS, seed=%7.1e, epsilon=%7.1e\n", seed, epsilon) ;
                0046     printf("=============================================\n") ;
                0047     dbad_phase = 99 ;
                0048   } else {
                0049     printf("DBAD_PHASE environment variable must be set to 1 or 2\n") ;
                0050     exit(0) ;
                0051   }
                0052 }
                0053
                0054 void adContextCpx_initReal8(char* varname, cdcmplx *indep) {
                0055   indep->di = dbad_ddeps*dbad_nextRandom() ;
                0056   if (dbad_phase==1)
                0057     indep->dr = indep->dr + indep->di ;
                0058   else if (dbad_phase==99)
                0059     printf("initReal8 of %s: %24.16e + i*%24.16e\n", varname, indep->dr, indep->di) ;
                0060 }
                0061
                0062 void adContextCpx_initReal8Array(char* varname, cdcmplx *indep, int length) {
                0063   int i ;
                0064   for (i=0 ; i<length ; ++i) {
                0065     indep[i].di = dbad_ddeps*dbad_nextRandom() ;
                0066   }
                0067   if (dbad_phase==1) {
                0068     for (i=0 ; i<length ; ++i) {
                0069       indep[i].dr = indep[i].dr + indep[i].di ;
                0070     }
                0071   } else if (dbad_phase==99) {
                0072     printf("initReal8Array of %s, length=%i:\n", varname, length) ;
                0073     for (i=0 ; i<length ; ++i)
                0074       printf("    %i:%24.16e + i*%24.16e",i,indep[i].dr,indep[i].di) ;
                0075     printf("\n") ;
                0076   }
                0077 }
                0078
                0079 void adContextCpx_initReal4(char* varname, ccmplx *indep) {
                0080   indep->i = dbad_ddeps*(float)dbad_nextRandom() ;
                0081   if (dbad_phase==1)
                0082     indep->r = indep->r + indep->i ;
                0083   else if (dbad_phase==99)
                0084     printf("initReal4 of %s: %24.16e + i*%24.16e\n", varname, indep->r, indep->i) ;
                0085 }
                0086
                0087 void adContextCpx_initReal4Array(char* varname, ccmplx *indep, int length) {
                0088   int i ;
                0089   for (i=0 ; i<length ; ++i) {
                0090     indep[i].i = dbad_ddeps*(float)dbad_nextRandom() ;
                0091   }
                0092   if (dbad_phase==1) {
                0093     for (i=0 ; i<length ; ++i) {
                0094       indep[i].r = indep[i].r + indep[i].i ;
                0095     }
                0096   } else if (dbad_phase==99) {
                0097     printf("initReal4Array of %s, length=%i:\n", varname, length) ;
                0098     for (i=0 ; i<length ; ++i)
                0099       printf("    %i:%24.16e + i*%24.16e",i,indep[i].r,indep[i].i) ;
                0100     printf("\n") ;
                0101   }
                0102 }
                0103
                0104 // NO adContextCpx_initComplex* PRIMITIVES, BECAUSE COMPLEX-STEP MODE REQUIRES NO COMPLEX IN PRIMAL CODE
                0105
                0106 void adContextCpx_startConclude() {
                0107   dbad_currentSeed= 0.0 ;
                0108   dbad_condensed_val = 0.0 ;
                0109   dbad_condensed_tgt = 0.0 ;
                0110 }
                0111
                0112 void adContextCpx_concludeReal8(char* varname, cdcmplx *dep) {
                0113   double depb = dbad_nextRandom() ;
                0114   dbad_condensed_val += depb*(dep->dr) ;
                0115   if (dbad_phase==2 || dbad_phase==1)
                0116     dbad_condensed_tgt += depb*(dep->di)/dbad_ddeps ;
                0117   else if (dbad_phase==99)
                0118     printf("concludeReal8 of %s [%24.16e *] %24.16e + i*%24.16e\n", varname, depb, dep->dr, dep->di) ;
                0119 }
                0120
                0121 void adContextCpx_concludeReal8Array(char* varname, cdcmplx *dep, int length) {
                0122   int i ;
                0123   double depb ;
                0124   if (dbad_phase==99) printf("concludeReal8Array of %s, length=%i:\n", varname, length) ;
                0125   for (i=0 ; i<length ; ++i) {
                0126     depb = dbad_nextRandom() ;
                0127     dbad_condensed_val += depb*dep[i].dr ;
                0128     if (dbad_phase==2 || dbad_phase==1) {
                0129        dbad_condensed_tgt += depb*dep[i].di/dbad_ddeps ;
                0130     } else if (dbad_phase==99) {
                0131       printf("    %i:[%24.16e *] %24.16e + i*%24.16e",i,depb,dep[i].dr,dep[i].di) ;
                0132     }
                0133   }
                0134   if (dbad_phase==99) printf("\n") ;
                0135 }
                0136
                0137 void adContextCpx_concludeReal4(char* varname, ccmplx *dep) {
                0138   float depb = (float)dbad_nextRandom() ;
                0139   dbad_condensed_val += depb*(dep->r) ;
                0140   if (dbad_phase==2 || dbad_phase==1)
                0141     dbad_condensed_tgt += depb*(dep->i)/dbad_ddeps ;
                0142   else if (dbad_phase==99)
                0143     printf("concludeReal4 of %s [%24.16e *] %24.16e + i*%24.16e\n", varname, depb, dep->r, dep->i) ;
                0144 }
                0145
                0146 void adContextCpx_concludeReal4Array(char* varname, ccmplx *dep, int length) {
                0147   int i ;
                0148   float depb ;
                0149   if (dbad_phase==99) printf("concludeReal4Array of %s, length=%i:\n", varname, length) ;
                0150   for (i=0 ; i<length ; ++i) {
                0151     depb = (float)dbad_nextRandom() ;
                0152     dbad_condensed_val += depb*dep[i].r ;
                0153     if (dbad_phase==2 || dbad_phase==1) {
                0154        dbad_condensed_tgt += depb*dep[i].i/dbad_ddeps ;
                0155     } else if (dbad_phase==99) {
                0156       printf("    %i:[%24.16e *] %24.16e + i*%24.16e",i,depb,dep[i].r,dep[i].i) ;
                0157     }
                0158   }
                0159   if (dbad_phase==99) printf("\n") ;
                0160 }
                0161
                0162 // NO adContextCpx_concludeComplex* PRIMITIVES, BECAUSE COMPLEX-STEP MODE REQUIRES NO COMPLEX IN PRIMAL CODE
                0163
                0164 void adContextCpx_conclude() {
                0165   if (dbad_phase==2) {
                0166     printf("[seed:%7.1e] Condensed result : %24.16e\n", dbad_seed, dbad_condensed_val) ;
                0167     printf("[seed:%7.1e] Condensed tangent: %24.16e\n", dbad_seed, dbad_condensed_tgt) ;
                0168   } else if (dbad_phase==1) {
                0169     printf("[seed:%7.1e] Condensed perturbed result : %24.16e (epsilon:%7.1e)\n",
                0170            dbad_seed, dbad_condensed_val, dbad_ddeps) ;
                0171     printf("[seed:%7.1e] Condensed perturbed tangent: %24.16e\n", dbad_seed, dbad_condensed_tgt) ;
                0172   }
                0173 }
                0174
                0175 // NO adContextAdj_* PRIMITIVES, BECAUSE COMPLEX-STEP MODE WORKS ONLY IN TANGENT MODE!
                0176
                0177 //############## INTERFACE PROCEDURES CALLED FROM FORTRAN ################
                0178
                0179 void adcontextcpx_init_(double *epsilon, double *seed) {
                0180   adContextCpx_init(*epsilon, *seed) ;
                0181 }
                0182
                0183 void adcontextcpx_initreal8_(char* varname, cdcmplx *indep) {
                0184   adContextCpx_initReal8(varname, indep) ;
                0185 }
                0186
                0187 void adcontextcpx_initreal8array_(char* varname, cdcmplx *indep, int *length) {
                0188   adContextCpx_initReal8Array(varname, indep, *length) ;
                0189 }
                0190
                0191 void adcontextcpx_initreal4_(char* varname, ccmplx *indep) {
                0192   adContextCpx_initReal4(varname, indep) ;
                0193 }
                0194
                0195 void adcontextcpx_initreal4array_(char* varname, ccmplx *indep, int *length) {
                0196   adContextCpx_initReal4Array(varname, indep, *length) ;
                0197 }
                0198
                0199 void adcontextcpx_startconclude_() {
                0200   adContextCpx_startConclude() ;
                0201 }
                0202
                0203 void adcontextcpx_concludereal8_(char* varname, cdcmplx *dep) {
                0204   if (dbad_phase==99)
                0205       printf("concludereal8_ of %s: \n", varname);
                0206   adContextCpx_concludeReal8(varname, dep) ;
                0207 }
                0208
                0209 void adcontextcpx_concludereal8array_(char* varname, cdcmplx *dep, int *length) {
                0210   if (dbad_phase==99)
                0211       printf("concludereal8array_ of %s: \n", varname);
                0212   adContextCpx_concludeReal8Array(varname, dep, *length) ;
                0213 }
                0214
                0215 void adcontextcpx_concludereal4_(char* varname, ccmplx *dep) {
                0216   adContextCpx_concludeReal4(varname, dep) ;
                0217 }
                0218
                0219 void adcontextcpx_concludereal4array_(char* varname, ccmplx *dep, int *length) {
                0220   adContextCpx_concludeReal4Array(varname, dep, *length) ;
                0221 }
                0222
                0223 void adcontextcpx_conclude_() {
                0224   adContextCpx_conclude() ;
                0225 }