ADMB Documentation  11.1.1890
 All Classes Files Functions Variables Typedefs Friends Defines
df1b2f21.cpp
Go to the documentation of this file.
00001 /*
00002  * $Id: df1b2f21.cpp 1110 2013-07-12 01:11:28Z johnoel $
00003  *
00004  * Author: David Fournier
00005  * Copyright (c) 2008-2012 Regents of the University of California
00006  */
00012 #include <df1b2fun.h>
00013 
00014 void ad_read_pass1_der_values(void);
00015 void read_pass1_1_dv(void);
00016 void read_pass1_2_dv(void);
00017 void read_pass1_3_dv(void);
00018 
00023 int df1b2_gradlist::write_pass1(const df1b2variable * _px,
00024   df1b2variable * pz, double df, double d2f, double d3f)
00025 {
00026   ADUNCONST(df1b2variable*,px)
00027   ncount++;
00028  /*
00029   if (28126<ncount && 28129> ncount)
00030     cout << "trap" << endl;
00031   */
00032 #if defined(CHECK_COUNT)
00033   if (ncount >= ncount_check)
00034     ncount_checker(ncount,ncount_check);
00035 #endif
00036   int nvar=df1b2variable::nvar;
00037 
00038   int total_bytes=sizeof(df1b2_header)+sizeof(df1b2_header)
00039    // +sizeof(char*)
00040     +4*sizeof(double)+nvar*sizeof(double);
00041 // string identifier debug stuff
00042 #if defined(SAFE_ALL)
00043   char ids[]="HU";
00044   int slen=strlen(ids);
00045   total_bytes+=slen;
00046 #endif
00047   list.check_buffer_size(total_bytes);
00048   void * tmpptr=list.bptr;
00049 #if defined(SAFE_ALL)
00050   memcpy(list,ids,slen);
00051 #endif
00052 // end of string identifier debug stuff
00053 
00054   memcpy(list,(df1b2_header*)(px),sizeof(df1b2_header));
00055   memcpy(list,(df1b2_header*)(pz),sizeof(df1b2_header));
00056 
00057   // replace this with the derivative values
00058   //memcpy(list,&pf,sizeof(char*));
00059 
00060   memcpy(list,&df,sizeof(double));
00061   memcpy(list,&d2f,sizeof(double));
00062   memcpy(list,&d3f,sizeof(double));
00063 
00064 
00065 
00066   memcpy(list,px->get_u(),sizeof(double));
00067   //list.bptr+=sizeof(double);
00068   memcpy(list,px->get_u_dot(),nvar*sizeof(double));
00069   //list.bptr+=nvar*sizeof(double);
00070   // ***** write  record size
00071   nlist.bptr->numbytes=adptr_diff(list.bptr,tmpptr);
00072   nlist.bptr->pf=(ADrfptr)(&ad_read_pass1_der_values);
00073   ++nlist;
00074   return 0;
00075 }
00076 
00081 void ad_read_pass1_der_values(void)
00082 {
00083   switch(df1b2variable::passnumber)
00084   {
00085   case 1:
00086     read_pass1_1_dv();
00087     break;
00088   case 2:
00089     read_pass1_2_dv();
00090     break;
00091   case 3:
00092     read_pass1_3_dv();
00093     break;
00094   default:
00095     cerr << "illegal value for df1b2variable::pass = "
00096          << df1b2variable::passnumber << endl;
00097     exit(1);
00098   }
00099 }
00100 
00105 void read_pass1_1_dv(void)
00106 {
00107   // We are going backword for bptr and forward for bptr2
00108   // the current entry+2 in bptr is the size of the record i.e
00109   // points to the next record
00110   int nvar=df1b2variable::nvar;
00111   fixed_smartlist & nlist=f1b2gradlist->nlist;
00112   test_smartlist& list=f1b2gradlist->list;
00113    // nlist-=sizeof(int);
00114   // get record size
00115   int num_bytes=nlist.bptr->numbytes;
00116   // backup the size of the record
00117   list-=num_bytes;
00118   //list.bptr-=num_bytes;
00119   list.saveposition(); // save pointer to beginning of record;
00120   // save the pointer to the beginning of the record
00121   double xu;
00122   double * xdot;
00123   //df1b2_header x,z;
00124   //df1b2function1 * pf;
00125 
00126   // get info from tape1
00127   // get info from tape1
00128 #if defined(SAFE_ARRAYS)
00129   checkidentiferstring("HU",list);
00130 #endif
00131   df1b2_header * px=(df1b2_header *) list.bptr;
00132   list.bptr+=sizeof(df1b2_header);
00133   df1b2_header * pz=(df1b2_header *) list.bptr;
00134   list.bptr+=sizeof(df1b2_header);
00135   double df=*(double*) list.bptr;
00136   list.bptr+=sizeof(double);
00137   double d2f=*(double*) list.bptr;
00138   list.bptr+=sizeof(double);
00139   //double d3f=*(double*) list.bptr;
00140   list.bptr+=sizeof(double);
00141   memcpy(&xu,list.bptr,sizeof(double));
00142   list.bptr+=sizeof(double);
00143   xdot=(double*)list.bptr;
00144   list.restoreposition(); // save pointer to beginning of record;
00145   int i;
00146 
00147   // Do first reverse paSS calculations
00148   // ****************************************************************
00149   // turn this off if no third derivatives are calculated
00150   // if (!no_third_derivatives)
00151   // {
00152   // save for second reverse pass
00153   // save identifier 1
00154      fixed_smartlist2& nlist2=f1b2gradlist->nlist2;
00155      test_smartlist& list2=f1b2gradlist->list2;
00156 
00157 
00158   int total_bytes=2*nvar*sizeof(double);
00159 // string identifier debug stuff
00160 #if defined(SAFE_ALL)
00161   char ids[]="JE";
00162   int slen=strlen(ids);
00163   total_bytes+=slen;
00164 #endif
00165 
00166   list2.check_buffer_size(total_bytes);
00167   void * tmpptr2=list2.bptr;
00168 
00169 #if defined(SAFE_ALL)
00170   memcpy(list2,ids,slen);
00171 #endif
00172 
00173    memcpy(list2,pz->get_u_bar(),nvar*sizeof(double));
00174    memcpy(list2,pz->get_u_dot_bar(),nvar*sizeof(double));
00175    *nlist2.bptr=adptr_diff(list2.bptr,tmpptr2);
00176    ++nlist2;
00177   // }
00178   //
00179   // ****************************************************************
00180 #if defined(PRINT_DERS)
00181  print_derivatives(pz,"z");
00182  print_derivatives(px,"x");
00183 #endif
00184 
00185   //double df=(pf->df)(xu);
00186   //double d2f=(pf->d2f)(xu);
00187   //double d3f=(pf->d3f)(xu);
00188 
00189   for (i=0;i<nvar;i++)
00190   {
00191     //px->u_bar[i]+=(pf->df)(xu)* pz->u_bar[i];
00192     px->u_bar[i]+=df * pz->u_bar[i];
00193   }
00194   for (i=0;i<nvar;i++)
00195   {
00196     //px->u_bar[i]+=(pf->d2f)(xu)*xdot[i]*pz->u_dot_bar[i];
00197     px->u_bar[i]+=d2f*xdot[i]*pz->u_dot_bar[i];
00198   }
00199   for (i=0;i<nvar;i++)
00200   {
00201     //px->u_dot_bar[i]+=(pf->df)(xu)*pz->u_dot_bar[i];
00202     px->u_dot_bar[i]+=df*pz->u_dot_bar[i];
00203   }
00204 
00205   // !!!!!!!!!!!!!!!!!!!!!!
00206   for (i=0;i<nvar;i++)
00207   {
00208     pz->u_bar[i]=0;
00209   }
00210   for (i=0;i<nvar;i++)
00211   {
00212     pz->u_dot_bar[i]=0;
00213   }
00214 
00215 #if defined(PRINT_DERS)
00216  print_derivatives(pz,"z");
00217  print_derivatives(px,"x");
00218 #endif
00219 }
00220 
00225 void read_pass1_2_dv(void)
00226 {
00227   //const int nlist_record_size=sizeof(int)+sizeof(char*);
00228   // We are going forward for bptr and backword for bptr2
00229   //
00230   // list 1
00231   //
00232   int nvar=df1b2variable::nvar;
00233   test_smartlist & list=f1b2gradlist->list;
00234 
00235   int total_bytes=sizeof(df1b2_header)+sizeof(df1b2_header)+sizeof(char*)
00236     +sizeof(double)+nvar*sizeof(double);
00237 // string identifier debug stuff
00238 #if defined(SAFE_ALL)
00239   char ids[]="HU";
00240   int slen=strlen(ids);
00241   total_bytes+=slen;
00242 #endif
00243   list.check_buffer_size(total_bytes);
00244 // end of string identifier debug stuff
00245 
00246   list.saveposition(); // save pointer to beginning of record;
00247   fixed_smartlist & nlist=f1b2gradlist->nlist;
00248    // nlist-=sizeof(int);
00249   // get record size
00250   int num_bytes=nlist.bptr->numbytes;
00251     // nlist+=nlist_record_size;
00252   //
00253   // list 2
00254   //
00255   test_smartlist & list2=f1b2gradlist->list2;
00256   fixed_smartlist2 & nlist2=f1b2gradlist->nlist2;
00257   // get record size
00258   int num_bytes2=*nlist2.bptr;
00259   --nlist2;
00260   // backup the size of the record
00261   list2-=num_bytes2;
00262   list2.saveposition();
00263   // save the pointer to the beginning of the record
00264   // bptr and bptr2 now both point to the beginning of their records
00265 
00266   double xu;
00267   double * xdot;
00268   //df1b2_header x,z;
00269   //df1b2function1 * pf;
00270 
00271   // get info from tape1
00272   // get info from tape1
00273 #if defined(SAFE_ARRAYS)
00274   checkidentiferstring("HU",list);
00275   checkidentiferstring("JE",list2);
00276 #endif
00277   df1b2_header * px=(df1b2_header *) list.bptr;
00278   list.bptr+=sizeof(df1b2_header);
00279   df1b2_header * pz=(df1b2_header *) list.bptr;
00280   list.bptr+=sizeof(df1b2_header);
00281 
00282   //pf=*(df1b2function1 **) list.bptr;
00283 
00284   double df=*(double*) list.bptr;
00285   list.bptr+=sizeof(double);
00286   double d2f=*(double*) list.bptr;
00287   list.bptr+=sizeof(double);
00288   double d3f=*(double*) list.bptr;
00289   list.bptr+=sizeof(double);
00290   //list.bptr+=sizeof(char*);
00291   memcpy(&xu,list.bptr,sizeof(double));
00292   list.bptr+=sizeof(double);
00293   xdot=(double*)list.bptr;
00294   list.restoreposition(num_bytes); // save pointer to beginning of record;
00295 
00296   double * zbar;
00297   double * zdotbar;
00298 
00299   zbar=(double*)list2.bptr;
00300   zdotbar=(double*)(list2.bptr+nvar*sizeof(double));
00301 
00302   double * x_tilde=px->get_u_tilde();
00303   double * x_dot_tilde=px->get_u_dot_tilde();
00304   double * x_bar_tilde=px->get_u_bar_tilde();
00305   double * x_dot_bar_tilde=px->get_u_dot_bar_tilde();
00306   double * z_bar_tilde=pz->get_u_bar_tilde();
00307   double * z_dot_bar_tilde=pz->get_u_dot_bar_tilde();
00308 #if defined(PRINT_DERS)
00309  print_derivatives(pz,"z");
00310  print_derivatives(px,"x");
00311 #endif
00312   // Do second "reverse-reverse" pass calculations
00313   int i;
00314 
00315   // !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00316   for (i=0;i<nvar;i++)
00317   {
00318     z_bar_tilde[i]=0;
00319     z_dot_bar_tilde[i]=0;
00320   }
00321 
00322   //double df=(pf->df)(xu);
00323   //double d2f=(pf->d2f)(xu);
00324   //double d3f=(pf->d3f)(xu);
00325   for (i=0;i<nvar;i++)
00326   {
00327     //*x_tilde+=(pf->d2f)(xu)*zbar[i]*x_bar_tilde[i];
00328     *x_tilde+=d2f*zbar[i]*x_bar_tilde[i];
00329     //z_bar_tilde[i]+=(pf->df)(xu)*x_bar_tilde[i];
00330     z_bar_tilde[i]+=df*x_bar_tilde[i];
00331   }
00332 
00333   for (i=0;i<nvar;i++)
00334   {
00335     //*x_tilde+=(pf->d2f)(xu)*zdotbar[i]*x_dot_bar_tilde[i];
00336     *x_tilde+=d2f*zdotbar[i]*x_dot_bar_tilde[i];
00337     //z_dot_bar_tilde[i]+=(pf->df)(xu)*x_dot_bar_tilde[i];
00338     z_dot_bar_tilde[i]+=df*x_dot_bar_tilde[i];
00339   }
00340 
00341   for (i=0;i<nvar;i++)
00342   {
00343     //x_dot_tilde[i]+=(pf->d2f)(xu)*zdotbar[i]*x_bar_tilde[i];
00344     //z_dot_bar_tilde[i]+=(pf->d2f)(xu)*xdot[i]*x_bar_tilde[i];
00345     //*x_tilde+=(pf->d3f)(xu)*xdot[i]*zdotbar[i]*x_bar_tilde[i];
00346     x_dot_tilde[i]+=d2f*zdotbar[i]*x_bar_tilde[i];
00347     z_dot_bar_tilde[i]+=d2f*xdot[i]*x_bar_tilde[i];
00348     *x_tilde+=d3f*xdot[i]*zdotbar[i]*x_bar_tilde[i];
00349   }
00350   list2.restoreposition();
00351 #if defined(PRINT_DERS)
00352  print_derivatives(pz,"z");
00353  print_derivatives(px,"x");
00354 #endif
00355 }
00356 
00361 void read_pass1_3_dv(void)
00362 {
00363   // We are going backword for bptr and forward for bptr2
00364   // the current entry+2 in bptr is the size of the record i.e
00365   // points to the next record
00366   int nvar=df1b2variable::nvar;
00367   fixed_smartlist & nlist=f1b2gradlist->nlist;
00368   test_smartlist& list=f1b2gradlist->list;
00369    // nlist-=sizeof(int);
00370   // get record size
00371   int num_bytes=nlist.bptr->numbytes;
00372   // backup the size of the record
00373   list-=num_bytes;
00374   list.saveposition(); // save pointer to beginning of record;
00375   // save the pointer to the beginning of the record
00376   double xu;
00377   double * xdot;
00378   //df1b2_header x,z;
00379   //df1b2function1 * pf;
00380 
00381   // get info from tape1
00382 #if defined(SAFE_ARRAYS)
00383   checkidentiferstring("HU",list);
00384 #endif
00385   df1b2_header * px=(df1b2_header *) list.bptr;
00386   list.bptr+=sizeof(df1b2_header);
00387   df1b2_header * pz=(df1b2_header *) list.bptr;
00388   list.bptr+=sizeof(df1b2_header);
00389   //pf=*(df1b2function1 **) list.bptr;
00390   //list.bptr+=sizeof(char*);
00391   double df=*(double*) list.bptr;
00392   list.bptr+=sizeof(double);
00393   double d2f=*(double*) list.bptr;
00394   list.bptr+=sizeof(double);
00395   //double d3f=*(double*) list.bptr;
00396   list.bptr+=sizeof(double);
00397   memcpy(&xu,list.bptr,sizeof(double));
00398   list.bptr+=sizeof(double);
00399   xdot=(double*)list.bptr;
00400   list.restoreposition(); // save pointer to beginning of record;
00401   int i;
00402 
00403 #if defined(PRINT_DERS)
00404  print_derivatives(pz,"z");
00405  print_derivatives(px,"x");
00406 #endif
00407 
00408   //double df=(pf->df)(xu);
00409   //double d2f=(pf->d2f)(xu);
00410   //*(px->u_tilde)+=(pf->df)(xu)* *(pz->u_tilde);
00411   *(px->u_tilde)+=df * *(pz->u_tilde);
00412   for (i=0;i<nvar;i++)
00413   {
00414     //*(px->u_tilde)+=(pf->d2f)(xu)*xdot[i]*pz->u_dot_tilde[i];
00415     *(px->u_tilde)+=d2f*xdot[i]*pz->u_dot_tilde[i];
00416   }
00417   for (i=0;i<nvar;i++)
00418   {
00419     //px->u_dot_tilde[i]+=(pf->df)(xu)*pz->u_dot_tilde[i];
00420     px->u_dot_tilde[i]+=df*pz->u_dot_tilde[i];
00421   }
00422   *(pz->u_tilde)=0;
00423   for (i=0;i<nvar;i++)
00424   {
00425     pz->u_dot_tilde[i]=0;
00426   }
00427 #if defined(PRINT_DERS)
00428  print_derivatives(pz,"z");
00429  print_derivatives(px,"x");
00430 #endif
00431 }