ADMB Documentation  11.6rc.3396
 All Classes Files Functions Variables Typedefs Friends Defines
df33fun1.cpp
Go to the documentation of this file.
00001 /*
00002  * $Id$
00003  *
00004  * Author: David Fournier
00005  * Copyright (c) 2008-2012 Regents of the University of California
00006  */
00011 #include <df1b2fun.h>
00012 #include "df33fun.h"
00013 
00014 #ifndef OPT_LIB
00015   #include <cassert>
00016   #include <climits>
00017 #endif
00018 
00019 void ad_read_pass2_dvdvdv(void);
00020 
00025  int df1b2_gradlist::write_pass1(const df1b2variable * _px,
00026    const df1b2variable * _py,const df1b2variable * pw,
00027    const df1b2variable * pz,
00028    double df_x, double df_y, double df_z,
00029    double df_xx, double df_xy,double df_xz,double df_yy,
00030    double df_yz,double df_zz,
00031    double df_xxx,
00032    double df_xxy,
00033    double df_xxz,
00034    double df_xyy,
00035    double df_xyz,
00036    double df_xzz,
00037    double df_yyy,
00038    double df_yyz,
00039    double df_yzz,
00040    double df_zzz)
00041  {
00042    ADUNCONST(df1b2variable*,px)
00043    ADUNCONST(df1b2variable*,py)
00044    ncount++;
00045 #if defined(CHECK_COUNT)
00046   if (ncount >= ncount_check)
00047     ncount_checker(ncount,ncount_check);
00048 #endif
00049   size_t nvar=df1b2variable::nvar;
00050 
00051    size_t total_bytes=4*sizeof(df1b2_header)+sizeof(char*)
00052      +(3*nvar+22)*sizeof(double);
00053 
00054 // string identifier debug stuff
00055 #if defined(SAFE_ALL)
00056   char ids[]="U8";
00057   int slen=strlen(ids);
00058   total_bytes+=slen;
00059 #endif
00060 
00061   list.check_buffer_size(total_bytes);
00062 
00063   void * tmpptr=list.bptr;
00064 #if defined(SAFE_ALL)
00065   memcpy(list,ids,slen);
00066 #endif
00067 // end of string identifier debug stuff
00068 
00069    memcpy(list,(df1b2_header*)(px),sizeof(df1b2_header));
00070    memcpy(list,(df1b2_header*)(py),sizeof(df1b2_header));
00071    memcpy(list,(df1b2_header*)(pw),sizeof(df1b2_header));
00072    memcpy(list,(df1b2_header*)(pz),sizeof(df1b2_header));
00073    //memcpy(list,&pf,sizeof(char *));
00074    //*(char**)(list.bptr)=(char*)pf;
00075 
00076   /*
00077    memcpy(list,&df_x,sizeof(double));
00078    memcpy(list,&df_y,sizeof(double));
00079    memcpy(list,&df_xx,sizeof(double));
00080    memcpy(list,&df_xy,sizeof(double));
00081    memcpy(list,&df_yy,sizeof(double));
00082    memcpy(list,&df_xxx,sizeof(double));
00083    memcpy(list,&df_xxy,sizeof(double));
00084    memcpy(list,&df_xyy,sizeof(double));
00085    memcpy(list,&df_yyy,sizeof(double));
00086   */
00087    size_t sizeofdouble = sizeof(double);
00088    memcpy(list,&df_x,sizeofdouble);
00089    memcpy(list,&df_y,sizeofdouble);
00090    memcpy(list,&df_z,sizeofdouble);
00091    memcpy(list,&df_xx,sizeofdouble);
00092    memcpy(list,&df_xy,sizeofdouble);
00093    memcpy(list,&df_xz,sizeofdouble);
00094    memcpy(list,&df_yy,sizeofdouble);
00095    memcpy(list,&df_yz,sizeofdouble);
00096    memcpy(list,&df_zz,sizeofdouble);
00097    memcpy(list,&df_xxx,sizeofdouble);
00098    memcpy(list,&df_xxy,sizeofdouble);
00099    memcpy(list,&df_xxz,sizeofdouble);
00100    memcpy(list,&df_xyy,sizeofdouble);
00101    memcpy(list,&df_xyz,sizeofdouble);
00102    memcpy(list,&df_xzz,sizeofdouble);
00103    memcpy(list,&df_yyy,sizeofdouble);
00104    memcpy(list,&df_yyz,sizeofdouble);
00105    memcpy(list,&df_yzz,sizeofdouble);
00106    memcpy(list,&df_zzz,sizeofdouble);
00107 
00108    memcpy(list,px->get_u(),sizeofdouble);
00109    memcpy(list,py->get_u(),sizeofdouble);
00110    memcpy(list,pw->get_u(),sizeofdouble);
00111    memcpy(list,px->get_u_dot(),nvar*sizeofdouble);
00112    memcpy(list,py->get_u_dot(),nvar*sizeofdouble);
00113    memcpy(list,pw->get_u_dot(),nvar*sizeofdouble);
00114    // ***** write  record size
00115    nlist.bptr->numbytes=adptr_diff(list.bptr,tmpptr);
00116    nlist.bptr->pf=(ADrfptr)(&ad_read_pass2_dvdvdv);
00117    ++nlist;
00118    return 0;
00119  }
00120 
00121 
00122 void read_pass2_1_dvdvdv(void);
00123 void read_pass2_2_dvdvdv(void);
00124 void read_pass2_3_dvdvdv(void);
00125 
00130 void ad_read_pass2_dvdvdv(void)
00131 {
00132   switch(df1b2variable::passnumber)
00133   {
00134   case 1:
00135     read_pass2_1_dvdvdv();
00136     break;
00137   case 2:
00138     read_pass2_2_dvdvdv();
00139     break;
00140   case 3:
00141     read_pass2_3_dvdvdv();
00142     break;
00143   default:
00144     cerr << "illegal value for df1b2variable::pass = "
00145          << df1b2variable::passnumber << endl;
00146     exit(1);
00147   }
00148 }
00149 
00154 void read_pass2_1_dvdvdv(void)
00155 {
00156   // We are going backword for bptr and nbptr
00157   // and  forward for bptr2 and nbptr2
00158   // the current entry+2 in bptr is the size of the record i.e
00159   // points to the next record
00160   //char * bptr=f1b2gradlist->bptr;
00161   //char * bptr2=f1b2gradlist2->bptr;
00162   size_t nvar = df1b2variable::nvar;
00163   test_smartlist& list=f1b2gradlist->list;
00164   //f1b2gradlist->nlist-=sizeof(int);
00165   int num_bytes=f1b2gradlist->nlist.bptr->numbytes;
00166   list-=num_bytes;
00167   list.saveposition(); // save pointer to beginning of record;
00168   double xu,yu,wu;
00169   //ad_dstar xdot,ydot;
00170   //df1b2function2 * pf;
00171 
00172   // get info from tape1
00173 #if defined(SAFE_ARRAYS)
00174   checkidentiferstring("U8",f1b2gradlist->list);
00175 #endif
00176   char * bptr=f1b2gradlist->list.bptr;
00177   df1b2_header * px=(df1b2_header *) bptr;
00178   bptr+=sizeof(df1b2_header);
00179   df1b2_header * py=(df1b2_header *) bptr;
00180   bptr+=sizeof(df1b2_header);
00181   df1b2_header * pw=(df1b2_header *) bptr;
00182   bptr+=sizeof(df1b2_header);
00183   df1b2_header * pz=(df1b2_header *) bptr;
00184   bptr+=sizeof(df1b2_header);
00185   //pf=*(df1b2function2 **) bptr;
00186   //bptr+=sizeof(char*);
00187 
00188   double df1=*(double*) bptr;
00189   bptr+=sizeof(double);
00190 
00191   double df2=*(double*) bptr;
00192   bptr+=sizeof(double);
00193 
00194   double df3=*(double*) bptr;
00195   bptr+=sizeof(double);
00196 
00197   double d2f11=*(double*) bptr;
00198   bptr+=sizeof(double);
00199 
00200   double d2f12=*(double*) bptr;
00201   bptr+=sizeof(double);
00202 
00203   double d2f13=*(double*) bptr;
00204   bptr+=sizeof(double);
00205 
00206   double d2f22=*(double*) bptr;
00207   bptr+=sizeof(double);
00208 
00209   double d2f23=*(double*) bptr;
00210   bptr+=sizeof(double);
00211 
00212   double d2f33=*(double*) bptr;
00213   bptr+=sizeof(double);
00214 
00215 #if defined(PRINT_DERS)
00216   double d3f111=*(double*) bptr;
00217 #endif
00218   bptr+=sizeof(double);
00219 
00220 #if defined(PRINT_DERS)
00221   double d3f112=*(double*) bptr;
00222 #endif
00223   bptr+=sizeof(double);
00224 
00225 #if defined(PRINT_DERS)
00226   double d3f113=*(double*) bptr;
00227 #endif
00228   bptr+=sizeof(double);
00229 
00230 #if defined(PRINT_DERS)
00231   double d3f122=*(double*) bptr;
00232 #endif
00233   bptr+=sizeof(double);
00234 
00235 #if defined(PRINT_DERS)
00236   double d3f123=*(double*) bptr;
00237 #endif
00238   bptr+=sizeof(double);
00239 
00240 #if defined(PRINT_DERS)
00241   double d3f133=*(double*) bptr;
00242 #endif
00243   bptr+=sizeof(double);
00244 
00245 #if defined(PRINT_DERS)
00246   double d3f222=*(double*) bptr;
00247 #endif
00248   bptr+=sizeof(double);
00249 
00250 #if defined(PRINT_DERS)
00251   double d3f223=*(double*) bptr;
00252 #endif
00253   bptr+=sizeof(double);
00254 
00255 #if defined(PRINT_DERS)
00256   double d3f233=*(double*) bptr;
00257 #endif
00258   bptr+=sizeof(double);
00259 
00260 #if defined(PRINT_DERS)
00261   double d3f333=*(double*) bptr;
00262 #endif
00263   bptr+=sizeof(double);
00264 
00265   memcpy(&xu,bptr,sizeof(double));
00266   bptr+=sizeof(double);
00267   memcpy(&yu,bptr,sizeof(double));
00268   bptr+=sizeof(double);
00269   memcpy(&wu,bptr,sizeof(double));
00270   bptr+=sizeof(double);
00271   double * xdot=(double*)bptr;
00272   bptr+=nvar*sizeof(double);
00273   double * ydot=(double*)bptr;
00274   bptr+=nvar*sizeof(double);
00275   double * wdot=(double*)bptr;
00276 
00277   list.restoreposition(); // save pointer to beginning of record;
00278 
00279   // ****************************************************************
00280   // turn this off if no third derivatives are calculated
00281   // if (!no_third_derivatives)
00282   // {
00283   // save for second reverse pass
00284   // save identifier 1
00285      test_smartlist & list2 = f1b2gradlist->list2;
00286 
00287 
00288    size_t total_bytes=2*nvar*sizeof(double);
00289 // string identifier debug stuff
00290 #if defined(SAFE_ALL)
00291   char ids[]="FX";
00292   int slen=strlen(ids);
00293   total_bytes+=slen;
00294 #endif
00295 
00296   list2.check_buffer_size(total_bytes);
00297   void * tmpptr=list2.bptr;
00298 #if defined(SAFE_ALL)
00299   memcpy(list2,ids,slen);
00300 #endif
00301 
00302      fixed_smartlist2 & nlist2 = f1b2gradlist->nlist2;
00303      size_t sizeofdouble = sizeof(double);
00304      memcpy(list2,pz->get_u_bar(),nvar*sizeofdouble);
00305      memcpy(list2,pz->get_u_dot_bar(),nvar*sizeofdouble);
00306      *nlist2.bptr=adptr_diff(list2.bptr,tmpptr);
00307      ++nlist2;
00308   // }
00309   //
00310   // ****************************************************************
00311 #if defined(PRINT_DERS)
00312  print_derivatives(funname,(f),(df1),
00313   (df2),(d2f11),(d2f12),(d2f22),
00314   (d3f111),(d3f112),(d3f122),
00315   (d3f222),1);
00316  print_derivatives(pz,"z");
00317  print_derivatives(px,"x");
00318  print_derivatives(py,"y");
00319 #endif
00320 #if defined(__DERCHECK__)
00321   if (derchecker)
00322   if (derchecker->node_number)
00323   {
00324     if (derchecker->counter == derchecker->node_number)
00325     {
00326       switch (derchecker->pass_number) // increment the variable of interest
00327       {
00328       case 2:
00329         switch(derchecker->vartype)
00330         {
00331         case 1:
00332           if (!derchecker->dotflag)
00333             px->u_bar[derchecker->index-1]+=derchecker->delta;
00334           else
00335             px->u_dot_bar[derchecker->index-1]+=derchecker->delta;
00336           break;
00337         case 2:
00338           if (!derchecker->dotflag)
00339             py->u_bar[derchecker->index-1]+=derchecker->delta;
00340           else
00341             py->u_dot_bar[derchecker->index-1]+=derchecker->delta;
00342           break;
00343         case 3:
00344           if (!derchecker->dotflag)
00345             pz->u_bar[derchecker->index-1]+=derchecker->delta;
00346           else
00347             pz->u_dot_bar[derchecker->index-1]+=derchecker->delta;
00348           break;
00349         default:
00350           cerr << "Invalid index value for dercheck_index was "
00351                << derchecker->index << endl;
00352           break;
00353         }
00354         break;
00355       case 3:
00356         switch(derchecker->vartype)
00357         {
00358         case 1:
00359           if (!derchecker->dotflag)
00360             px->u_bar[derchecker->index-1]-=derchecker->delta;
00361           else
00362             px->u_dot_bar[derchecker->index-1]-=derchecker->delta;
00363           break;
00364         case 2:
00365           if (!derchecker->dotflag)
00366             py->u_bar[derchecker->index-1]-=derchecker->delta;
00367           else
00368             py->u_dot_bar[derchecker->index-1]-=derchecker->delta;
00369           break;
00370         case 3:
00371           if (!derchecker->dotflag)
00372             pz->u_bar[derchecker->index-1]-=derchecker->delta;
00373           else
00374             pz->u_dot_bar[derchecker->index-1]-=derchecker->delta;
00375           break;
00376         default:
00377           cerr << "Invalid index value for dercheck_index was "
00378                << derchecker->index << endl;
00379           break;
00380         }
00381         break;
00382       }
00383     }
00384   }
00385 #endif
00386 
00387   // Do first reverse pass calculations
00388   for (size_t i=0;i<nvar;i++)
00389   {
00390     px->u_bar[i]+=(df1)*pz->u_bar[i];
00391   }
00392   for (size_t i=0;i<nvar;i++)
00393   {
00394     py->u_bar[i]+=(df2)*pz->u_bar[i];
00395   }
00396   for (size_t i=0;i<nvar;i++)
00397   {
00398     pw->u_bar[i]+=(df3)*pz->u_bar[i];
00399   }
00400   for (size_t i=0;i<nvar;i++)
00401   {
00402     px->u_bar[i]+=(d2f11)*xdot[i]*pz->u_dot_bar[i];
00403     px->u_bar[i]+=(d2f12)*ydot[i]*pz->u_dot_bar[i];
00404     px->u_bar[i]+=(d2f13)*wdot[i]*pz->u_dot_bar[i];
00405   }
00406   for (size_t i=0;i<nvar;i++)
00407   {
00408     py->u_bar[i]+=(d2f12)*xdot[i]*pz->u_dot_bar[i];
00409     py->u_bar[i]+=(d2f22)*ydot[i]*pz->u_dot_bar[i];
00410     //py->u_bar[i]+=(d2f13)*wdot[i]*pz->u_dot_bar[i];
00411     py->u_bar[i]+=(d2f23)*wdot[i]*pz->u_dot_bar[i];
00412   }
00413   for (size_t i=0;i<nvar;i++)
00414   {
00415     pw->u_bar[i]+=(d2f13)*xdot[i]*pz->u_dot_bar[i];
00416     pw->u_bar[i]+=(d2f23)*ydot[i]*pz->u_dot_bar[i];
00417     pw->u_bar[i]+=(d2f33)*wdot[i]*pz->u_dot_bar[i];
00418   }
00419   for (size_t i=0;i<nvar;i++)
00420   {
00421     px->u_dot_bar[i]+=(df1)*pz->u_dot_bar[i];
00422   }
00423   for (size_t i=0;i<nvar;i++)
00424   {
00425     py->u_dot_bar[i]+=(df2)*pz->u_dot_bar[i];
00426   }
00427   for (size_t i=0;i<nvar;i++)
00428   {
00429     pw->u_dot_bar[i]+=(df3)*pz->u_dot_bar[i];
00430   }
00431 
00432   for (size_t i=0;i<nvar;i++)
00433   {
00434     pz->u_bar[i]=0;
00435   }
00436   for (size_t i=0;i<nvar;i++)
00437   {
00438     pz->u_dot_bar[i]=0;
00439   }
00440 }
00441 
00446 void read_pass2_2_dvdvdv(void)
00447 {
00448   //const int nlist_record_size=sizeof(int)+sizeof(char*);
00449   // We are going forward for bptr and backword for bptr2
00450   //
00451   // list 1
00452   //
00453   unsigned int nvar=df1b2variable::nvar;
00454   test_smartlist & list=f1b2gradlist->list;
00455 
00456   size_t total_bytes=3*sizeof(df1b2_header)+sizeof(char*)
00457     +2*(nvar+1)*sizeof(double);
00458 // string identifier debug stuff
00459 #if defined(SAFE_ALL)
00460   char ids[]="U8";
00461   int slen=strlen(ids);
00462   total_bytes+=slen;
00463 #endif
00464 
00465   list.check_buffer_size(total_bytes);
00466 // end of string identifier debug stuff
00467 
00468   list.saveposition(); // save pointer to beginning of record;
00469   fixed_smartlist & nlist=f1b2gradlist->nlist;
00470    // nlist-=sizeof(int);
00471   // get record size
00472   int num_bytes=nlist.bptr->numbytes;
00473     // nlist+=nlist_record_size;
00474   //
00475   // list 2
00476   //
00477   test_smartlist & list2=f1b2gradlist->list2;
00478   fixed_smartlist2 & nlist2=f1b2gradlist->nlist2;
00479   // get record size
00480   int num_bytes2=*nlist2.bptr;
00481   --nlist2;
00482   // backup the size of the record
00483   list2-=num_bytes2;
00484   list2.saveposition(); // save pointer to beginning of record;
00485   // save the pointer to the beginning of the record
00486   // bptr and bptr2 now both point to the beginning of their records
00487 
00488   double xu,yu,wu;
00489   //df1b2_header x,z;
00490   //df1b2function2 * pf;
00491 
00492   // get info from tape1
00493   // get info from tape1
00494 #if defined(SAFE_ARRAYS)
00495   checkidentiferstring("U8",list);
00496   checkidentiferstring("FX",list2);
00497 #endif
00498   /*
00499   df1b2_header * px=(df1b2_header *) list.bptr;
00500   list.bptr+=sizeof(df1b2_header);
00501   df1b2_header * py=(df1b2_header *) list.bptr;
00502   list.bptr+=sizeof(df1b2_header);
00503   df1b2_header * pz=(df1b2_header *) list.bptr;
00504   list.bptr+=sizeof(df1b2_header);
00505   pf=*(df1b2function2 **) list.bptr;
00506   list.bptr+=sizeof(char*);
00507   memcpy(&xu,list.bptr,sizeof(double));
00508   list.bptr+=sizeof(double);
00509   memcpy(&yu,list.bptr,sizeof(double));
00510   list.bptr+=sizeof(double);
00511   xdot=(double*)list.bptr;
00512   list.bptr+=nvar*sizeof(double);
00513   ydot=(double*)list.bptr;
00514   */
00515   //char * bptr=f1b2gradlist->list.bptr;
00516   df1b2_header * px=(df1b2_header *) list.bptr;
00517   list.bptr+=sizeof(df1b2_header);
00518   df1b2_header * py=(df1b2_header *) list.bptr;
00519   list.bptr+=sizeof(df1b2_header);
00520   df1b2_header * pw=(df1b2_header *) list.bptr;
00521   list.bptr+=sizeof(df1b2_header);
00522   df1b2_header * pz=(df1b2_header *) list.bptr;
00523   list.bptr+=sizeof(df1b2_header);
00524 
00525   double df1=*(double*) list.bptr;
00526   list.bptr+=sizeof(double);
00527 
00528   double df2=*(double*) list.bptr;
00529   list.bptr+=sizeof(double);
00530 
00531   double df3=*(double*) list.bptr;
00532   list.bptr+=sizeof(double);
00533 
00534   double d2f11=*(double*) list.bptr;
00535   list.bptr+=sizeof(double);
00536 
00537   double d2f12=*(double*) list.bptr;
00538   list.bptr+=sizeof(double);
00539 
00540   double d2f13=*(double*) list.bptr;
00541   list.bptr+=sizeof(double);
00542 
00543   double d2f22=*(double*) list.bptr;
00544   list.bptr+=sizeof(double);
00545 
00546   double d2f23=*(double*) list.bptr;
00547   list.bptr+=sizeof(double);
00548 
00549   double d2f33=*(double*) list.bptr;
00550   list.bptr+=sizeof(double);
00551 
00552   double d3f111=*(double*) list.bptr;
00553   list.bptr+=sizeof(double);
00554 
00555   double d3f112=*(double*) list.bptr;
00556   list.bptr+=sizeof(double);
00557 
00558   double d3f113=*(double*) list.bptr;
00559   list.bptr+=sizeof(double);
00560 
00561   double d3f122=*(double*) list.bptr;
00562   list.bptr+=sizeof(double);
00563 
00564   double d3f123=*(double*) list.bptr;
00565   list.bptr+=sizeof(double);
00566 
00567   double d3f133=*(double*) list.bptr;
00568   list.bptr+=sizeof(double);
00569 
00570   double d3f222=*(double*) list.bptr;
00571   list.bptr+=sizeof(double);
00572 
00573   double d3f223=*(double*) list.bptr;
00574   list.bptr+=sizeof(double);
00575 
00576   double d3f233=*(double*) list.bptr;
00577   list.bptr+=sizeof(double);
00578 
00579   double d3f333=*(double*) list.bptr;
00580   list.bptr+=sizeof(double);
00581 
00582   memcpy(&xu,list.bptr,sizeof(double));
00583   list.bptr+=sizeof(double);
00584   memcpy(&yu,list.bptr,sizeof(double));
00585   list.bptr+=sizeof(double);
00586   memcpy(&wu,list.bptr,sizeof(double));
00587   list.bptr+=sizeof(double);
00588   double * xdot=(double*)list.bptr;
00589   list.bptr+=nvar*sizeof(double);
00590   double * ydot=(double*)list.bptr;
00591   list.bptr+=nvar*sizeof(double);
00592   double * wdot=(double*)list.bptr;
00593 
00594   list.restoreposition(num_bytes); // save pointer to beginning of record;
00595 
00596   double * zbar;
00597   double * zdotbar;
00598 
00599 
00600   zbar=(double*)list2.bptr;
00601   zdotbar=(double*)(list2.bptr+nvar*sizeof(double));
00602   list2.restoreposition(); // save pointer to beginning of record;
00603 
00604   double * x_tilde=px->get_u_tilde();
00605   double * x_dot_tilde=px->get_u_dot_tilde();
00606   double * x_bar_tilde=px->get_u_bar_tilde();
00607   double * x_dot_bar_tilde=px->get_u_dot_bar_tilde();
00608   double * y_tilde=py->get_u_tilde();
00609   double * y_dot_tilde=py->get_u_dot_tilde();
00610   double * y_bar_tilde=py->get_u_bar_tilde();
00611   double * y_dot_bar_tilde=py->get_u_dot_bar_tilde();
00612   double * w_tilde=pw->get_u_tilde();
00613   double * w_dot_tilde=pw->get_u_dot_tilde();
00614   double * w_bar_tilde=pw->get_u_bar_tilde();
00615   double * w_dot_bar_tilde=pw->get_u_dot_bar_tilde();
00616   double * z_bar_tilde=pz->get_u_bar_tilde();
00617   double * z_dot_bar_tilde=pz->get_u_dot_bar_tilde();
00618   // Do second "reverse-reverse" pass calculations
00619 #if defined(PRINT_DERS)
00620  print_derivatives(funname,(f),(df1),
00621   (df2),(d2f11),(d2f12),(d2f22),
00622   (d3f111),(d3f112),(d3f122),
00623   (d3f222),1);
00624  print_derivatives(pz,"z");
00625  print_derivatives(px,"x");
00626  print_derivatives(pw,"x");
00627  print_derivatives(py,"y");
00628 #endif
00629 
00630   for (size_t i=0;i<nvar;i++)
00631   {
00632     z_bar_tilde[i]=0;
00633     z_dot_bar_tilde[i]=0;
00634   }
00635 
00636   for (size_t i=0;i<nvar;i++)
00637   {
00638     *x_tilde+=(d2f11)*zbar[i]*x_bar_tilde[i];
00639     *y_tilde+=(d2f12)*zbar[i]*x_bar_tilde[i];
00640     *w_tilde+=(d2f13)*zbar[i]*x_bar_tilde[i];
00641 
00642     *x_tilde+=(d2f12)*zbar[i]*y_bar_tilde[i];
00643     *y_tilde+=(d2f22)*zbar[i]*y_bar_tilde[i];
00644     *w_tilde+=(d2f23)*zbar[i]*y_bar_tilde[i];
00645 
00646     *x_tilde+=(d2f13)*zbar[i]*w_bar_tilde[i];
00647     *y_tilde+=(d2f23)*zbar[i]*w_bar_tilde[i];
00648     *w_tilde+=(d2f33)*zbar[i]*w_bar_tilde[i];
00649   }
00650 
00651   for (size_t i=0;i<nvar;i++)
00652   {
00653     *x_tilde+=(d2f11)*zdotbar[i]*x_dot_bar_tilde[i];
00654     *y_tilde+=(d2f12)*zdotbar[i]*x_dot_bar_tilde[i];
00655     *w_tilde+=(d2f13)*zdotbar[i]*x_dot_bar_tilde[i];
00656 
00657     *x_tilde+=(d2f12)*zdotbar[i]*y_dot_bar_tilde[i];
00658     *y_tilde+=(d2f22)*zdotbar[i]*y_dot_bar_tilde[i];
00659     *w_tilde+=(d2f23)*zdotbar[i]*y_dot_bar_tilde[i];
00660 
00661     *x_tilde+=(d2f13)*zdotbar[i]*w_dot_bar_tilde[i];
00662     *y_tilde+=(d2f23)*zdotbar[i]*w_dot_bar_tilde[i];
00663     *w_tilde+=(d2f33)*zdotbar[i]*w_dot_bar_tilde[i];
00664   }
00665 
00666   for (size_t i=0;i<nvar;i++)
00667   {
00668     *x_tilde+=(d3f111)*xdot[i]*zdotbar[i]*x_bar_tilde[i];
00669     *y_tilde+=(d3f112)*xdot[i]*zdotbar[i]*x_bar_tilde[i];
00670     *w_tilde+=(d3f113)*xdot[i]*zdotbar[i]*x_bar_tilde[i];
00671 
00672     *x_tilde+=(d3f112)*xdot[i]*zdotbar[i]*y_bar_tilde[i];
00673     *y_tilde+=(d3f122)*xdot[i]*zdotbar[i]*y_bar_tilde[i];
00674     *w_tilde+=(d3f123)*xdot[i]*zdotbar[i]*y_bar_tilde[i];
00675 
00676     *x_tilde+=(d3f113)*xdot[i]*zdotbar[i]*w_bar_tilde[i];
00677     *y_tilde+=(d3f123)*xdot[i]*zdotbar[i]*w_bar_tilde[i];
00678     *w_tilde+=(d3f133)*xdot[i]*zdotbar[i]*w_bar_tilde[i];
00679 
00680     *x_tilde+=(d3f112)*ydot[i]*zdotbar[i]*x_bar_tilde[i];
00681     *y_tilde+=(d3f122)*ydot[i]*zdotbar[i]*x_bar_tilde[i];
00682     *w_tilde+=(d3f123)*ydot[i]*zdotbar[i]*x_bar_tilde[i];
00683 
00684     *x_tilde+=(d3f122)*ydot[i]*zdotbar[i]*y_bar_tilde[i];
00685     *y_tilde+=(d3f222)*ydot[i]*zdotbar[i]*y_bar_tilde[i];
00686     *w_tilde+=(d3f223)*ydot[i]*zdotbar[i]*y_bar_tilde[i];
00687 
00688     *x_tilde+=(d3f123)*ydot[i]*zdotbar[i]*w_bar_tilde[i];
00689     *y_tilde+=(d3f223)*ydot[i]*zdotbar[i]*w_bar_tilde[i];
00690     *w_tilde+=(d3f233)*ydot[i]*zdotbar[i]*w_bar_tilde[i];
00691 
00692     *x_tilde+=(d3f113)*wdot[i]*zdotbar[i]*x_bar_tilde[i];
00693     *y_tilde+=(d3f123)*wdot[i]*zdotbar[i]*x_bar_tilde[i];
00694     *w_tilde+=(d3f133)*wdot[i]*zdotbar[i]*x_bar_tilde[i];
00695 
00696     *x_tilde+=(d3f123)*wdot[i]*zdotbar[i]*y_bar_tilde[i];
00697     *y_tilde+=(d3f223)*wdot[i]*zdotbar[i]*y_bar_tilde[i];
00698     *w_tilde+=(d3f233)*wdot[i]*zdotbar[i]*y_bar_tilde[i];
00699 
00700     *x_tilde+=(d3f133)*wdot[i]*zdotbar[i]*w_bar_tilde[i];
00701     *y_tilde+=(d3f233)*wdot[i]*zdotbar[i]*w_bar_tilde[i];
00702     *w_tilde+=(d3f333)*wdot[i]*zdotbar[i]*w_bar_tilde[i];
00703   }
00704 
00705   for (size_t i=0;i<nvar;i++)
00706   {
00707     z_bar_tilde[i]+=(df1)*x_bar_tilde[i];
00708     z_bar_tilde[i]+=(df2)*y_bar_tilde[i];
00709     z_bar_tilde[i]+=(df3)*w_bar_tilde[i];
00710   }
00711 
00712   for (size_t i=0;i<nvar;i++)
00713   {
00714     z_dot_bar_tilde[i]+=(df1)*x_dot_bar_tilde[i];
00715     z_dot_bar_tilde[i]+=(df2)*y_dot_bar_tilde[i];
00716     z_dot_bar_tilde[i]+=(df3)*w_dot_bar_tilde[i];
00717   }
00718 
00719   for (size_t i=0;i<nvar;i++)
00720   {
00721     x_dot_tilde[i]+=(d2f11)*zdotbar[i]*x_bar_tilde[i];
00722     y_dot_tilde[i]+=(d2f12)*zdotbar[i]*x_bar_tilde[i];
00723     w_dot_tilde[i]+=(d2f13)*zdotbar[i]*x_bar_tilde[i];
00724 
00725     x_dot_tilde[i]+=(d2f12)*zdotbar[i]*y_bar_tilde[i];
00726     y_dot_tilde[i]+=(d2f22)*zdotbar[i]*y_bar_tilde[i];
00727     w_dot_tilde[i]+=(d2f23)*zdotbar[i]*y_bar_tilde[i];
00728 
00729     x_dot_tilde[i]+=(d2f13)*zdotbar[i]*w_bar_tilde[i];
00730     y_dot_tilde[i]+=(d2f23)*zdotbar[i]*w_bar_tilde[i];
00731     w_dot_tilde[i]+=(d2f33)*zdotbar[i]*w_bar_tilde[i];
00732   }
00733 
00734   for (size_t i=0;i<nvar;i++)
00735   {
00736     z_dot_bar_tilde[i]+=(d2f11)*xdot[i]*x_bar_tilde[i];
00737     z_dot_bar_tilde[i]+=(d2f12)*xdot[i]*y_bar_tilde[i];
00738     z_dot_bar_tilde[i]+=(d2f13)*xdot[i]*w_bar_tilde[i];
00739 
00740     z_dot_bar_tilde[i]+=(d2f12)*ydot[i]*x_bar_tilde[i];
00741     z_dot_bar_tilde[i]+=(d2f22)*ydot[i]*y_bar_tilde[i];
00742     z_dot_bar_tilde[i]+=(d2f23)*ydot[i]*w_bar_tilde[i];
00743 
00744     z_dot_bar_tilde[i]+=(d2f13)*wdot[i]*x_bar_tilde[i];
00745     z_dot_bar_tilde[i]+=(d2f23)*wdot[i]*y_bar_tilde[i];
00746     z_dot_bar_tilde[i]+=(d2f33)*wdot[i]*w_bar_tilde[i];
00747   }
00748 
00749 
00750 #if defined(__DERCHECK__)
00751   if (derchecker->node_number)
00752   {
00753     if (derchecker->counter == derchecker->node_number)
00754     {
00755       if (derchecker->pass_number==1) // increment the variable of interest
00756       {
00757         switch(derchecker->vartype)
00758         {
00759         case 1:
00760           if (!derchecker->dotflag)
00761             derchecker->der_value=
00762               px->u_bar_tilde[derchecker->index-1];
00763           else
00764             derchecker->der_value=
00765               px->u_dot_bar_tilde[derchecker->index-1];
00766           break;
00767         case 2:
00768           if (!derchecker->dotflag)
00769             derchecker->der_value=
00770               py->u_bar_tilde[derchecker->index-1];
00771           else
00772             derchecker->der_value=
00773               py->u_dot_bar_tilde[derchecker->index-1];
00774           break;
00775         case 3:
00776           if (!derchecker->dotflag)
00777             derchecker->der_value=
00778               pz->u_bar_tilde[derchecker->index-1];
00779           else
00780             derchecker->der_value=
00781               pz->u_dot_bar_tilde[derchecker->index-1];
00782           break;
00783         default:
00784           cerr << "Invalid index value for dercheck_index was "
00785                << derchecker->index << endl;
00786         }
00787       }
00788     }
00789   }
00790 #endif
00791 #if defined(PRINT_DERS)
00792  print_derivatives(pz,"z");
00793  print_derivatives(px,"x");
00794  print_derivatives(py,"y");
00795 #endif
00796 }
00797 
00802 void read_pass2_3_dvdvdv(void)
00803 {
00804   // We are going backword for bptr and forward for bptr2
00805   // the current entry+2 in bptr is the size of the record i.e
00806   // points to the next record
00807   unsigned int nvar=df1b2variable::nvar;
00808   fixed_smartlist & nlist=f1b2gradlist->nlist;
00809   test_smartlist& list=f1b2gradlist->list;
00810    // nlist-=sizeof(int);
00811   // get record size
00812   int num_bytes=nlist.bptr->numbytes;
00813   // backup the size of the record
00814   list-=num_bytes;
00815   list.saveposition(); // save pointer to beginning of record;
00816   // save the pointer to the beginning of the record
00817   double xu;
00818   double yu;
00819   double wu;
00820   //df1b2_header x,z;
00821   //df1b2function2 * pf;
00822 
00823   // get info from tape1
00824   // get info from tape1
00825 #if defined(SAFE_ARRAYS)
00826   checkidentiferstring("U8",list);
00827 #endif
00828 
00829   df1b2_header * px=(df1b2_header *) list.bptr;
00830   list.bptr+=sizeof(df1b2_header);
00831   df1b2_header * py=(df1b2_header *) list.bptr;
00832   list.bptr+=sizeof(df1b2_header);
00833   df1b2_header * pw=(df1b2_header *) list.bptr;
00834   list.bptr+=sizeof(df1b2_header);
00835   df1b2_header * pz=(df1b2_header *) list.bptr;
00836   list.bptr+=sizeof(df1b2_header);
00837 
00838   double df1=*(double*) list.bptr;
00839   list.bptr+=sizeof(double);
00840 
00841   double df2=*(double*) list.bptr;
00842   list.bptr+=sizeof(double);
00843 
00844   double df3=*(double*) list.bptr;
00845   list.bptr+=sizeof(double);
00846 
00847   double d2f11=*(double*) list.bptr;
00848   list.bptr+=sizeof(double);
00849 
00850   double d2f12=*(double*) list.bptr;
00851   list.bptr+=sizeof(double);
00852 
00853   double d2f13=*(double*) list.bptr;
00854   list.bptr+=sizeof(double);
00855 
00856   double d2f22=*(double*) list.bptr;
00857   list.bptr+=sizeof(double);
00858 
00859   double d2f23=*(double*) list.bptr;
00860   list.bptr+=sizeof(double);
00861 
00862   double d2f33=*(double*) list.bptr;
00863   list.bptr+=sizeof(double);
00864 
00865 #if defined(PRINT_DERS)
00866   double d3f111=*(double*) list.bptr;
00867 #endif
00868   list.bptr+=sizeof(double);
00869 
00870 #if defined(PRINT_DERS)
00871   double d3f112=*(double*) list.bptr;
00872 #endif
00873   list.bptr+=sizeof(double);
00874 
00875 #if defined(PRINT_DERS)
00876   double d3f113=*(double*) list.bptr;
00877 #endif
00878   list.bptr+=sizeof(double);
00879 
00880 #if defined(PRINT_DERS)
00881   double d3f122=*(double*) list.bptr;
00882 #endif
00883   list.bptr+=sizeof(double);
00884 
00885 #if defined(PRINT_DERS)
00886   double d3f123=*(double*) list.bptr;
00887 #endif
00888   list.bptr+=sizeof(double);
00889 
00890 #if defined(PRINT_DERS)
00891   double d3f133=*(double*) list.bptr;
00892 #endif
00893   list.bptr+=sizeof(double);
00894 
00895 #if defined(PRINT_DERS)
00896   double d3f222=*(double*) list.bptr;
00897 #endif
00898   list.bptr+=sizeof(double);
00899 
00900 #if defined(PRINT_DERS)
00901   double d3f223=*(double*) list.bptr;
00902 #endif
00903   list.bptr+=sizeof(double);
00904 
00905 #if defined(PRINT_DERS)
00906   double d3f233=*(double*) list.bptr;
00907 #endif
00908   list.bptr+=sizeof(double);
00909 
00910 #if defined(PRINT_DERS)
00911   double d3f333=*(double*) list.bptr;
00912 #endif
00913   list.bptr+=sizeof(double);
00914 
00915   memcpy(&xu,list.bptr,sizeof(double));
00916   list.bptr+=sizeof(double);
00917   memcpy(&yu,list.bptr,sizeof(double));
00918   list.bptr+=sizeof(double);
00919   memcpy(&wu,list.bptr,sizeof(double));
00920   list.bptr+=sizeof(double);
00921   double * xdot=(double*)list.bptr;
00922   list.bptr+=nvar*sizeof(double);
00923   double * ydot=(double*)list.bptr;
00924   list.bptr+=nvar*sizeof(double);
00925   double * wdot=(double*)list.bptr;
00926 
00927   list.restoreposition(); // save pointer to beginning of record;
00928 #if defined(PRINT_DERS)
00929  print_derivatives(funname,(f),(df1),
00930   (df2),(d2f11),(d2f12),(d2f22),
00931   (d3f111),(d3f112),(d3f122),
00932   (d3f222),1);
00933  print_derivatives(pz,"z");
00934  print_derivatives(px,"x");
00935  print_derivatives(py,"y");
00936  print_derivatives(pw,"y");
00937 #endif
00938 
00939   *(px->u_tilde)+=(df1)* *(pz->u_tilde);
00940   *(py->u_tilde)+=(df2)* *(pz->u_tilde);
00941   *(pw->u_tilde)+=(df3)* *(pz->u_tilde);
00942 
00943   for (size_t i=0;i<nvar;i++)
00944   {
00945     *(px->u_tilde)+=(d2f11)*xdot[i]*pz->u_dot_tilde[i];
00946     *(py->u_tilde)+=(d2f12)*xdot[i]*pz->u_dot_tilde[i];
00947     *(pw->u_tilde)+=(d2f13)*xdot[i]*pz->u_dot_tilde[i];
00948 
00949     *(px->u_tilde)+=(d2f12)*ydot[i]*pz->u_dot_tilde[i];
00950     *(py->u_tilde)+=(d2f22)*ydot[i]*pz->u_dot_tilde[i];
00951     *(pw->u_tilde)+=(d2f23)*ydot[i]*pz->u_dot_tilde[i];
00952 
00953     *(px->u_tilde)+=(d2f13)*wdot[i]*pz->u_dot_tilde[i];
00954     *(py->u_tilde)+=(d2f23)*wdot[i]*pz->u_dot_tilde[i];
00955     *(pw->u_tilde)+=(d2f33)*wdot[i]*pz->u_dot_tilde[i];
00956   }
00957   for (size_t i=0;i<nvar;i++)
00958   {
00959     px->u_dot_tilde[i]+=(df1)*pz->u_dot_tilde[i];
00960     py->u_dot_tilde[i]+=(df2)*pz->u_dot_tilde[i];
00961     pw->u_dot_tilde[i]+=(df3)*pz->u_dot_tilde[i];
00962   }
00963   *(pz->u_tilde)=0;
00964   for (size_t i=0;i<nvar;i++)
00965   {
00966     pz->u_dot_tilde[i]=0;
00967   }
00968 
00969 
00970 #if defined(PRINT_DERS)
00971  print_derivatives(pz,"z");
00972  print_derivatives(px,"x");
00973  print_derivatives(py,"y");
00974  print_derivatives(pw,"y");
00975 #endif
00976 }