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