ADMB Documentation  11.1.2495
 All Classes Files Functions Variables Typedefs Friends Defines
df1b2prb.cpp
Go to the documentation of this file.
00001 /*
00002  * $Id: df1b2prb.cpp 2351 2014-09-15 23:40:15Z johnoel $
00003  *
00004  * Author: David Fournier
00005  * Copyright (c) 2008-2012 Regents of the University of California
00006  */
00011 #include <df1b2fun.h>
00012 
00013 #ifndef OPT_LIB
00014   #include <cassert>
00015   #include <climits>
00016 #endif
00017 
00022   df1b2variable operator * (const df1b2variable& _x, double y)
00023   {
00024     ADUNCONST(df1b2variable,x)
00025     df1b2variable z;
00026     double * xd=x.get_u_dot();
00027     double * zd=z.get_u_dot();
00028     double xu=*x.get_u();
00029 
00030     *z.get_u()=y*xu;
00031 
00032     for (int i=0;i<df1b2variable::nvar;i++)
00033     {
00034       *zd++ = y * *xd++;
00035     }
00036 
00037     // WRITE WHATEVER ON TAPE
00038     if (!df1b2_gradlist::no_derivatives)
00039       f1b2gradlist->write_pass1_prod(&x,y,&z);
00040     return z;
00041   }
00042 
00043 void ad_read_pass2_prodc2(void);
00044 
00049  int df1b2_gradlist::write_pass1_prod(const df1b2variable * _px,double y,
00050    df1b2variable * pz)
00051  {
00052    ADUNCONST(df1b2variable*,px)
00053    ncount++;
00054 #if defined(CHECK_COUNT)
00055   if (ncount >= ncount_check)
00056     cout << ncount << endl;
00057 #endif
00058    int nvar=df1b2variable::nvar;
00059 
00060    //int total_bytes=3*sizeof(df1b2_header)+sizeof(char*)
00061    //  +2*(nvar+1)*sizeof(double);
00062    size_t total_bytes=2*sizeof(df1b2_header)
00063      +(nvar+2)*sizeof(double);
00064 // string identifier debug stuff
00065 #if defined(SAFE_ALL)
00066   char ids[]="DL";
00067   size_t slen=strlen(ids);
00068   total_bytes+=slen;
00069 #endif
00070 
00071 #ifndef OPT_LIB
00072   assert(total_bytes <= INT_MAX);
00073 #endif
00074   list.check_buffer_size((int)total_bytes);
00075   void * tmpptr=list.bptr;
00076 #if defined(SAFE_ALL)
00077   memcpy(list,ids,slen);
00078 #endif
00079 // end of string identifier debug stuff
00080 
00081    memcpy(list,(df1b2_header*)(px),sizeof(df1b2_header));
00082    memcpy(list,&y,sizeof(double));
00083    memcpy(list,(df1b2_header*)(pz),sizeof(df1b2_header));
00084   const int sizeofdouble = sizeof(double);
00085    memcpy(list,px->get_u(),sizeofdouble);
00086    memcpy(list,px->get_u_dot(),nvar*sizeofdouble);
00087    // ***** write  record size
00088    nlist.bptr->numbytes=adptr_diff(list.bptr,tmpptr);
00089    if ((int)total_bytes != nlist.bptr->numbytes)
00090    {
00091      cerr << "error in byte calculation in "
00092        " write_pass1_prod" << endl;
00093      ad_exit(1);
00094    }
00095    nlist.bptr->pf=(ADrfptr)(&ad_read_pass2_prodc2);
00096       ++nlist;
00097    return 0;
00098  }
00099 
00100 
00101 void read_pass2_1_prodc2(void);
00102 void read_pass2_2_prodc2(void);
00103 void read_pass2_3_prodc2(void);
00104 
00109 void ad_read_pass2_prodc2(void)
00110 {
00111   switch(df1b2variable::passnumber)
00112   {
00113   case 1:
00114     read_pass2_1_prodc2();
00115     break;
00116   case 2:
00117     read_pass2_2_prodc2();
00118     break;
00119   case 3:
00120     read_pass2_3_prodc2();
00121     break;
00122   default:
00123     cerr << "illegal value for df1b2variable::pass = "
00124          << df1b2variable::passnumber << endl;
00125     exit(1);
00126   }
00127 }
00128 
00133 void read_pass2_1_prodc2(void)
00134 {
00135   // We are going backword for bptr and nbptr
00136   // and  forward for bptr2 and nbptr2
00137   // the current entry+2 in bptr is the size of the record i.e
00138   // points to the next record
00139   //char * bptr=f1b2gradlist->bptr;
00140   //char * bptr2=f1b2gradlist2->bptr;
00141   int nvar=df1b2variable::nvar;
00142   test_smartlist& list=f1b2gradlist->list;
00143   //f1b2gradlist->nlist-=sizeof(int);
00144   int num_bytes=f1b2gradlist->nlist.bptr->numbytes;
00145   list-=num_bytes;
00146   list.saveposition(); // save pointer to beginning of record;
00147   double xu;
00148 
00149   // get info from tape1
00150 #if defined(SAFE_ALL)
00151   checkidentiferstring("DL",f1b2gradlist->list);
00152 #endif
00153   char * bptr=f1b2gradlist->list.bptr;
00154   df1b2_header * px=(df1b2_header *) bptr;
00155   bptr+=sizeof(df1b2_header);
00156   double yu=*(double *) bptr;
00157   bptr+=sizeof(double);
00158   df1b2_header * pz=(df1b2_header *) bptr;
00159   bptr+=sizeof(df1b2_header);
00160   memcpy(&xu,bptr,sizeof(double));
00161   //bptr+=sizeof(double); double* xdot=(double*)bptr;
00162 
00163   list.restoreposition(); // save pointer to beginning of record;
00164 
00165   // ****************************************************************
00166   // turn this off if no third derivatives are calculated
00167   // if (!no_third_derivatives)
00168   // {
00169   // save for second reverse pass
00170   // save identifier 1
00171      test_smartlist & list2 = f1b2gradlist->list2;
00172 
00173   size_t total_bytes=2*nvar*sizeof(double);
00174 // string identifier debug stuff
00175 #if defined(SAFE_ALL)
00176   char ids[]="QK";
00177   size_t slen=strlen(ids);
00178   total_bytes+=slen;
00179 #endif
00180 
00181 #ifndef OPT_LIB
00182   assert(total_bytes <= INT_MAX);
00183 #endif
00184 
00185   list2.check_buffer_size((int)total_bytes);
00186   void * tmpptr=list2.bptr;
00187 #if defined(SAFE_ALL)
00188   memcpy(list2,ids,slen);
00189 #endif
00190 
00191   fixed_smartlist2 & nlist2 = f1b2gradlist->nlist2;
00192   const int sizeofdouble = sizeof(double);
00193   memcpy(list2,pz->get_u_bar(),nvar*sizeofdouble);
00194   memcpy(list2,pz->get_u_dot_bar(),nvar*sizeofdouble);
00195   *nlist2.bptr=adptr_diff(list2.bptr,tmpptr);
00196   ++nlist2;
00197 
00198   // Do first reverse pass calculations
00199   for (int i=0;i<nvar;i++)
00200   {
00201     px->u_bar[i]+=yu*pz->u_bar[i];
00202   }
00203 
00204   for (int i=0;i<nvar;i++)
00205   {
00206     px->u_dot_bar[i]+=yu*pz->u_dot_bar[i];
00207   }
00208 
00209   // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00210   for (int i=0;i<nvar;i++)
00211   {
00212     pz->u_bar[i]=0;
00213   }
00214   for (int i=0;i<nvar;i++)
00215   {
00216     pz->u_dot_bar[i]=0;
00217   }
00218 }
00219 
00224 void read_pass2_2_prodc2(void)
00225 {
00226   //const int nlist_record_size=sizeof(int)+sizeof(char*);
00227   // We are going forward for bptr and backword for bptr2
00228   //
00229   // list 1
00230   //
00231   int nvar=df1b2variable::nvar;
00232   test_smartlist & list=f1b2gradlist->list;
00233 
00234   size_t total_bytes=2*sizeof(df1b2_header)
00235      +(nvar+2)*sizeof(double);
00236 // string identifier debug stuff
00237 #if defined(SAFE_ALL)
00238   char ids[]="DL";
00239   size_t slen=strlen(ids);
00240   total_bytes+=slen;
00241 #endif
00242 
00243 #ifndef OPT_LIB
00244   assert(total_bytes <= INT_MAX);
00245 #endif
00246   list.check_buffer_size((int)total_bytes);
00247 // end of string identifier debug stuff
00248 
00249   list.saveposition(); // save pointer to beginning of record;
00250   fixed_smartlist & nlist=f1b2gradlist->nlist;
00251    // nlist-=sizeof(int);
00252   // get record size
00253   int num_bytes=nlist.bptr->numbytes;
00254   //
00255   // list 2
00256   //
00257   test_smartlist & list2=f1b2gradlist->list2;
00258   fixed_smartlist2 & nlist2=f1b2gradlist->nlist2;
00259   // get record size
00260   int num_bytes2=*nlist2.bptr;
00261   --nlist2;
00262   // backup the size of the record
00263   list2-=num_bytes2;
00264   list2.saveposition(); // save pointer to beginning of record;
00265   // save the pointer to the beginning of the record
00266   // bptr and bptr2 now both point to the beginning of their records
00267 
00268   double xu;
00269   //df1b2_header x,z;
00270   //df1b2function2 * pf;
00271 
00272   // get info from tape1
00273   // get info from tape1
00274 #if defined(SAFE_ALL)
00275   checkidentiferstring("DL",list);
00276   checkidentiferstring("QK",list2);
00277 #endif
00278   df1b2_header * px=(df1b2_header *) list.bptr;
00279   list.bptr+=sizeof(df1b2_header);
00280   double yu=*(double *) list.bptr;
00281   list.bptr+=sizeof(double);
00282   df1b2_header * pz=(df1b2_header *) list.bptr;
00283   list.bptr+=sizeof(df1b2_header);
00284   memcpy(&xu,list.bptr,sizeof(double));
00285   list.bptr+=sizeof(double);
00286   //double* xdot=(double*)list.bptr;
00287   list.restoreposition(num_bytes); // save pointer to beginning of record;
00288 
00289   //double* zbar=(double*)list2.bptr;
00290   //double* zdotbar=(double*)(list2.bptr+nvar*sizeof(double));
00291 
00292   list2.restoreposition(); // save pointer to beginning of record;
00293 
00294   //double * x_tilde=px->get_u_tilde();
00295   //double * x_dot_tilde=px->get_u_dot_tilde();
00296   double * x_bar_tilde=px->get_u_bar_tilde();
00297   double * x_dot_bar_tilde=px->get_u_dot_bar_tilde();
00298   double * z_bar_tilde=pz->get_u_bar_tilde();
00299   double * z_dot_bar_tilde=pz->get_u_dot_bar_tilde();
00300   // Do second "reverse-reverse" pass calculations
00301 
00302   int i;
00303 
00304   for (i=0;i<nvar;i++)
00305   {
00306     z_bar_tilde[i]=0;
00307     z_dot_bar_tilde[i]=0;
00308   }
00309 
00310   // start with y and add x
00311   for (i=0;i<nvar;i++)
00312   {
00313     z_bar_tilde[i]+=yu*x_bar_tilde[i];
00314   }
00315 
00316   for (i=0;i<nvar;i++)
00317   {
00318     z_dot_bar_tilde[i]+=yu*x_dot_bar_tilde[i];
00319   }
00320 }
00321 
00326 void read_pass2_3_prodc2(void)
00327 {
00328   // We are going backword for bptr and forward for bptr2
00329   // the current entry+2 in bptr is the size of the record i.e
00330   // points to the next record
00331   int nvar=df1b2variable::nvar;
00332   fixed_smartlist & nlist=f1b2gradlist->nlist;
00333   test_smartlist& list=f1b2gradlist->list;
00334    // nlist-=sizeof(int);
00335   // get record size
00336   int num_bytes=nlist.bptr->numbytes;
00337   // backup the size of the record
00338   list-=num_bytes;
00339   list.saveposition(); // save pointer to beginning of record;
00340   // save the pointer to the beginning of the record
00341   double xu;
00342   //df1b2_header z;
00343   //df1b2function2 * pf;
00344 
00345   // get info from tape1
00346   // get info from tape1
00347 #if defined(SAFE_ALL)
00348   checkidentiferstring("DL",list);
00349 #endif
00350   df1b2_header * px=(df1b2_header *) list.bptr;
00351   list.bptr+=sizeof(df1b2_header);
00352   double yu=*(double*) list.bptr;
00353   list.bptr+=sizeof(double);
00354   df1b2_header * pz=(df1b2_header *) list.bptr;
00355   list.bptr+=sizeof(df1b2_header);
00356   memcpy(&xu,list.bptr,sizeof(double));
00357   list.bptr+=sizeof(double);
00358   //double* xdot=(double*)list.bptr;
00359   list.restoreposition(); // save pointer to beginning of record;
00360   int i;
00361 
00362   *(px->u_tilde)+=yu* *(pz->u_tilde);
00363   for (i=0;i<nvar;i++)
00364   {
00365     px->u_dot_tilde[i]+=yu*pz->u_dot_tilde[i];
00366   }
00367   *(pz->u_tilde)=0;
00368   for (i=0;i<nvar;i++)
00369   {
00370     pz->u_dot_tilde[i]=0;
00371   }
00372 }