ADMB Documentation  11.2.2828
 All Classes Files Functions Variables Typedefs Friends Defines
jacob3.cpp
Go to the documentation of this file.
00001 /*
00002  * $Id: jacob3.cpp 2643 2014-11-13 21:00:02Z johnoel $
00003  *
00004  * Author: David Fournier
00005  * Copyright (c) 2008-2012 Regents of the University of California
00006  */
00011 #include "fvar.hpp"
00012 
00013 #include <sys/stat.h>
00014 #include <fcntl.h>
00015 #include <string.h>
00016 
00017 #ifdef __TURBOC__
00018   #pragma hdrstop
00019   #include <iostream.h>
00020 #endif
00021 
00022 #include <stdio.h>
00023 #include <stdlib.h>
00024 
00025 #ifdef __SUN__
00026   #include <iostream.h>
00027   #include <sys/stat.h>
00028   #include <sys/types.h>
00029   #include <unistd.h>
00030   #ifdef _MSC_VER
00031     #define lseek _lseek
00032     #define  read _read
00033     #define write _write
00034     #define open _open
00035     #define close _close
00036   #endif
00037 #endif
00038 
00039 #if !defined(_MSC_VER)
00040   #include <iostream>
00041   using namespace std;
00042   #include <sys/stat.h>
00043   #include <sys/types.h>
00044   #include <unistd.h>
00045 #endif
00046 
00047 #if defined(__NDPX__ )
00048   extern "C" {
00049     int lseek(int, int, int);
00050     int read(int, char*, int);
00051   };
00052 #endif
00053 
00054 #include <math.h>
00055 
00060 void jacobcalc(int nvar, const uostream& ofs)
00061 {
00062   gradient_structure::jacobcalc(nvar,ofs);
00063 }
00064 
00069 void gradient_structure::jacobcalc(int nvar, const uostream& ofs)
00070 {
00071   dvector jac(1,nvar);
00072   off_t lpos;
00073   int depvar_count=DEPVARS_INFO->depvar_count;
00074 
00075   int& _GRADFILE_PTR=GRAD_STACK1->_GRADFILE_PTR;
00076   // check to see if anything has been written into the file
00077   off_t last_gpos=lseek(_GRADFILE_PTR,0L,SEEK_CUR);
00078 
00079   //save current contents of the buffer so we can get them later
00080   if (last_gpos)
00081   {
00082     GRAD_STACK1->write_grad_stack_buffer();
00083   }
00084 
00085   // check to see if anything has been written into the file
00086   off_t last_cpos=lseek(fp->file_ptr,0L,SEEK_CUR);
00087 
00088   //save current contents of the buffer so we can get them later
00089   if (last_cpos)
00090   {
00091     fp->write_cmpdif_stack_buffer();
00092   }
00093 
00094   // save variable values if desired
00095   if (save_var_flag)
00096   {
00097     save_arrays();
00098     save_variables();
00099   }
00100   // now evalueate the jacobian
00101   for (int ijac=1;ijac<=depvar_count;ijac++)
00102   {
00103     dvector& g=jac;
00104     //max_num_dependent_variables=ndv;
00105     if (depvar_count>DEPVARS_INFO->max_num_dependent_variables)
00106     {
00107       cout << "maximum number of depdendent variables of "
00108          << DEPVARS_INFO->max_num_dependent_variables << " exceeded "
00109          << endl
00110          << "use gradient_structure::set_NUM_DEPENDENT_VARIABLES(int i);"
00111          << endl << "to increase the number of dependent variables"
00112          << endl;
00113       ad_exit(1);
00114     }
00115 
00116     fp->offset=DEPVARS_INFO->cmpdif_buffer_position(ijac);
00117     fp->toffset=fp->offset;
00118     _GRADFILE_PTR=DEPVARS_INFO->grad_file_count(ijac);
00119     fp->file_ptr=DEPVARS_INFO->cmpdif_file_count(ijac);
00120     lpos=DEPVARS_INFO->grad_file_position(ijac);
00121     // position the cmpdif file correctly;
00122     if (last_cpos)
00123     {
00124       off_t cmp_lpos=DEPVARS_INFO->cmpdif_file_position(ijac);
00125       lseek(fp->file_ptr,cmp_lpos,SEEK_SET);
00126       fp->read_cmpdif_stack_buffer(cmp_lpos);
00127     }
00128     GRAD_STACK1->_GRADFILE_PTR = GRAD_STACK1->gradfile_handle();
00129 
00130     if (last_gpos)
00131     {
00132       // just use the end of the buffer
00133       GRAD_STACK1->set_gbuffer_pointers();
00134 
00135       // check to sere if buffer was written into the beginning of
00136       // the next file
00137       if ( (GRAD_STACK1->_GRADFILE_PTR == GRAD_STACK1->_GRADFILE_PTR1)
00138          && (lpos == GRAD_STACK1->end_pos1) && (lpos>0) )
00139       {
00140         // get the next file
00141         GRAD_STACK1->increment_current_gradfile_ptr();
00142         lpos=0;
00143       }
00144       // and position the file to the begginig of the buffer image
00145       lseek(_GRADFILE_PTR,lpos,SEEK_SET);
00146       // now fill the buffer with the appropriate stuff
00147       GRAD_STACK1->read_grad_stack_buffer(lpos);
00148       // now reposition the grad_buffer pointer
00149     }
00150     GRAD_STACK1->ptr=
00151          (grad_stack_entry *)DEPVARS_INFO->grad_buffer_position(ijac);
00152 
00153     if(GRAD_STACK1->ptr <= GRAD_STACK1->ptr_first)
00154     {
00155 #ifdef SAFE_ALL
00156         cerr << "warning -- calling gradcalc when no calculations generating"
00157          << endl << "derivative information have occurred" << endl;
00158 #endif
00159       g.initialize();
00160       return;
00161     }    // current is one past the end so -- it
00162 
00163     gradient_structure::GRAD_STACK1->ptr--;
00164 
00165     for (unsigned int i=0; i< GRAD_LIST->nlinks; i++)
00166     {
00167       * (double*) (GRAD_LIST->dlink_addresses[i]) = 0;
00168     }
00169 
00170     double_and_int* tmp =
00171       (double_and_int*)gradient_structure::ARRAY_MEMBLOCK_BASE;
00172 
00173     unsigned long int max_last_offset
00174                = gradient_structure::ARR_LIST1->get_max_last_offset();
00175 
00176     unsigned int size = sizeof(double_and_int );
00177 
00178     for (unsigned int i = 0; i < (max_last_offset/size); i++)
00179     {
00180       tmp->x = 0;
00181 #if defined (__ZTC__)
00182       tmp = (double_and_int*)_farptr_norm((void*)(++tmp));
00183 #else
00184       tmp++;
00185 #endif
00186     }
00187 
00188     * gradient_structure::GRAD_STACK1->ptr->dep_addr  = 1;
00189     //double* zptr = gradient_structure::GRAD_STACK1->ptr->dep_addr;
00190 
00191     int break_flag=1;
00192     do
00193     {
00194       gradient_structure::GRAD_STACK1->ptr++;
00195       #ifdef FAST_ASSEMBLER
00196         gradloop();
00197       #else
00198         //int counter=0;
00199       while (gradient_structure::GRAD_STACK1->ptr-- >
00200              gradient_structure::GRAD_STACK1->ptr_first)
00201       {
00202         //grad_stack_entry* grad_ptr =
00203         //gradient_structure::GRAD_STACK1->ptr;
00204         {
00205           (* gradient_structure::GRAD_STACK1->ptr->func)();
00206         }
00207       }
00208       #endif
00209 
00210       // \todo Need test
00211       // back up the file one buffer size and read forward
00212       off_t offset = (off_t)(sizeof(grad_stack_entry)
00213         * gradient_structure::GRAD_STACK1->length);
00214       lpos = lseek(gradient_structure::GRAD_STACK1->_GRADFILE_PTR,
00215         -offset, SEEK_CUR);
00216 
00217        break_flag=gradient_structure::
00218                   GRAD_STACK1->read_grad_stack_buffer(lpos);
00219     }  while (break_flag); // do
00220 
00221     int mindx = g.indexmin();
00222     for (int i=0; i < nvar; i++)
00223     {
00224       g[i+mindx] =  * gradient_structure::INDVAR_LIST->get_address(i);
00225     }
00226     gradient_structure::GRAD_STACK1->ptr =
00227          gradient_structure::GRAD_STACK1->ptr_first;
00228     //ofs << setprecision(10) << g << endl;
00229     ofs << g ;
00230   }// loop over dep vars
00231   DEPVARS_INFO->depvar_count=0;
00232   if (gradient_structure::save_var_flag)
00233   {
00234     gradient_structure::restore_arrays();
00235     gradient_structure::restore_variables();
00236   }
00237 }