ADMB Documentation  11.1x.2711
 All Classes Files Functions Variables Typedefs Friends Defines
df1b2prb.cpp
Go to the documentation of this file.
00001 /*
00002  * $Id: df1b2prb.cpp 2601 2014-11-08 22:34:58Z 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 (unsigned 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    unsigned 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   list.check_buffer_size(total_bytes);
00072   void * tmpptr=list.bptr;
00073 #if defined(SAFE_ALL)
00074   memcpy(list,ids,slen);
00075 #endif
00076 // end of string identifier debug stuff
00077 
00078    memcpy(list,(df1b2_header*)(px),sizeof(df1b2_header));
00079    memcpy(list,&y,sizeof(double));
00080    memcpy(list,(df1b2_header*)(pz),sizeof(df1b2_header));
00081   const int sizeofdouble = sizeof(double);
00082    memcpy(list,px->get_u(),sizeofdouble);
00083    memcpy(list,px->get_u_dot(),nvar*sizeofdouble);
00084    // ***** write  record size
00085    nlist.bptr->numbytes=adptr_diff(list.bptr,tmpptr);
00086    if ((int)total_bytes != nlist.bptr->numbytes)
00087    {
00088      cerr << "error in byte calculation in "
00089        " write_pass1_prod" << endl;
00090      ad_exit(1);
00091    }
00092    nlist.bptr->pf=(ADrfptr)(&ad_read_pass2_prodc2);
00093       ++nlist;
00094    return 0;
00095  }
00096 
00097 
00098 void read_pass2_1_prodc2(void);
00099 void read_pass2_2_prodc2(void);
00100 void read_pass2_3_prodc2(void);
00101 
00106 void ad_read_pass2_prodc2(void)
00107 {
00108   switch(df1b2variable::passnumber)
00109   {
00110   case 1:
00111     read_pass2_1_prodc2();
00112     break;
00113   case 2:
00114     read_pass2_2_prodc2();
00115     break;
00116   case 3:
00117     read_pass2_3_prodc2();
00118     break;
00119   default:
00120     cerr << "illegal value for df1b2variable::pass = "
00121          << df1b2variable::passnumber << endl;
00122     exit(1);
00123   }
00124 }
00125 
00130 void read_pass2_1_prodc2(void)
00131 {
00132   // We are going backword for bptr and nbptr
00133   // and  forward for bptr2 and nbptr2
00134   // the current entry+2 in bptr is the size of the record i.e
00135   // points to the next record
00136   //char * bptr=f1b2gradlist->bptr;
00137   //char * bptr2=f1b2gradlist2->bptr;
00138   unsigned int nvar=df1b2variable::nvar;
00139   test_smartlist& list=f1b2gradlist->list;
00140   //f1b2gradlist->nlist-=sizeof(int);
00141   int num_bytes=f1b2gradlist->nlist.bptr->numbytes;
00142   list-=num_bytes;
00143   list.saveposition(); // save pointer to beginning of record;
00144   double xu;
00145 
00146   // get info from tape1
00147 #if defined(SAFE_ALL)
00148   checkidentiferstring("DL",f1b2gradlist->list);
00149 #endif
00150   char * bptr=f1b2gradlist->list.bptr;
00151   df1b2_header * px=(df1b2_header *) bptr;
00152   bptr+=sizeof(df1b2_header);
00153   double yu=*(double *) bptr;
00154   bptr+=sizeof(double);
00155   df1b2_header * pz=(df1b2_header *) bptr;
00156   bptr+=sizeof(df1b2_header);
00157   memcpy(&xu,bptr,sizeof(double));
00158   //bptr+=sizeof(double); double* xdot=(double*)bptr;
00159 
00160   list.restoreposition(); // save pointer to beginning of record;
00161 
00162   // ****************************************************************
00163   // turn this off if no third derivatives are calculated
00164   // if (!no_third_derivatives)
00165   // {
00166   // save for second reverse pass
00167   // save identifier 1
00168      test_smartlist & list2 = f1b2gradlist->list2;
00169 
00170   size_t total_bytes=2*nvar*sizeof(double);
00171 // string identifier debug stuff
00172 #if defined(SAFE_ALL)
00173   char ids[]="QK";
00174   size_t slen=strlen(ids);
00175   total_bytes+=slen;
00176 #endif
00177 
00178   list2.check_buffer_size(total_bytes);
00179   void * tmpptr=list2.bptr;
00180 #if defined(SAFE_ALL)
00181   memcpy(list2,ids,slen);
00182 #endif
00183 
00184   fixed_smartlist2 & nlist2 = f1b2gradlist->nlist2;
00185   const int sizeofdouble = sizeof(double);
00186   memcpy(list2,pz->get_u_bar(),nvar*sizeofdouble);
00187   memcpy(list2,pz->get_u_dot_bar(),nvar*sizeofdouble);
00188   *nlist2.bptr=adptr_diff(list2.bptr,tmpptr);
00189   ++nlist2;
00190 
00191   // Do first reverse pass calculations
00192   for (size_t i=0;i<nvar;i++)
00193   {
00194     px->u_bar[i]+=yu*pz->u_bar[i];
00195   }
00196 
00197   for (size_t i=0;i<nvar;i++)
00198   {
00199     px->u_dot_bar[i]+=yu*pz->u_dot_bar[i];
00200   }
00201 
00202   // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00203   for (size_t i=0;i<nvar;i++)
00204   {
00205     pz->u_bar[i]=0;
00206   }
00207   for (size_t i=0;i<nvar;i++)
00208   {
00209     pz->u_dot_bar[i]=0;
00210   }
00211 }
00212 
00217 void read_pass2_2_prodc2(void)
00218 {
00219   //const int nlist_record_size=sizeof(int)+sizeof(char*);
00220   // We are going forward for bptr and backword for bptr2
00221   //
00222   // list 1
00223   //
00224   unsigned int nvar=df1b2variable::nvar;
00225   test_smartlist & list=f1b2gradlist->list;
00226 
00227   size_t total_bytes=2*sizeof(df1b2_header)
00228      +(nvar+2)*sizeof(double);
00229 // string identifier debug stuff
00230 #if defined(SAFE_ALL)
00231   char ids[]="DL";
00232   size_t slen=strlen(ids);
00233   total_bytes+=slen;
00234 #endif
00235 
00236   list.check_buffer_size(total_bytes);
00237 // end of string identifier debug stuff
00238 
00239   list.saveposition(); // save pointer to beginning of record;
00240   fixed_smartlist & nlist=f1b2gradlist->nlist;
00241    // nlist-=sizeof(int);
00242   // get record size
00243   int num_bytes=nlist.bptr->numbytes;
00244   //
00245   // list 2
00246   //
00247   test_smartlist & list2=f1b2gradlist->list2;
00248   fixed_smartlist2 & nlist2=f1b2gradlist->nlist2;
00249   // get record size
00250   int num_bytes2=*nlist2.bptr;
00251   --nlist2;
00252   // backup the size of the record
00253   list2-=num_bytes2;
00254   list2.saveposition(); // save pointer to beginning of record;
00255   // save the pointer to the beginning of the record
00256   // bptr and bptr2 now both point to the beginning of their records
00257 
00258   double xu;
00259   //df1b2_header x,z;
00260   //df1b2function2 * pf;
00261 
00262   // get info from tape1
00263   // get info from tape1
00264 #if defined(SAFE_ALL)
00265   checkidentiferstring("DL",list);
00266   checkidentiferstring("QK",list2);
00267 #endif
00268   df1b2_header * px=(df1b2_header *) list.bptr;
00269   list.bptr+=sizeof(df1b2_header);
00270   double yu=*(double *) list.bptr;
00271   list.bptr+=sizeof(double);
00272   df1b2_header * pz=(df1b2_header *) list.bptr;
00273   list.bptr+=sizeof(df1b2_header);
00274   memcpy(&xu,list.bptr,sizeof(double));
00275   list.bptr+=sizeof(double);
00276   //double* xdot=(double*)list.bptr;
00277   list.restoreposition(num_bytes); // save pointer to beginning of record;
00278 
00279   //double* zbar=(double*)list2.bptr;
00280   //double* zdotbar=(double*)(list2.bptr+nvar*sizeof(double));
00281 
00282   list2.restoreposition(); // save pointer to beginning of record;
00283 
00284   //double * x_tilde=px->get_u_tilde();
00285   //double * x_dot_tilde=px->get_u_dot_tilde();
00286   double * x_bar_tilde=px->get_u_bar_tilde();
00287   double * x_dot_bar_tilde=px->get_u_dot_bar_tilde();
00288   double * z_bar_tilde=pz->get_u_bar_tilde();
00289   double * z_dot_bar_tilde=pz->get_u_dot_bar_tilde();
00290   // Do second "reverse-reverse" pass calculations
00291 
00292   for (size_t i=0;i<nvar;i++)
00293   {
00294     z_bar_tilde[i]=0;
00295     z_dot_bar_tilde[i]=0;
00296   }
00297 
00298   // start with y and add x
00299   for (size_t i=0;i<nvar;i++)
00300   {
00301     z_bar_tilde[i]+=yu*x_bar_tilde[i];
00302   }
00303 
00304   for (size_t i=0;i<nvar;i++)
00305   {
00306     z_dot_bar_tilde[i]+=yu*x_dot_bar_tilde[i];
00307   }
00308 }
00309 
00314 void read_pass2_3_prodc2(void)
00315 {
00316   // We are going backword for bptr and forward for bptr2
00317   // the current entry+2 in bptr is the size of the record i.e
00318   // points to the next record
00319   unsigned int nvar=df1b2variable::nvar;
00320   fixed_smartlist & nlist=f1b2gradlist->nlist;
00321   test_smartlist& list=f1b2gradlist->list;
00322    // nlist-=sizeof(int);
00323   // get record size
00324   int num_bytes=nlist.bptr->numbytes;
00325   // backup the size of the record
00326   list-=num_bytes;
00327   list.saveposition(); // save pointer to beginning of record;
00328   // save the pointer to the beginning of the record
00329   double xu;
00330   //df1b2_header z;
00331   //df1b2function2 * pf;
00332 
00333   // get info from tape1
00334   // get info from tape1
00335 #if defined(SAFE_ALL)
00336   checkidentiferstring("DL",list);
00337 #endif
00338   df1b2_header * px=(df1b2_header *) list.bptr;
00339   list.bptr+=sizeof(df1b2_header);
00340   double yu=*(double*) list.bptr;
00341   list.bptr+=sizeof(double);
00342   df1b2_header * pz=(df1b2_header *) list.bptr;
00343   list.bptr+=sizeof(df1b2_header);
00344   memcpy(&xu,list.bptr,sizeof(double));
00345   list.bptr+=sizeof(double);
00346   //double* xdot=(double*)list.bptr;
00347   list.restoreposition(); // save pointer to beginning of record;
00348 
00349   *(px->u_tilde)+=yu* *(pz->u_tilde);
00350   for (size_t i=0;i<nvar;i++)
00351   {
00352     px->u_dot_tilde[i]+=yu*pz->u_dot_tilde[i];
00353   }
00354   *(pz->u_tilde)=0;
00355   for (size_t i=0;i<nvar;i++)
00356   {
00357     pz->u_dot_tilde[i]=0;
00358   }
00359 }