ADMB Documentation  11.2.2817
 All Classes Files Functions Variables Typedefs Friends Defines
df1b2prc.cpp
Go to the documentation of this file.
00001 /*
00002  * $Id: df1b2prc.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 #include <cassert>
00014 #include <climits>
00015 
00016 /*
00017 //int debugcounter=0;
00018   df1b2variable operator / (double x,const df1b2variable& y)
00019   {
00020     return x*inv(y);
00021   }
00022 */
00023 
00028   df1b2variable operator / (const df1b2variable& x,double y)
00029   {
00030     return x*(1.0/y);
00031   }
00032 
00037   df1b2variable operator * (double x,const df1b2variable& _y)
00038   {
00039     ADUNCONST(df1b2variable,y)
00040     df1b2variable z;
00041     double * yd=y.get_u_dot();
00042     double * zd=z.get_u_dot();
00043     double yu=*y.get_u();
00044 
00045     *z.get_u()=x*yu;
00046 
00047     for (unsigned int i=0;i<df1b2variable::nvar;i++)
00048     {
00049       *zd++ = x * *yd++;
00050     }
00051 
00052     // WRITE WHATEVER ON TAPE
00053     if (!df1b2_gradlist::no_derivatives)
00054       f1b2gradlist->write_pass1_prod(x,&y,&z);
00055     return z;
00056   }
00057 
00058 
00059 void ad_read_pass2_prodc1(void);
00060 
00065  int df1b2_gradlist::write_pass1_prod(double x,const df1b2variable * _py,
00066    df1b2variable * pz)
00067  {
00068    ADUNCONST(df1b2variable*,py)
00069    ncount++;
00070 #if defined(CHECK_COUNT)
00071   if (ncount >= ncount_check)
00072     cout << ncount << endl;
00073 #endif
00074    unsigned int nvar=df1b2variable::nvar;
00075 
00076    //int total_bytes=3*sizeof(df1b2_header)+sizeof(char*)
00077    //  +2*(nvar+1)*sizeof(double);
00078    size_t total_bytes=2*sizeof(df1b2_header)
00079      +(nvar+2)*sizeof(double);
00080 // string identifier debug stuff
00081 #if defined(SAFE_ALL)
00082   char ids[]="DL";
00083   size_t slen=strlen(ids);
00084   total_bytes+=slen;
00085 #endif
00086 
00087   list.check_buffer_size(total_bytes);
00088   void * tmpptr=list.bptr;
00089 #if defined(SAFE_ALL)
00090   memcpy(list,ids,slen);
00091 #endif
00092 // end of string identifier debug stuff
00093 
00094    memcpy(list,&x,sizeof(double));
00095    memcpy(list,(df1b2_header*)(py),sizeof(df1b2_header));
00096    memcpy(list,(df1b2_header*)(pz),sizeof(df1b2_header));
00097    const int sizeofdouble = sizeof(double);
00098    memcpy(list,py->get_u(),sizeofdouble);
00099    memcpy(list,py->get_u_dot(),nvar*sizeofdouble);
00100    // ***** write  record size
00101    nlist.bptr->numbytes=adptr_diff(list.bptr,tmpptr);
00102    if ((int)total_bytes != nlist.bptr->numbytes)
00103    {
00104      cerr << "error in byte calculation in "
00105        " write_pass1_prod" << endl;
00106      ad_exit(1);
00107    }
00108    nlist.bptr->pf=(ADrfptr)(&ad_read_pass2_prodc1);
00109       ++nlist;
00110    return 0;
00111  }
00112 
00113 
00114 void read_pass2_1_prodc1(void);
00115 void read_pass2_2_prodc1(void);
00116 void read_pass2_3_prodc1(void);
00117 
00122 void ad_read_pass2_prodc1(void)
00123 {
00124   switch(df1b2variable::passnumber)
00125   {
00126   case 1:
00127     read_pass2_1_prodc1();
00128     break;
00129   case 2:
00130     read_pass2_2_prodc1();
00131     break;
00132   case 3:
00133     read_pass2_3_prodc1();
00134     break;
00135   default:
00136     cerr << "illegal value for df1b2variable::pass = "
00137          << df1b2variable::passnumber << endl;
00138     exit(1);
00139   }
00140 }
00141 
00146 void read_pass2_1_prodc1(void)
00147 {
00148   //  vmon_begin();
00149   // We are going backword for bptr and nbptr
00150   // and  forward for bptr2 and nbptr2
00151   // the current entry+2 in bptr is the size of the record i.e
00152   // points to the next record
00153   //char * bptr=f1b2gradlist->bptr;
00154   //char * bptr2=f1b2gradlist2->bptr;
00155   unsigned int nvar=df1b2variable::nvar;
00156   test_smartlist& list=f1b2gradlist->list;
00157   //f1b2gradlist->nlist-=sizeof(int);
00158   int num_bytes=f1b2gradlist->nlist.bptr->numbytes;
00159   list-=num_bytes;
00160   list.saveposition(); // save pointer to beginning of record;
00161   double yu;
00162 
00163   // get info from tape1
00164 #if defined(SAFE_ALL)
00165   checkidentiferstring("DL",f1b2gradlist->list);
00166 #endif
00167   char * bptr=f1b2gradlist->list.bptr;
00168   double xu=*(double *) bptr;
00169   bptr+=sizeof(double);
00170   df1b2_header * py=(df1b2_header *) bptr;
00171   bptr+=sizeof(df1b2_header);
00172   df1b2_header * pz=(df1b2_header *) bptr;
00173   bptr+=sizeof(df1b2_header);
00174   memcpy(&yu,bptr,sizeof(double));
00175   //bptr+=sizeof(double); double* ydot=(double*)bptr;
00176 
00177   list.restoreposition(); // save pointer to beginning of record;
00178 
00179   // ****************************************************************
00180   // turn this off if no third derivatives are calculated
00181   // if (!no_third_derivatives)
00182   // {
00183   // save for second reverse pass
00184   // save identifier 1
00185      test_smartlist & list2 = f1b2gradlist->list2;
00186 
00187   size_t total_bytes=2*nvar*sizeof(double);
00188 // string identifier debug stuff
00189 #if defined(SAFE_ALL)
00190   char ids[]="QK";
00191   size_t slen=strlen(ids);
00192   total_bytes+=slen;
00193 #endif
00194 
00195   list2.check_buffer_size(total_bytes);
00196   void * tmpptr=list2.bptr;
00197 #if defined(SAFE_ALL)
00198   memcpy(list2,ids,slen);
00199 #endif
00200 
00201   fixed_smartlist2 & nlist2 = f1b2gradlist->nlist2;
00202   /*
00203   if (debugcounter++>569)
00204   {
00205     cout << debugcounter++ << endl;
00206     memcpy(list2,pz->get_u_bar(),2*nvar*sizeof(double));
00207     memcpy(list2,pz->get_u_bar(),nvar*sizeof(double));
00208     memcpy(list2,pz->get_u_dot_bar(),nvar*sizeof(double));
00209   }
00210   else
00211   */
00212   {
00213     const int sizeofdouble = sizeof(double);
00214     memcpy(list2,pz->get_u_bar(),nvar*sizeofdouble);
00215     memcpy(list2,pz->get_u_dot_bar(),nvar*sizeofdouble);
00216   }
00217   *nlist2.bptr=adptr_diff(list2.bptr,tmpptr);
00218   ++nlist2;
00219 
00220   // Do first reverse pass calculations
00221   for (size_t i=0;i<nvar;i++)
00222   {
00223     py->u_bar[i]+=xu*pz->u_bar[i];
00224   }
00225 
00226   for (size_t i=0;i<nvar;i++)
00227   {
00228     py->u_dot_bar[i]+=xu*pz->u_dot_bar[i];
00229   }
00230 
00231   // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00232   for (size_t i=0;i<nvar;i++)
00233   {
00234     pz->u_bar[i]=0;
00235   }
00236   for (size_t i=0;i<nvar;i++)
00237   {
00238     pz->u_dot_bar[i]=0;
00239   }
00240 }
00241 
00246 void read_pass2_2_prodc1(void)
00247 {
00248   //const int nlist_record_size=sizeof(int)+sizeof(char*);
00249   // We are going forward for bptr and backword for bptr2
00250   //
00251   // list 1
00252   //
00253   unsigned int nvar=df1b2variable::nvar;
00254   test_smartlist & list=f1b2gradlist->list;
00255 
00256   size_t total_bytes=2*sizeof(df1b2_header)
00257      +(nvar+2)*sizeof(double);
00258 // string identifier debug stuff
00259 #if defined(SAFE_ALL)
00260   char ids[]="DL";
00261   size_t slen=strlen(ids);
00262   total_bytes+=slen;
00263 #endif
00264 
00265   list.check_buffer_size(total_bytes);
00266 // end of string identifier debug stuff
00267 
00268   list.saveposition(); // save pointer to beginning of record;
00269   fixed_smartlist & nlist=f1b2gradlist->nlist;
00270    // nlist-=sizeof(int);
00271   // get record size
00272   int num_bytes=nlist.bptr->numbytes;
00273   //
00274   // list 2
00275   //
00276   test_smartlist & list2=f1b2gradlist->list2;
00277   fixed_smartlist2 & nlist2=f1b2gradlist->nlist2;
00278   // get record size
00279   int num_bytes2=*nlist2.bptr;
00280   --nlist2;
00281   // backup the size of the record
00282   list2-=num_bytes2;
00283   list2.saveposition(); // save pointer to beginning of record;
00284   // save the pointer to the beginning of the record
00285   // bptr and bptr2 now both point to the beginning of their records
00286 
00287   double yu;
00288   //df1b2_header x,z;
00289   //df1b2function2 * pf;
00290 
00291   // get info from tape1
00292   // get info from tape1
00293 #if defined(SAFE_ALL)
00294   checkidentiferstring("DL",list);
00295   checkidentiferstring("QK",list2);
00296 #endif
00297   double xu=*(double *) list.bptr;
00298   list.bptr+=sizeof(double);
00299   df1b2_header * py=(df1b2_header *) list.bptr;
00300   list.bptr+=sizeof(df1b2_header);
00301   df1b2_header * pz=(df1b2_header *) list.bptr;
00302   list.bptr+=sizeof(df1b2_header);
00303   //pf=*(df1b2function2 **) list.bptr;
00304   //list.bptr+=sizeof(char*);
00305   memcpy(&yu,list.bptr,sizeof(double));
00306   list.bptr+=sizeof(double);
00307 
00308   //double* ydot=(double*)list.bptr;
00309   list.restoreposition(num_bytes); // save pointer to beginning of record;
00310 
00311   //double* zbar=(double*)list2.bptr;
00312   //double* zdotbar=(double*)(list2.bptr+nvar*sizeof(double));
00313 
00314   list2.restoreposition(); // save pointer to beginning of record;
00315 
00316   //double * y_tilde=py->get_u_tilde();
00317   //double * y_dot_tilde=py->get_u_dot_tilde();
00318   double * y_bar_tilde=py->get_u_bar_tilde();
00319   double * y_dot_bar_tilde=py->get_u_dot_bar_tilde();
00320   double * z_bar_tilde=pz->get_u_bar_tilde();
00321   double * z_dot_bar_tilde=pz->get_u_dot_bar_tilde();
00322   // Do second "reverse-reverse" pass calculations
00323 
00324   for (size_t i=0;i<nvar;i++)
00325   {
00326     z_bar_tilde[i]=0;
00327     z_dot_bar_tilde[i]=0;
00328   }
00329 
00330   // start with y and add x
00331   for (size_t i=0;i<nvar;i++)
00332   {
00333     z_bar_tilde[i]+=xu*y_bar_tilde[i];
00334   }
00335 
00336   for (size_t i=0;i<nvar;i++)
00337   {
00338     z_dot_bar_tilde[i]+=xu*y_dot_bar_tilde[i];
00339   }
00340 }
00341 
00346 void read_pass2_3_prodc1(void)
00347 {
00348   // We are going backword for bptr and forward for bptr2
00349   // the current entry+2 in bptr is the size of the record i.e
00350   // points to the next record
00351   unsigned int nvar=df1b2variable::nvar;
00352   fixed_smartlist & nlist=f1b2gradlist->nlist;
00353   test_smartlist& list=f1b2gradlist->list;
00354    // nlist-=sizeof(int);
00355   // get record size
00356   int num_bytes=nlist.bptr->numbytes;
00357   // backup the size of the record
00358   list-=num_bytes;
00359   list.saveposition(); // save pointer to beginning of record;
00360   // save the pointer to the beginning of the record
00361   double yu;
00362   //df1b2_header z;
00363   //df1b2function2 * pf;
00364 
00365   // get info from tape1
00366   // get info from tape1
00367 #if defined(SAFE_ALL)
00368   checkidentiferstring("DL",list);
00369 #endif
00370   double xu=*(double*) list.bptr;
00371   list.bptr+=sizeof(double);
00372   df1b2_header * py=(df1b2_header *) list.bptr;
00373   list.bptr+=sizeof(df1b2_header);
00374   df1b2_header * pz=(df1b2_header *) list.bptr;
00375   list.bptr+=sizeof(df1b2_header);
00376   memcpy(&yu,list.bptr,sizeof(double));
00377   list.bptr+=sizeof(double);
00378 
00379   //double * ydot = (double*)list.bptr;
00380 
00381   list.restoreposition(); // save pointer to beginning of record;
00382 
00383   *(py->u_tilde)+=xu* *(pz->u_tilde);
00384   for (size_t i=0;i<nvar;i++)
00385   {
00386     py->u_dot_tilde[i]+=xu*pz->u_dot_tilde[i];
00387   }
00388   *(pz->u_tilde)=0;
00389   for (size_t i=0;i<nvar;i++)
00390   {
00391     pz->u_dot_tilde[i]=0;
00392   }
00393 }