ADMB Documentation  11.1.2192
 All Classes Files Functions Variables Typedefs Friends Defines
sgradclc.cpp
Go to the documentation of this file.
00001 /*
00002  * $Id: sgradclc.cpp 1941 2014-04-28 22:18:29Z johnoel $
00003  *
00004  * Author: David Fournier
00005  * Copyright (c) 2008-2012 Regents of the University of California
00006  */
00011 #if !defined(DOS386)
00012   #define DOS386
00013 #endif
00014 
00015 #include "fvar.hpp"
00016 
00017 #include <sys/stat.h>
00018 #include <fcntl.h>
00019 #include <string.h>
00020 
00021 #ifdef __TURBOC__
00022   #pragma hdrstop
00023   #include <iostream.h>
00024 #endif
00025 
00026 #ifdef __ZTC__
00027   #include <iostream.hpp>
00028 #endif
00029 
00030 #if defined (__WAT32__)
00031 #  include <io.h>
00032 #endif
00033 
00034 #include <stdio.h>
00035 #include <stdlib.h>
00036 
00037 #ifdef _MSC_VER
00038   #define lseek _lseek
00039   #define  read _read
00040   #define write _write
00041 #else
00042   #include <iostream>
00043   using namespace std;
00044   #include <sys/stat.h>
00045   #include <sys/types.h>
00046   #include <unistd.h>
00047 #endif
00048 
00049 #ifdef __SUN__
00050   #include <iostream.h>
00051   #include <sys/stat.h>
00052   #include <sys/types.h>
00053   #include <unistd.h>
00054 #endif
00055 
00056 #if defined(__NDPX__ )
00057   extern "C" {
00058     int lseek(int, int, int);
00059     int read(int, char*, int);
00060   };
00061 #endif
00062 
00063 #include <math.h>
00064 
00065 #if defined(__ZTC__)
00066   void _far * _cdecl _farptr_norm(void _far *);
00067   void _far * _cdecl _farptr_fromlong(unsigned long);
00068   long _cdecl _farptr_tolong(void _far *);
00069 #endif
00070 
00077 void gradcalc(int nvar, const dvector& _g)
00078 {
00079   if (nvar!=0)
00080   {
00081     if (nvar != gradient_structure::NVAR)
00082     {
00083       cerr << "nvar != gradient_structure::NVAR in gradcalc" << endl;
00084       cerr << "  nvar = " << nvar << endl;
00085       cerr << "  gradient_structure::NVAR = " << gradient_structure::NVAR
00086            << endl;
00087       cerr << "  in " __FILE__ << endl;
00088       ad_exit(1);
00089     }
00090   }
00091   dvector& g= (dvector&) _g;
00092   gradient_structure::TOTAL_BYTES = 0;
00093   gradient_structure::PREVIOUS_TOTAL_BYTES=0;
00094   unsigned int i;
00095   if(!gradient_structure::instances)
00096   {
00097     g.initialize();
00098     return;
00099   }
00100 
00101   if (g.size() < nvar)
00102   {
00103     cerr  << "gradient vector size is less than the number of variables.\n";
00104     ad_exit(1);
00105   }
00106 
00107    gradient_structure::GRAD_STACK1->_GRADFILE_PTR =
00108      gradient_structure::GRAD_STACK1->gradfile_handle();
00109 
00110   int& _GRADFILE_PTR=gradient_structure::GRAD_STACK1->_GRADFILE_PTR;
00111 
00112   long int lpos = lseek(_GRADFILE_PTR,0L,SEEK_CUR);
00113 
00114   if (gradient_structure::GRAD_STACK1->ptr <=
00115         gradient_structure::GRAD_STACK1->ptr_first)
00116   {
00117 #ifdef SAFE_ALL
00118     cerr << "warning -- calling gradcalc when no calculations generating"
00119          << endl << "derivative information have occurred" << endl;
00120 #endif
00121     g.initialize();
00122     return;
00123   }    // current is one past the end so -- it
00124 
00125   if (gradient_structure::save_var_flag)
00126   {
00127     gradient_structure::save_arrays();
00128     gradient_structure::save_variables();
00129   }
00130 
00131   gradient_structure::GRAD_STACK1->ptr--;
00132 
00133   for (i = 0; i < gradient_structure::GRAD_LIST->nlinks; i++)
00134   {
00135     *(double*)(gradient_structure::GRAD_LIST->dlink_addresses[i]) = 0;
00136   }
00137 
00138   double_and_int* tmp =
00139     (double_and_int*)gradient_structure::ARRAY_MEMBLOCK_BASE;
00140 
00141   unsigned long int max_last_offset =
00142     gradient_structure::ARR_LIST1->get_max_last_offset();
00143   unsigned int size = sizeof(double_and_int);
00144   for (i = 0; i < (max_last_offset/size); i++)
00145   {
00146      tmp->x = 0;
00147 #if defined (__ZTC__)
00148      tmp = (double_and_int*)_farptr_norm((void*)(++tmp));
00149 #else
00150      tmp++;
00151 #endif
00152   }
00153 
00154   *gradient_structure::GRAD_STACK1->ptr->dep_addr = 1;
00155 
00156   //int icount=0;
00157   int break_flag=1;
00158   do
00159   {
00160     gradient_structure::GRAD_STACK1->ptr++;
00161 #ifdef FAST_ASSEMBLER
00162     gradloop();
00163 #else
00164     grad_stack_entry* grad_ptr_first =
00165       gradient_structure::GRAD_STACK1->ptr_first;
00166     while (gradient_structure::GRAD_STACK1->ptr-- >
00167              grad_ptr_first)
00168     {
00169       (*(gradient_structure::GRAD_STACK1->ptr->func))();
00170 /*
00171       icount++;
00172       if (icount%1000==0)
00173       {
00174         //cout << "icount = " << icount << endl;
00175       }
00176 */
00177     }
00178 #endif
00179 
00180    // back up the file one buffer size and read forward
00181    lpos = lseek(gradient_structure::GRAD_STACK1->_GRADFILE_PTR,
00182       -((long int)(sizeof(grad_stack_entry)*gradient_structure::
00183         GRAD_STACK1->length)),SEEK_CUR);
00184 
00185     break_flag=gradient_structure::GRAD_STACK1->read_grad_stack_buffer(lpos);
00186   } while (break_flag);
00187 
00188   {
00189 #ifdef GRAD_DIAG
00190     long int ttmp =
00191 #endif
00192     lseek(gradient_structure::GRAD_STACK1->_GRADFILE_PTR, 0,SEEK_CUR);
00193 #ifdef GRAD_DIAG
00194     cout << "Offset in file at end of gradcalc is " << ttmp
00195          << " bytes from the beginning\n";
00196 #endif
00197   }
00198 
00199   int mindx = g.indexmin();
00200   for (i=0; i < (unsigned int)nvar; i++)
00201   {
00202     g[i + mindx] = *gradient_structure::INDVAR_LIST->get_address(i);
00203   }
00204 
00205   gradient_structure::GRAD_STACK1->ptr =
00206     gradient_structure::GRAD_STACK1->ptr_first;
00207 
00208   if (gradient_structure::save_var_flag)
00209   {
00210     gradient_structure::restore_arrays();
00211     gradient_structure::restore_variables();
00212   }
00213 }
00222 double gradcalc(int nvar, const dvector& _g, dvariable& f)
00223 {
00224   double v = value(f);
00225   gradcalc(nvar, _g);
00226   return v;
00227 }
00230 void gradient_structure::save_arrays()
00231 {
00232   void * temp_ptr;
00233   long bytes_needed=min(gradient_structure::ARR_LIST1->get_last_offset()+1,
00234     ARRAY_MEMBLOCK_SIZE);
00235   gradient_structure::save_var_file_flag=0;
00236 #ifdef __ZTC__
00237    if ( (temp_ptr = farmalloc(bytes_needed) ) == 0)
00238 #else
00239    //if ( (temp_ptr = malloc(bytes_needed) ) == 0)
00240    if ((temp_ptr = (void *)malloc(bytes_needed)) == 0)
00241   #define __USE_IOSTREAM__
00242 #endif
00243    {
00244      gradient_structure::save_var_file_flag=1;
00245      cerr << "insufficient memory to allocate space for ARRAY_MEMBLOCK"
00246           << " save buffer " << endl;
00247    }
00248    if (gradient_structure::save_var_file_flag==0)
00249    {
00250      ARRAY_MEMBLOCK_SAVE = temp_ptr;
00251 #if defined(DOS386)
00252   #ifndef USE_ASSEMBLER
00253          memcpy((char*)ARRAY_MEMBLOCK_SAVE,(char*)ARRAY_MEMBLOCK_BASE,
00254            bytes_needed);
00255   #else
00256          dw_block_move((double*)ARRAY_MEMBLOCK_SAVE,
00257            (double*)ARRAY_MEMBLOCK_BASE,bytes_needed/8);
00258   #endif
00259 #else
00260      unsigned long int max_move=50000;
00261      unsigned long int left_to_move=bytes_needed;
00262      humungous_pointer dest = ARRAY_MEMBLOCK_SAVE;
00263      humungous_pointer src = ARRAY_MEMBLOCK_BASE;
00264      while(left_to_move > max_move)
00265      {
00266        memcpy((char*)dest,(char*)src,max_move);
00267        left_to_move-=max_move;
00268        src+=max_move;
00269        dest+=max_move;
00270      }
00271      memcpy((char*)dest,(char*)src,left_to_move);
00272 #endif
00273   }
00274   else
00275   {
00276      humungous_pointer src = ARRAY_MEMBLOCK_BASE;
00277      lseek(gradient_structure::GRAD_STACK1->_VARSSAV_PTR,0L,SEEK_SET);
00278 #if defined(DOS386)
00279        write(gradient_structure::GRAD_STACK1->_VARSSAV_PTR,
00280          (char*)src,bytes_needed);
00281 #else
00282      unsigned long int max_move=500;
00283      unsigned long int left_to_move=bytes_needed;
00284      while(left_to_move > max_move)
00285      {
00286        write(_VARSSAV_PTR,(char*)src,max_move);
00287        left_to_move-=max_move;
00288        src+=max_move;
00289      }
00290      write(gradient_structure::GRAD_STACK1->_VARSSAV_PTR,(char*)src,
00291        left_to_move);
00292 #endif
00293   }
00294 }
00295 
00298 void gradient_structure::restore_arrays()
00299 {
00300   long bytes_needed=min(gradient_structure::ARR_LIST1->get_last_offset()+1,
00301     ARRAY_MEMBLOCK_SIZE);
00302   if (gradient_structure::save_var_file_flag==0)
00303   {
00304 #if defined(DOS386)
00305   #ifndef USE_ASSEMBLER
00306         memcpy((char*)ARRAY_MEMBLOCK_BASE,(char*)ARRAY_MEMBLOCK_SAVE,
00307           bytes_needed);
00308   #else
00309          dw_block_move((double*)ARRAY_MEMBLOCK_BASE,
00310            (double*)ARRAY_MEMBLOCK_SAVE,bytes_needed/8);
00311   #endif
00312 #else
00313      unsigned long max_move=50000;
00314 
00315      long int left_to_move=bytes_needed;
00316      humungous_pointer src = ARRAY_MEMBLOCK_SAVE;
00317      humungous_pointer dest = ARRAY_MEMBLOCK_BASE;
00318      while(left_to_move > max_move)
00319      {
00320        memcpy((char*)dest,(char*)src,max_move);
00321        left_to_move-=max_move;
00322        src+=max_move;
00323        dest+=max_move;
00324      }
00325      memcpy((char*)dest,(char*)src,left_to_move);
00326 #endif
00327     ARRAY_MEMBLOCK_SAVE.free();
00328   }
00329   else
00330   {
00331     humungous_pointer dest = ARRAY_MEMBLOCK_BASE;
00332     lseek(gradient_structure::GRAD_STACK1->_VARSSAV_PTR,0L,SEEK_SET);
00333 #if defined(DOS386)
00334       read(gradient_structure::GRAD_STACK1->_VARSSAV_PTR,
00335        (char*)dest,bytes_needed);
00336 #else
00337      unsigned long int max_move=50000;
00338 
00339      long int left_to_move=bytes_needed;
00340      while(left_to_move > max_move)
00341      {
00342        read(gradient_structure::GRAD_STACK1->_VARSSAV_PTR,
00343          (char*)dest,max_move);
00344        left_to_move-=max_move;
00345        dest+=max_move;
00346      }
00347      read(gradient_structure::GRAD_STACK1->_VARSSAV_PTR,
00348        (char*)dest,left_to_move);
00349 #endif
00350   }
00351 }
00355 void gradient_structure::save_variables()
00356 {
00357   if ((variables_save = new double[gradient_structure::MAX_DLINKS])==NULL)
00358   {
00359     //_VARSSAV_PTR=open(var_store_file_name, O_RDWR | O_CREAT | O_TRUNC ,
00360     //               S_IREAD | S_IWRITE);
00361     cerr << "insufficient memory to allocate space for dvariables"
00362          << " save buffer " << endl;
00363     ad_exit(1);
00364   }
00365   for (unsigned int i=0; i<gradient_structure::GRAD_LIST->nlinks; i++)
00366   {
00367     //variables_save[i]= * (double*)
00368     //    (gradient_structure::GRAD_LIST->dlink_addresses[i]);
00369     memcpy(&(variables_save[i]),
00370        gradient_structure::GRAD_LIST->dlink_addresses[i],sizeof(double));
00371   }
00372 }
00376 void gradient_structure::restore_variables()
00377 {
00378   for (unsigned int i=0; i<gradient_structure::GRAD_LIST->nlinks; i++)
00379   {
00380     //* (double*)(gradient_structure::GRAD_LIST->dlink_addresses[i])
00381     //  = variables_save[i];
00382     memcpy(gradient_structure::GRAD_LIST->dlink_addresses[i],
00383        &(variables_save[i]),sizeof(double));
00384   }
00385   delete [] variables_save;
00386   variables_save = 0;
00387 }
00391 void reset_gradient_stack(void)
00392 {
00393   gradient_structure::GRAD_STACK1->ptr =
00394     gradient_structure::GRAD_STACK1->ptr_first;
00395 
00396   int& _GRADFILE_PTR=gradient_structure::GRAD_STACK1->_GRADFILE_PTR;
00397 
00398   lseek(_GRADFILE_PTR,0L,SEEK_SET);
00399 }
00409 void grad_stack::set_gradient_stack1(void (* func)(void),
00410   double* dep_addr, double* ind_addr1)
00411 {
00412 #ifdef NO_DERIVS
00413   if (!gradient_structure::no_derivatives)
00414   {
00415 #endif
00416     if (ptr > ptr_last)
00417     {
00418        // current buffer is full -- write it to disk and reset pointer
00419        // and counter
00420        this->write_grad_stack_buffer();
00421     }
00422     //test_the_pointer();
00423     ptr->func = func;
00424     ptr->dep_addr = dep_addr;
00425     ptr->ind_addr1 = ind_addr1;
00426     ptr++;
00427 #ifdef NO_DERIVS
00428   }
00429 #endif
00430 }
00431 void test_the_pointer(void)
00432 {
00433 }