ADMB Documentation  11.4.2891
 All Classes Files Functions Variables Typedefs Friends Defines
admodel.h
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  *
00007  * ADModelbuilder and associated libraries and documentations are
00008  * provided under the general terms of the "BSD" license.
00009  *
00010  * License:
00011  *
00012  * Redistribution and use in source and binary forms, with or without
00013  * modification, are permitted provided that the following conditions are
00014  * met:
00015  *
00016  * 1. Redistributions of source code must retain the above copyright
00017  * notice, this list of conditions and the following disclaimer.
00018  *
00019  * 2.  Redistributions in binary form must reproduce the above copyright
00020  * notice, this list of conditions and the following disclaimer in the
00021  * documentation and/or other materials provided with the distribution.
00022  *
00023  * 3.  Neither the name of the  University of California, Otter Research,
00024  * nor the ADMB Foundation nor the names of its contributors may be used
00025  * to endorse or promote products derived from this software without
00026  * specific prior written permission.
00027  *
00028  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
00029  * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
00030  * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
00031  * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
00032  * OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
00033  * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
00034  * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
00035  * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
00036  * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
00037  * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
00038  * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00039  *
00040  */
00045 #ifndef _ADMODEL_H_
00046 #define _ADMODEL_H_
00047 
00048 #define USE_SHARE_FLAGS
00049 //#define DO_PROFILE
00050 #define __MINI_MAX__
00051 #if defined(__GNUC__) && (__GNUC__ < 3)
00052   #pragma interface
00053 #endif
00054 
00055 #define BIG_INIT_PARAMS
00056 
00057   class laplace_approximation_calculator;
00058   void cleanup_laplace_stuff(laplace_approximation_calculator *);
00059 
00060 #include <fvar.hpp>
00061 
00062 //#include <d4arr.hpp>
00063 #include <cifstrem.h>
00064 
00065 #include <adstring.hpp>
00066 class init_xml_doc;
00067 
00068 //#define _ADSTD_ std::
00069 #define param_matrix named_dvar_matrix
00070 #define param_vector named_dvar_vector
00071 #define param_number named_dvariable
00072 #define param_3array named_dvar3_array
00073 #define param_4array named_dvar4_array
00074 #define param_5array named_dvar5_array
00075 #define param_6array named_dvar6_array
00076 #define param_7array named_dvar7_array
00077 
00078 #define SPparam_matrix SPnamed_dvar_matrix
00079 #define SPparam_vector SPnamed_dvar_vector
00080 #define SPparam_number SPnamed_dvariable
00081 #define SPparam_3array SPnamed_dvar3_array
00082 #define SPparam_4array SPnamed_dvar4_array
00083   double mfexp(const double);
00084   dvariable mfexp(const prevariable& v1);
00085   dvar_vector mfexp(const dvar_vector&);
00086   dvector mfexp(const dvector&);
00087 
00088   class param_init_bounded_number_vector;
00089   class model_parameters;
00090   extern int AD_gaussflag;
00091 
00092   extern function_minimizer * pfm;
00093   extern int traceflag;
00094   extern int ADqd_flag;
00095 /*
00096  void set_value_inv(const dvariable& x, const dvector& v, const int& ii);
00097  void set_value_inv(const dvar_matrix& x, const dvector& v, const int& ii);
00098  void set_value_inv(const dvar_vector&, const dvector&, const int &ii);
00099  void set_value_inv(const dvariable& u, const dvector& x, const int& ii,
00100    const double fmin, const double fmax);
00101  void set_value_inv(const dvar_matrix& u, const dvector& x, const int& ii,
00102   const double fmin,const double fmax);
00103  void set_value_inv(const dvar3_array& u, const dvector& x, const int& ii,
00104   const double fmin,const double fmax);
00105  void set_value_inv(const dvar3_array& u, const dvector& x, const int& ii);
00106 */
00107 
00108   void copy_value_to_vector(const prevariable& x, const dvector& v,
00109     const int& ii);
00110   void copy_value_to_vector(const dvar_vector& x, const dvector& v,
00111     const int& ii);
00112   void copy_value_to_vector(const dvar_matrix& x, const dvector& v,
00113     const int& ii);
00114   void copy_value_to_vector(const dvar3_array& x, const dvector& v,
00115     const int& ii);
00116 
00117   void restore_value_from_vector(const prevariable& x, const dvector& v,
00118     const int& ii);
00119   void restore_value_from_vector(const dvar_vector& x, const dvector& v,
00120     const int& ii);
00121   void restore_value_from_vector(const dvar_matrix& x, const dvector& v,
00122     const int& ii);
00123   void restore_value_from_vector(dvar3_array& x, const dvector& v,
00124     const int& ii);
00125 
00130 class AD_matherror
00131 {
00132 public:
00133 #if defined (_BORLANDC_)
00134   exception * err;
00135   AD_matherror(exception * _err) : err(_err) {;}
00136 #endif
00137 #if defined (_MSC_VER)
00138   _exception * err;
00139   AD_matherror(_exception * _err) : err(_err) {;}
00140 #endif
00141 };
00142 
00143 class model_data;
00144 
00149 class label_class
00150 {
00151   const char * name;
00152   friend ostream& operator<<(const ostream& s, const label_class& lc);
00153   friend class model_name_tag;
00154 public:
00155   const char * mychar(void) { return name;}
00156   label_class(const char * s){name=s;}
00157 };
00158 
00163 class model_name_tag
00164 {
00165 protected:
00166   adstring name;
00167   //friend ostream& operator<<(const ostream& os, const model_name_tag& mnt);
00168 public:
00169   model_name_tag(void){;}
00170   void allocate(const char * s);
00171   label_class label(void){return (char*)(name);}
00172   const char * get_name(void) { return name;}
00173 };
00174 
00179 class named_dvar_vector : public dvar_vector, public model_name_tag
00180 {
00181 protected:
00182   named_dvar_vector(void) : dvar_vector(), model_name_tag() {;}
00183   named_dvar_vector& operator=(const dvar_vector& m);
00184   named_dvar_vector& operator=(const dvector& m);
00185   named_dvar_vector& operator=(const double m);
00186   named_dvar_vector& operator=(const prevariable& m);
00187   friend class model_parameters;
00188   void allocate(int mmin,int mmax,const char * s);
00189   void allocate(const char * s);
00190 };
00191 
00196 class equality_constraint_vector : public named_dvar_vector
00197 {
00198 protected:
00199   equality_constraint_vector(void);
00200   equality_constraint_vector& operator=(const dvar_vector& m);
00201   equality_constraint_vector& operator=(const dvector& m);
00202   equality_constraint_vector& operator=(const double m);
00203   equality_constraint_vector& operator=(const prevariable& m);
00204   friend class model_parameters;
00205   void allocate(int mmin,int mmax,const char * s);
00206 };
00207 
00212 class inequality_constraint_vector : public named_dvar_vector
00213 {
00214 protected:
00215   inequality_constraint_vector(void);
00216   inequality_constraint_vector& operator=(const dvar_vector& m);
00217   inequality_constraint_vector& operator=(const dvector& m);
00218   inequality_constraint_vector& operator=(const double m);
00219   inequality_constraint_vector& operator=(const prevariable& m);
00220   friend class model_parameters;
00221   void allocate(int mmin,int mmax,const char * s);
00222 };
00223 
00228 class dll_param_vector : public named_dvar_vector
00229 {
00230   double * pd;
00231 public:
00232   ~dll_param_vector();
00233   void allocate(double *_pd,int mmin,int mmax,const char * s);
00234   dll_param_vector& operator=(const dvar_vector& m);
00235   dll_param_vector& operator=(const dvector& m);
00236   dll_param_vector& operator=(const double m);
00237   dll_param_vector& operator=(const prevariable& m);
00238 };
00239 
00244 class named_dvariable : public dvariable, public model_name_tag
00245 {
00246   //named_dvariable& operator=(const dvariable& m);
00247 protected:
00248   named_dvariable(void) : dvariable(), model_name_tag() {;}
00249   void allocate(const char * s);
00250   named_dvariable& operator=(const prevariable& m);
00251   named_dvariable& operator=(const double m);
00252   friend class model_parameters;
00253 };
00254 
00259 class dll_param_number : public named_dvariable
00260 {
00261   double * pd;
00262 protected:
00263   //named_dvariable(void) : dvariable(), model_name_tag() {;}
00264   void allocate(double *_d,const char * s);
00265   dll_param_number& operator=(const prevariable& m);
00266   dll_param_number& operator=(const double m);
00267   virtual ~dll_param_number();
00268   friend class model_parameters;
00269 };
00270 
00275 class named_dvar_matrix : public dvar_matrix, public model_name_tag
00276 {
00277 protected:
00278   named_dvar_matrix(void) : dvar_matrix(), model_name_tag() {;}
00279   void allocate(int rmin,int rmax,int cmin,int cmax,const char * s);
00280   void allocate(int rmin,int rmax,const char * s);
00281   void allocate(const char * s);
00282   //void allocate(int rmin,int rmax,int,const ivector&, const char * s);
00283   void allocate(
00284     int rmin, int rmax,
00285     const index_type& cmin, const index_type& cmax,
00286     const char* s);
00287   void allocate(
00288     int rmin, int rmax,
00289     const data_int& cmin, const data_int& cmax,
00290     const char* s)
00291   {
00292     allocate(rmin, rmax, index_type(cmin), index_type(cmax), s);
00293   }
00294 public:
00295   named_dvar_matrix& operator=(const double m);
00296   named_dvar_matrix& operator=(const dmatrix& m);
00297   named_dvar_matrix& operator=(const dvar_matrix& m);
00298   named_dvar_matrix& operator=(const dvariable& m);
00299 
00300   friend class model_parameters;
00301 };
00302 
00307 class dll_param_matrix : public named_dvar_matrix
00308 {
00309   double * pd;
00310 public:
00311   //named_dvar_matrix(void) : dvar_matrix(), model_name_tag() {;}
00312   void allocate(double * ,int rmin,int rmax,int cmin,int cmax,const char * s);
00313   void allocate(double * ,int rmin,int rmax,const index_type&
00314      ,const index_type& ,const char * s);
00315   dll_param_matrix& operator=(const dvariable& m);
00316   dll_param_matrix & operator = (const double m);
00317   dll_param_matrix& operator=(const dmatrix& m);
00318   dll_param_matrix& operator=(const dvar_matrix& m);
00319   virtual ~dll_param_matrix();
00320 };
00321 
00326 class named_dvar3_array : public dvar3_array, public model_name_tag
00327 {
00328 protected:
00329   named_dvar3_array(void) : dvar3_array(), model_name_tag() {;}
00330  void allocate(const ad_integer& sl,const ad_integer& sh,
00331   const index_type& nrl,const index_type& nrh,const index_type& ncl,
00332     const index_type& nch,const char * s="UNNAMED");
00333   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
00334     const char * s="UNNAMED");
00335 
00336   void allocate(int hsl,int hsu,int rmin,int rmax,const char * s="UNNAMED");
00337   void allocate(int hsl,int hsu,const index_type& rmin,const index_type& rmax,
00338     const char * s="UNNAMED");
00339   void allocate(int hsl,int hsu,const char * s="UNNAMED");
00340   void allocate(const char * s="UNNAMED");
00341 
00342 
00343   void allocate(int hsl, int hsu, const ivector& rmin, int rmax,
00344                 int cmin, int cmax, const char *s = "UNNAMED");
00345   void allocate(int hsl, int hsu, int rmin, const ivector& rmax,
00346     int cmin,int cmax,const char * s="UNNAMED");
00347   void allocate(int hsl, int hsu, const ivector& rmin, const ivector& rmax,
00348     int cmin,int cmax,const char * s="UNNAMED");
00349   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
00350     int cmax,const char * s="UNNAMED");
00351   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
00352                 const ivector& cmax, const char *s = "UNNAMED");
00353   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,
00354                 const ivector& cmax, const char *s = "UNNAMED");
00355   named_dvar3_array& operator=(const dvar3_array& m);
00356   named_dvar3_array& operator=(const d3_array& m);
00357   friend class model_parameters;
00358 };
00359 
00364 class named_dvar4_array : public dvar4_array, public model_name_tag
00365 {
00366 protected:
00367   named_dvar4_array(void) : dvar4_array(), model_name_tag() {;}
00368   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,int rmax,
00369   int cmin,int cmax,const char * s);
00370   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00371     const index_type&,const index_type&,const index_type&,const index_type&,
00372     const char *s);
00373 
00374   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00375     const index_type&,const index_type&,const char *s);
00376 
00377   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00378     const char *s);
00379   void allocate(ad_integer,ad_integer,const char *s);
00380   void allocate(const char *s);
00381 
00382   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,int rmax,
00383     const char * s);
00384   void allocate(int hhsl,int hhsu,int hsl,int hsu,const char * s);
00385   void allocate(int hhsl,int hhsu,const char * s);
00386 
00387   named_dvar4_array& operator=(const dvar4_array& m);
00388   named_dvar4_array& operator=(const d4_array& m);
00389   friend class model_parameters;
00390 };
00391 
00396 class named_dvar5_array : public dvar5_array, public model_name_tag
00397 {
00398 protected:
00399   named_dvar5_array(void) : dvar5_array(), model_name_tag() {;}
00400   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,int rmax,
00401   int cmin,int cmax,
00402   int l5,int u5,
00403   const char * s);
00404 
00405   void allocate(const ad_integer& hhsl,const ad_integer& hhsu,
00406     const index_type& hsl,const index_type& hsu,
00407     const index_type& sl,const index_type& sh,
00408     const index_type& nrl,const index_type& nrh,
00409     const index_type& ncl,const index_type& nch,
00410     const char * s);
00411 
00412   named_dvar5_array& operator=(const dvar5_array& m);
00413   named_dvar5_array& operator=(const d5_array& m);
00414   friend class model_parameters;
00415 };
00416 
00421 class named_dvar6_array : public dvar6_array, public model_name_tag
00422 {
00423 protected:
00424   named_dvar6_array(void) : dvar6_array(), model_name_tag() {;}
00425   void allocate(int hhsl,int hhsu,
00426     int hsl,int hsu,
00427     int rmin,int rmax,
00428     int cmin,int cmax,
00429     int l5,int u5,
00430     int l6,int u6,
00431     const char * s);
00432 
00433   void allocate(const ad_integer& hhsl,const ad_integer& hhsu,
00434     const index_type& hsl,const index_type& hsu,
00435     const index_type& sl,const index_type& sh,
00436     const index_type& nrl,const index_type& nrh,
00437     const index_type& ncl,const index_type& nch,
00438     const index_type& l5,const index_type& u5,
00439     const char * s);
00440 
00441   named_dvar6_array& operator=(const dvar6_array& m);
00442   named_dvar6_array& operator=(const d6_array& m);
00443   friend class model_parameters;
00444 };
00445 
00450 class named_dvar7_array : public dvar7_array, public model_name_tag
00451 {
00452 protected:
00453   named_dvar7_array(void) : dvar7_array(), model_name_tag() {;}
00454   void allocate(
00455     int hhsl,int hhsu,
00456     int hsl,int hsu,
00457     int rmin,int rmax,
00458     int cmin,int cmax,
00459     int l5,int u5,
00460     int l6,int u6,
00461     int l7,int u7,
00462     const char * s);
00463 
00464   void allocate(const ad_integer& hhsl,const ad_integer& hhsu,
00465     const index_type& hsl,const index_type& hsu,
00466     const index_type& sl,const index_type& sh,
00467     const index_type& nrl,const index_type& nrh,
00468     const index_type& ncl,const index_type& nch,
00469     const index_type& l5,const index_type& u5,
00470     const index_type& l6,const index_type& u6,
00471     const char * s);
00472 
00473   named_dvar7_array& operator=(const dvar7_array& m);
00474   named_dvar7_array& operator=(const d7_array& m);
00475   friend class model_parameters;
00476 };
00477 
00482 class named_dvector : public dvector, public model_name_tag
00483 {
00484 protected:
00485   named_dvector(void) : dvector(), model_name_tag() {;}
00486   void allocate(int mmin,int mmax,const char * s);
00487   void allocate(const char * s);
00488   void allocate(int mmin, const ivector& mmax, const char *s);
00489   named_dvector& operator=(const dvector& m);
00490   named_dvector& operator=(const double m);
00491   friend class model_data;
00492 };
00493 
00498 class named_ivector : public ivector, public model_name_tag
00499 {
00500 protected:
00501   named_ivector(void) : ivector(), model_name_tag() {;}
00502   void allocate(int mmin,int mmax,const char * s);
00503 };
00504 
00509 class named_dmatrix : public dmatrix, public model_name_tag
00510 {
00511 protected:
00512   named_dmatrix(void) : dmatrix(), model_name_tag() {;}
00513   void allocate(int rmin,int rmax,const char * s);
00514   void allocate(const char * s);
00515   void allocate(int rmin,int rmax,int cmin,int cmax,const char * s);
00516   void allocate(int rmin,int rmax,const ivector& cmin,int cmax,const char * s);
00517   void allocate(int rmin,int rmax,const ivector& cmin,const ivector& cmax,
00518     const char * s);
00519   void allocate(int rmin,int rmax,int cmin,const ivector& cmax,const char * s);
00520   named_dmatrix& operator=(const dmatrix& m);
00521   named_dmatrix& operator=(const double m);
00522 };
00523 
00528 class named_imatrix : public imatrix, public model_name_tag
00529 {
00530 protected:
00531   named_imatrix(void) : imatrix(), model_name_tag() {;}
00532   void allocate(int rmin,int rmax,int cmin,int cmax,const char * s);
00533   void allocate(int rmin, int rmax, const index_type& cmin,
00534                 const index_type& cmax, const char *s);
00535   named_imatrix& operator=(const imatrix& m);
00536   named_imatrix& operator=(const int& m);
00537 };
00538 
00543 class named_d3_array : public d3_array, public model_name_tag
00544 {
00545 protected:
00546   named_d3_array(void) : d3_array(), model_name_tag() {;}
00547   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
00548     const char * s);
00549   void allocate(int hsl, int hsu, const index_type& rmin,
00550                 const index_type& rmax, const index_type& cmin,
00551                 const index_type& cmax, const char *s);
00552   void allocate(int hsl, int hsu, const ivector& rmin, int rmax,
00553     int cmin,int cmax,const char * s);
00554   void allocate(int hsl, int hsu, int rmin, const ivector& rmax,
00555     int cmin,int cmax,const char * s);
00556   void allocate(int hsl, int hsu, const ivector& rmin, const ivector& rmax,
00557     int cmin,int cmax,const char * s);
00558   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
00559     int cmax, const char * s);
00560   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,
00561                 const ivector& cmax, const char *s);
00562   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
00563                 const ivector& cmax, const char *s);
00564   named_d3_array& operator=(const d3_array& m);
00565 
00566   void allocate(int hsl,int hsu,const index_type& rmin,const index_type& rmax,
00567     const char * s);
00568   void allocate(int hsl,int hsu,int rmin,int rmax,const char * s);
00569   void allocate(int hsl,int hsu,const char * s);
00570   void allocate(const char * s);
00571 };
00572 
00577 class named_i3_array : public i3_array, public model_name_tag
00578 {
00579 protected:
00580   named_i3_array(void) : i3_array(), model_name_tag() {;}
00581   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
00582     const char * s);
00583   void allocate(int hsl, int hsu, const index_type& rmin,
00584                 const index_type& rmax, const index_type& cmin,
00585                 const index_type& cmax, const char *s);
00586   named_i3_array& operator=(const i3_array& m);
00587 };
00588 
00593 class named_d4_array : public d4_array, public model_name_tag
00594 {
00595 protected:
00596   named_d4_array(void) : d4_array(), model_name_tag() {;}
00597   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00598     const index_type&,const index_type&,const index_type&,const index_type&,
00599     const char * s);
00600   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,
00601     int rmax,int cmin,int cmax,const char * s);
00602 
00603   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,
00604     int rmax,const char * s);
00605   void allocate(int hhsl,int hhsu,int hsl,int hsu,const char * s);
00606   void allocate(int hhsl,int hhsu,const char * s);
00607   void allocate(const char * s);
00608 
00609   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00610     const index_type&,const index_type&,const char * s);
00611   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00612     const char * s);
00613   void allocate(ad_integer,ad_integer,const char * s);
00614 
00615 
00616   named_d4_array& operator=(const d4_array& m);
00617 };
00618 class named_i4_array : public i4_array, public model_name_tag
00619 {
00620 protected:
00621   named_i4_array(void) : i4_array(), model_name_tag() {;}
00622   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00623     const index_type&,const index_type&,const index_type&,const index_type&,
00624     const char * s);
00625   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,
00626     int rmax,int cmin,int cmax,const char * s);
00627 
00628   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,
00629     int rmax,const char * s);
00630   void allocate(int hhsl,int hhsu,int hsl,int hsu,const char * s);
00631   void allocate(int hhsl,int hhsu,const char * s);
00632   void allocate(const char * s);
00633 
00634   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00635     const index_type&,const index_type&,const char * s);
00636   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00637     const char * s);
00638   void allocate(ad_integer,ad_integer,const char * s);
00639 
00640 
00641   named_i4_array& operator=(const i4_array& m);
00642 };
00643 
00648 class named_d5_array : public d5_array, public model_name_tag
00649 {
00650 protected:
00651   named_d5_array(void) : d5_array(), model_name_tag() {;}
00652   void allocate(int l5,int u5,int hhsl,int hhsu,int hsl,int hsu,int rmin,
00653     int rmax,int cmin,int cmax,const char * s);
00654   void allocate(const ad_integer& hhsl,const ad_integer& hhsu,
00655     const index_type& hsl,const index_type& hsu, const index_type& sl,
00656     const index_type& sh,const index_type& nrl,const index_type& nrh,
00657     const index_type& ncl,const index_type& nch,const char * s);
00658 
00659   named_d5_array& operator=(const d5_array& m);
00660 };
00661 
00666 class named_d6_array : public d6_array, public model_name_tag
00667 {
00668 protected:
00669   named_d6_array(void) : d6_array(), model_name_tag() {;}
00670   void allocate(int l6,int u6,int l5,int u5,int hhsl,int hhsu,int hsl,
00671     int hsu,int rmin,int rmax,int cmin,int cmax,const char * s);
00672   void allocate(const ad_integer& l6,const ad_integer& u6,
00673     const index_type& l5,const index_type& u5,
00674     const index_type& hhsl,const index_type& hhsu,
00675     const index_type& hsl,const index_type& hsu,
00676     const index_type& sl,const index_type& sh,
00677     const index_type& nrl,const index_type& nrh,
00678     const char * s);
00679 
00680   named_d6_array& operator=(const d6_array& m);
00681 };
00682 
00687 class named_d7_array : public d7_array, public model_name_tag
00688 {
00689 protected:
00690   named_d7_array(void) : d7_array(), model_name_tag() {;}
00691   void allocate(int l7,int u7,int l6,int u6,int l5,int u5,int hhsl,
00692     int hhsu,int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
00693     const char * s);
00694   void allocate(const ad_integer& l7,const ad_integer& u7,
00695     const index_type& l6,const index_type& u6,
00696     const index_type& l5,const index_type& u5,
00697     const index_type& hhsl,const index_type& hhsu,
00698     const index_type& hsl,const index_type& hsu,
00699     const index_type& sl,const index_type& sh,
00700     const index_type& nrl,const index_type& nrh,
00701     const char * s);
00702 
00703   named_d7_array& operator=(const d7_array& m);
00704 };
00705 
00706 
00707 class function_minimizer;
00708 
00709 #if defined(USE_ADPVM)
00710 
00714 class pvm_params
00715 {
00716   static pvm_params * varsptr[]; // this should be a resizeable array
00717   static int num_pvm_params;
00718   static const int maxnum_pvm_params;
00719   void add_to_list(void);
00720   virtual void send_to_slaves(void)=0;
00721   virtual void get_from_master(void)=0;
00722 public:
00723   static void pvm_params::send_all_to_slaves(void);
00724   static void pvm_params::get_all_from_master(void);
00725   void allocate(const char *);
00726   void allocate(void);
00727 };
00728 
00733 class pvm_number : public pvm_params
00734 {
00735 public:
00736   virtual void send_to_slaves(void);
00737   virtual void get_from_master(void);
00738   dvector v;
00739   double d;
00740   operator double();
00741   void assign(const dvector&);
00742   void assign(double);
00743 };
00744 
00749 class pvm_int : public pvm_params
00750 {
00751 public:
00752   virtual void send_to_slaves(void);
00753   virtual void get_from_master(void);
00754   ivector v;
00755   int d;
00756   operator int();
00757   void assign(const ivector&);
00758   void assign(int);
00759 };
00760 #endif // #if defined(USE_ADPVM)
00761 
00762 class initial_params;
00763 typedef initial_params* pinitial_params;
00764 typedef void* ptovoid;
00768 class adlist_ptr
00769 {
00770   ptovoid * ptr;
00771   int current_size;
00772   int current;
00773   void resize(void);
00774   void add_to_list(void* p);
00775 public:
00776   adlist_ptr(int init_size);
00777   ~adlist_ptr();
00778 
00779   void initialize();
00780 
00781   pinitial_params& operator[](int i);
00782 
00783   friend class initial_params;
00784 };
00785 
00786 #if defined(USE_SHARE_FLAGS)
00787 
00792   class shareinfo
00793   {
00794     index_type * sflags;
00795     index_type * original_sflags;
00796     index_type * aflags;
00797     index_type * invflags;
00798     i3_array * bmap;
00799     int dimension;
00800     int maxshare;
00801     int current_phase;
00802   public:
00803     void get_inv_matrix_shared( int _cf);
00804     void get_inv_vector_shared( int _cf);
00805     int & get_maxshare(void) { return maxshare; }
00806     i3_array & get_bmap(void);
00807     int & get_dimension(void){ return dimension;}
00808     index_type * get_original_shareflags(void);
00809     index_type * get_shareflags(void);
00810     int& get_current_phase(void);
00811     index_type * get_activeflags(void);
00812     index_type * get_invflags(void);
00813     void set_shareflags(const index_type& sf);
00814     void set_original_shareflags(const index_type& sf);
00815     void reset_shareflags(const index_type& sf);
00816     void set_activeflags(const index_type& af);
00817     void set_bmap(const i3_array& af);
00818     void reset_bmap(const i3_array& af);
00819     void set_invflags(const index_type& af);
00820     shareinfo(const index_type& sf,const index_type& af);
00821     ~shareinfo();
00822   };
00823 #endif
00824 
00829 class initial_params
00830 {
00831 protected:
00832 #if defined(USE_SHARE_FLAGS)
00833   shareinfo * share_flags;
00834 #endif
00835   virtual ~initial_params();
00836   int active_flag;
00837   int initial_value_flag;
00838   double initial_value;
00839   double scalefactor;
00840 public:
00841 #if defined(USE_SHARE_FLAGS)
00842   virtual void setshare(const index_type& sf,const index_type& af);
00843   virtual void shared_set_value_inv(const dvector&,const int&);
00844   virtual void shared_set_value(const dvar_vector&,const int&,
00845     const dvariable& pen);
00846   virtual int shared_size_count(void); // get the number of active parameters
00847 #endif
00848   double get_scalefactor();
00849   void set_scalefactor(const double);
00850 #if !defined(BIG_INIT_PARAMS)
00851   static initial_params  varsptr[]; // this should be a resizeable array
00852 #else
00853   static adlist_ptr varsptr; // this should be a resizeable array
00854 #endif
00855   static int num_initial_params;
00856   static const int max_num_initial_params;
00857   static int straight_through_flag;
00858   static int num_active_initial_params;
00859   static int max_number_phases;
00860   static int current_phase;
00861   static int restart_phase;
00862   static int sd_phase;
00863   static int mc_phase;
00864   static int mceval_phase;
00865   int phase_start;
00866   int phase_save;
00867   int phase_stop;
00868   virtual void set_random_effects_active();
00869   void restore_phase_start(void);
00870   virtual void set_random_effects_inactive();
00871   virtual void set_only_random_effects_active();
00872   virtual void set_only_random_effects_inactive();
00873   virtual void set_value(const dvar_vector&, const int&,
00874     const dvariable& pen) = 0;
00875   virtual void dev_correction(const dmatrix&, const int&) = 0;
00876   void set_initial_value(double x);
00877   double get_initial_value(void);
00878   void set_phase_start(int x);
00879   int get_phase_start(void);
00880   static void set_all_simulation_bounds(const dmatrix& symbds);
00881   static void set_all_simulation_bounds(const dmatrix& symbds, const dvector&);
00882   static void get_jacobian_value(const dvector& y, const dvector& jac);
00883   static int correct_for_dev_objects(const dmatrix &H);
00884   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii) = 0;
00885   //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
00886   //  const int& ii) = 0;
00887   virtual void set_value_inv(const dvector&, const int&) = 0;
00888   virtual void add_value(const dvector&, const int&) = 0;
00889   virtual void add_value(const dvector&, const dvector&, const int&,
00890     const double&, const dvector&) = 0;
00891   virtual void get_jacobian(const dvector&, const dvector&, const int&) = 0;
00892   //virtual void check_tightness(const dvector&, const int&) = 0;
00893   virtual void copy_value_to_vector(const dvector&, const int&) = 0;
00894   virtual void restore_value_from_vector(const dvector&, const int&) = 0;
00895   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii) = 0;
00896   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
00897     const int& ii)=0;
00898   virtual void mc_scale(const dvector& d, const dvector& x, const int& ii) = 0;
00899   virtual void curv_scale(const dvector& d, const dvector& x,
00900     const int& ii) = 0;
00901   virtual void hess_scale(const dvector& d, const dvector& x,
00902     const int& ii) = 0;
00903   virtual int size_count(void)=0; // get the number of active parameters
00904 
00905   // save the objects value in to a text file
00906   virtual void save_value(ofstream& ofs) = 0;
00907 
00908   // save the objects value in a binary file
00909   virtual void bsave_value(uostream& uos) = 0;
00910 
00911   virtual void save_value(const ofstream& ofs, int prec) = 0;
00912   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
00913     int& offset)=0;
00914   virtual const char* label() = 0;
00915   void allocate(int _phase_start);
00916   void set_active_flag(void);
00917   void set_inactive_flag(void);
00918   friend int active(const initial_params& ip);
00919   static adstring get_reportfile_name(void);
00920   initial_params(void);
00921   static void xinit(const dvector& x); // get the number of active parameters
00922   // get the number of active parameters
00923   static void xinit_all(const dvector& x);
00924   static void save_all(const ofstream& _ofs,int prec,const dvector&g);
00925   // get the number of active parameters
00926   static void set_active_random_effects(void);
00927   // get the number of active parameters
00928   static void set_active_only_random_effects(void);
00929   // get the number of active parameters
00930   static void set_inactive_only_random_effects(void);
00931   // get the number of active parameters
00932   static void set_inactive_random_effects(void);
00933   // get the number of active parameters
00934   static void restore_start_phase(void);
00935   static void xinit1(const dvector& x, const dvector& g);
00936   //save all initial parameter values in a vector
00937   static void copy_all_values(const dvector& x, const int& ii);
00938   //get ivalues for all active parameters from a vector
00939   static void restore_all_values(const dvector& x, const int& ii);
00940   // get the number of active parameters
00941   static dvariable reset(const dvar_vector& x);
00942   // get the number of active parameters
00943   static dvariable reset(const dvector& x);
00944   static dvariable reset1(const dvar_vector& x, const dvector& g);
00945   // get the number of active parameters
00946   static dvariable reset(const dvar_vector& x, const dvector& pen);
00947   // get the number of active parameters
00948   static dvariable reset_all(const dvar_vector& x,const dvector& pen);
00949   static int nvarcalc();
00950   static int nvarcalc_all(void);
00951   static int num_active_calc(void);
00952   static int stddev_scale(const dvector& d, const dvector& x);
00953   static int stddev_vscale(const dvar_vector& d,const dvar_vector& x);
00954   static int montecarlo_scale(const dvector& d, const dvector& x);
00955   static int stddev_curvscale(const dvector& d, const dvector& x);
00956   static void read(void);
00957   static void save(void);
00958   static void save(const ofstream& ofs, int prec);
00959   static void restore(const ifstream& ifs);
00960   static void add_random_vector(const dvector& x);
00961   static void add_random_vector(const dvector& y, const dvector& x,
00962     const double& ll, const dvector& diag);
00963   virtual void restore_value(const ifstream& ifs) = 0;
00964   virtual void add_to_list(void);
00965 #if defined(USE_ADPVM)
00966   virtual void pvm_pack(void)=0;
00967   virtual void pvm_unpack(void)=0;
00968 #endif
00969 
00970   friend class function_minimizer;
00971 };
00972 
00973 void pvm_pack(const dvar_vector&);
00974 void pvm_unpack(const dvar_vector&);
00975 void pvm_pack(const prevariable&);
00976 void pvm_unpack(const prevariable&);
00977 void pvm_pack(const dvar_matrix&);
00978 void pvm_unpack(const dvar_matrix&);
00979 
00984 class param_init_vector: public named_dvar_vector, public initial_params
00985 {
00986 public:
00987 #if defined(USE_SHARE_FLAGS)
00988   virtual void setshare(const index_type& sf,const index_type& af);
00989   virtual void shared_set_value_inv(const dvector&,const int&);
00990   virtual void shared_set_value(const dvar_vector&,const int&,
00991     const dvariable& pen);
00992   virtual int shared_size_count(void); // get the number of active parameters
00993 #endif
00994   virtual void dev_correction(const dmatrix& H, const int& ii);
00995   virtual const char * label(void);
00996   param_init_vector();
00997 #if defined(USE_ADPVM)
00998   void pvm_pack(void){::pvm_pack(*this);}
00999   void pvm_unpack(void){::pvm_unpack(*this);}
01000 #endif
01001 private:
01002   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01003   //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01004   //  const int& ii);
01005   virtual void add_value(const dvector&, const dvector&, const int&,
01006     const double&, const dvector&);
01007   virtual void add_value(const dvector&, const int&);
01008   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01009   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01010     const int& ii);
01011 
01012   virtual void set_value(const dvar_vector& x, const int& ii,
01013     const dvariable& pen);
01014   virtual void set_value_inv(const dvector& x, const int& ii);
01015   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01016   virtual void restore_value_from_vector(const dvector&, const int&);
01017   virtual int size_count(void);
01018   virtual void sd_scale(const dvector& d, const dvector& x, const int&);
01019   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01020   virtual void curv_scale(const dvector& d, const dvector& x,const int&);
01021   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01022   virtual void save_value(ofstream& ofs);
01023   virtual void bsave_value(uostream& uos);
01024   virtual void save_value(const ofstream& ofs, int prec);
01025   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
01026     int& offset);
01027   virtual void restore_value(const ifstream& ifs);
01028   void report_value(void);
01029   //virtual void read_value(void);
01030   void allocate(int imin,int imax,int phasestart=1,const char * s="UNNAMED");
01031   void allocate(const ad_integer& imin,const ad_integer& imax,
01032     const ad_integer& phasestart=1,const char * s="UNNAMED");
01033   void allocate(int imin,int imax,const char * s="UNNAMED");
01034 #if defined(USE_SHARE_FLAGS)
01035   void allocate(int imin,int imax,const ivector& ishare,
01036     const char * s="UNNAMED");
01037 #endif
01038   friend class model_parameters;
01039   friend class param_init_vector_vector;
01040 public:
01041   param_init_vector& operator = (const dvector&);
01042   param_init_vector& operator = (const dvar_vector&);
01043   param_init_vector& operator = (const prevariable&);
01044   param_init_vector& operator = (const double&);
01045 };
01046 
01051 class dll_param_init_vector: public param_init_vector
01052 {
01053   double * pd;
01054 public:
01055   dll_param_init_vector& operator = (const dvector&);
01056   dll_param_init_vector& operator = (const dvar_vector&);
01057   dll_param_init_vector& operator = (const prevariable&);
01058   dll_param_init_vector& operator = (const double&);
01059   void allocate(double * _pd,int imin,int imax,
01060     int phasestart=1,const char * s="UNNAMED");
01061   void allocate(double * _pd,int imin,int imax,
01062     const char * s="UNNAMED");
01063 
01064   virtual ~dll_param_init_vector();
01065   friend class model_parameters;
01066 };
01067 
01072 class param_init_bounded_vector: public named_dvar_vector, public initial_params
01073 {
01074   virtual void* parent_this(void){return this;}
01075 public:
01076   double get_minb(void);
01077   void set_minb(double b);
01078   double get_maxb(void);
01079   void set_maxb(double b);
01080 protected:
01081   double minb;
01082   double maxb;
01083   param_init_bounded_vector();
01084 private:
01085   virtual void dev_correction(const dmatrix&, const int&);
01086   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01087   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01088     const int& ii);
01089   //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01090   //  const int& ii);
01091   virtual void add_value(const dvector&, const dvector&, const int&,
01092     const double&, const dvector&);
01093   virtual void add_value(const dvector&, const int&);
01094   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01095   virtual void set_value(const dvar_vector& x, const int& ii,
01096     const dvariable& pen);
01097   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01098   virtual void restore_value_from_vector(const dvector&, const int&);
01099   virtual void set_value_inv(const dvector& x, const int& ii);
01100   virtual int size_count(void);
01101   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01102   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01103   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01104   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01105   void allocate(int imin,int imax,double _minb,double _maxb,
01106     int phasestart=1, const char * name="UNNAMED");
01107   void allocate(int imin,int imax,double _minb,double _maxb,
01108     const char * name="UNNAMED");
01109  //void param_init_bounded_vector::allocate(const ad_integer& imin,
01110    //const ad_integer& imax,const ad_double& _minb,const ad_double& _maxb,
01111    //const ad_integer& phase_start,const char * s);
01112   friend class model_parameters;
01113   friend class param_init_bounded_vector_vector;
01114   virtual const char * label(void);
01115   virtual void save_value(const ofstream& ofs, int prec);
01116   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
01117     int& offset);
01118   virtual void restore_value(const ifstream& ifs);
01119   virtual void save_value(ofstream& ofs);
01120   virtual void bsave_value(uostream& uos);
01121   void report_value(void);
01122   //virtual void read_value(void);
01123 public:
01124   param_init_bounded_vector& operator = (const dvector&);
01125   param_init_bounded_vector& operator = (const dvar_vector&);
01126   param_init_bounded_vector& operator = (const prevariable&);
01127   param_init_bounded_vector& operator = (const double&);
01128 #if defined(USE_ADPVM)
01129   void pvm_pack(void){::pvm_pack(*this);}
01130   void pvm_unpack(void){::pvm_unpack(*this);}
01131 #endif
01132 };
01133 
01138 class dll_param_init_bounded_vector: public param_init_bounded_vector
01139 {
01140   double * pd;
01141 public:
01142   void allocate(double * _pd,int imin,int imax,double _minb,double _maxb,
01143     int phasestart=1, const char * name="UNNAMED");
01144   void allocate(double * _pd,int imin,int imax,double _minb,double _maxb,
01145     const char * name="UNNAMED");
01146   ~dll_param_init_bounded_vector();
01147 };
01148 
01153 class param_init_bounded_dev_vector: public param_init_bounded_vector
01154 {
01155   virtual void set_value(const dvar_vector& x, const int& ii,
01156     const dvariable& pen);
01157   virtual void dev_correction(const dmatrix& H, const int& ii);
01158 public:
01159   param_init_bounded_dev_vector& operator = (const dvar_vector& m);
01160   param_init_bounded_dev_vector& operator = (const dvector& m);
01161   param_init_bounded_dev_vector& operator = (const prevariable& m);
01162   param_init_bounded_dev_vector& operator = (const double& m);
01163 };
01164 
01169 class param_init_number: public named_dvariable, public initial_params
01170 {
01171   virtual void dev_correction(const dmatrix&, const int&);
01172   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01173 
01174 #if defined(USE_ADPVM)
01175   void pvm_pack(void){::pvm_pack(*this);}
01176   void pvm_unpack(void){::pvm_unpack(*this);}
01177 #endif
01178 
01179 //  virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01180 //    const int& ii);
01181   virtual void add_value(const dvector&, const dvector&, const int&,
01182     const double&, const dvector&);
01183   virtual void add_value(const dvector&, const int&);
01184   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01185   virtual void set_value(const dvar_vector& x, const int& ii,
01186     const dvariable& pen);
01187   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01188   virtual void restore_value_from_vector(const dvector&, const int&);
01189   virtual void set_value_inv(const dvector& x, const int& ii);
01190   virtual int size_count(void);
01191   virtual void save_value(const ofstream& ofs, int prec);
01192   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
01193     int& offset);
01194   virtual void restore_value(const ifstream& ifs);
01195   virtual void save_value(ofstream& ofs);
01196   virtual void bsave_value(uostream& uos);
01197   void report_value(void);
01198   virtual const char * label(void);
01199   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01200   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01201   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01202   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01203   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01204     const int& ii);
01205   //virtual void read_value(void);
01206 protected:
01207   void allocate(int phase_start=1,const char *s="UNNAMED");
01208   void allocate(const char *s="UNNAMED");
01209   void allocate(init_xml_doc&, const char *s="UNNAMED");
01210   friend class model_parameters;
01211   friend class param_init_number_vector;
01212   param_init_number();
01213   param_init_number& operator=(const double m);
01214   param_init_number& operator=(const prevariable& m);
01215 };
01216 
01221 class dll_param_init_number: public param_init_number
01222 {
01223   double * pd;
01224 public:
01225   void allocate(double * pd,int phase_start=1,const char *s="UNNAMED");
01226   void allocate(double *pd,const char *s="UNNAMED");
01227   virtual ~dll_param_init_number();
01228   dll_param_init_number& operator=(const double m);
01229   dll_param_init_number& operator=(const prevariable& m);
01230 };
01231 
01232 //forward declaration
01233 class data_vector;
01240 class param_init_bounded_number: public param_init_number
01241 {
01242 public:
01243   double get_minb(void);
01244   void set_minb(double b);
01245   double get_maxb(void);
01246   void set_maxb(double b);
01247 protected:
01248   double minb;
01249   double maxb;
01250   void allocate(double _minb,double _maxb,int phase_start=1,
01251     const char * s="UNNAMED");
01252   void allocate(double _minb,double _maxb,const char * s="UNNAMED");
01253   // Added by Steve Martell for using input data for allocation.
01254   void allocate(const data_vector& v,const char * s="UNNAMED");
01255   void allocate(init_xml_doc&, const char * s="UNNAMED");
01256 
01257 public:
01258 #if defined(USE_ADPVM)
01259   void pvm_pack(void){::pvm_pack(*this);}
01260   void pvm_unpack(void){::pvm_unpack(*this);}
01261 #endif
01262   param_init_bounded_number();
01263 private:
01264   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01265 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01266 //  const int& ii);
01267   virtual void add_value(const dvector&, const dvector&, const int&,
01268     const double&, const dvector&);
01269   virtual void add_value(const dvector&, const int&);
01270   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01271   virtual void set_value(const dvar_vector& x, const int& ii,
01272     const dvariable& pen);
01273   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01274   virtual void restore_value_from_vector(const dvector&, const int&);
01275   virtual void set_value_inv(const dvector& x, const int& ii);
01276   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01277   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01278   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01279   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01280   virtual const char * label(void);
01281   void report_value(void);
01282   param_init_bounded_number& operator=(const double m);
01283   param_init_bounded_number& operator=(const prevariable& m);
01284   friend class model_parameters;
01285   friend class param_init_bounded_number_vector;
01286   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01287     const int& ii);
01288 };
01289 
01294 class dll_param_init_bounded_number: public param_init_bounded_number
01295 {
01296   double * pd;
01297 public:
01298   void allocate(double * _pd,double _minb,double _maxb,int phase_start=1,
01299     const char * s="UNNAMED");
01300   void allocate(double* _pd,double _minb, double _maxb,const char* s="UNNAMED");
01301 public:
01302   virtual ~dll_param_init_bounded_number();
01303   void report_value(void);
01304 };
01305 
01310 class param_init_matrix: public named_dvar_matrix,public initial_params
01311 {
01312 #if defined(USE_SHARE_FLAGS)
01313   virtual void shared_set_value_inv(const dvector& x,const int& ii);
01314   virtual int shared_size_count(void); // get the number of active parameters
01315   virtual void shared_set_value(const dvar_vector&,const int&,
01316     const dvariable& pen);
01317 #endif
01318   virtual void dev_correction(const dmatrix&, const int&);
01319   virtual void set_simulation_bounds(const dmatrix& symbds,const int& ii);
01320 #if defined(USE_ADPVM)
01321   void pvm_pack(void){::pvm_pack(*this);}
01322   void pvm_unpack(void){::pvm_unpack(*this);}
01323 #endif
01324 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01325 //  const int& ii);
01326   virtual void add_value(const dvector&, const dvector&, const int&,
01327     const double&, const dvector&);
01328   virtual void add_value(const dvector&, const int&);
01329   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01330 public:
01331   virtual void set_value(const dvar_vector& x, const int& ii,
01332     const dvariable& pen);
01333   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01334   virtual void restore_value_from_vector(const dvector&, const int&);
01335   virtual void set_value_inv(const dvector& x, const int& ii);
01336   virtual int size_count(void);
01337   virtual void save_value(ofstream& ofs);
01338   virtual void bsave_value(uostream& uos);
01339   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
01340     int& offset);
01341   virtual void save_value(const ofstream& ofs, int prec);
01342   virtual void restore_value(const ifstream& ifs);
01343   void report_value(void);
01344   //virtual void read_value(void);
01345   virtual const char * label(void);
01346   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01347   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01348   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01349   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01350   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01351     const int& ii);
01352 
01353 public:
01354   void allocate(int rmin,int rmax,int cmin,int cmax,
01355     int phase_start=1,const char * = "UNNAMED");
01356 #if defined(USE_SHARE_FLAGS)
01357   //void allocate(int rmin,int rmax,int cmin,int cmax, const imatrix& jshare,
01358   // const char * = "UNNAMED");
01359   virtual void setshare(const index_type& sf,const index_type& af);
01360 #endif
01361   void allocate(int rmin,int rmax,int cmin,int cmax,
01362     const char * = "UNNAMED");
01363  void allocate(const ad_integer& imin,const ad_integer&imax,
01364    const index_type& imin2,const index_type& imax2,
01365    const ad_integer& phase_start, const char * s);
01366   void allocate(const ad_integer& rmin, const ad_integer& rmax,
01367                 const index_type& cmin, const index_type& cmax,
01368                 const char * = "UNNAMED");
01369   void allocate(const ad_integer& rmin, const ad_integer& rmax,
01370                 const index_type& cmin, const index_type& cmax,
01371                 int phase_start = 1, const char * = "UNNAMED");
01372   param_init_matrix(void);
01373   param_init_matrix& operator = (const dmatrix& m);
01374   param_init_matrix& operator = (const dvar_matrix& m);
01375   param_init_matrix& operator = (const dvariable& m);
01376   param_init_matrix& operator = (const double& m);
01377 };
01378 
01383 class dll_param_init_matrix: public param_init_matrix
01384 {
01385   double * d;
01386 public:
01387   void allocate(double* _d,int rmin,int rmax,int cmin,int cmax,
01388     int phase_start=1,const char * = "UNNAMED");
01389   void allocate(double * _d,int rmin,int rmax,int cmin,int cmax,
01390     const char * = "UNNAMED");
01391   virtual ~dll_param_init_matrix();
01392   dll_param_init_matrix(){d=NULL;}
01393   dll_param_init_matrix& operator = (const dmatrix& m);
01394   dll_param_init_matrix& operator = (const dvar_matrix& m);
01395   dll_param_init_matrix& operator = (const dvariable& m);
01396   dll_param_init_matrix& operator = (const double& m);
01397 };
01398 
01403 class param_init_bounded_matrix: public param_init_matrix
01404 {
01405 public:
01406   double get_minb(void);
01407   void set_minb(double b);
01408   double get_maxb(void);
01409   void set_maxb(double b);
01410 protected:
01411   double minb;
01412   double maxb;
01413   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01414 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01415 //  const int& ii);
01416   virtual void add_value(const dvector&, const dvector&, const int&,
01417     const double&, const dvector&);
01418   virtual void add_value(const dvector&, const int&);
01419   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01420 public:
01421 #if defined(USE_ADPVM)
01422   void pvm_pack(void){::pvm_pack(*this);}
01423   void pvm_unpack(void){::pvm_unpack(*this);}
01424 #endif
01425   virtual void set_value(const dvar_vector& x, const int& ii,
01426     const dvariable& pen);
01427 #if defined(USE_SHARE_FLAGS)
01428   virtual void shared_set_value_inv(const dvector&,const int&);
01429   virtual void shared_set_value(const dvar_vector&,const int&,
01430     const dvariable& pen);
01431 #endif
01432   virtual void set_value_inv(const dvector& x, const int& ii);
01433   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01434   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01435   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01436   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01437   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01438     const int& ii);
01439 
01440 public:
01441  void allocate(const ad_integer& imin,
01442    const ad_integer& imax, const ad_integer& imin2,
01443    const ad_integer& imax2, const ad_double& _bmin,
01444    const ad_double& _bmax, const ad_integer& phase_start,
01445    const char * s);
01446 
01447   param_init_bounded_matrix(void);
01448   void allocate(int rmin,int rmax,int cmin,int cmax,
01449     double _minb,double _maxb,
01450     int phase_start=1,const char * = "UNNAMED");
01451   void allocate(int rmin,int rmax,int cmin,int cmax,
01452     double _minb,double _maxb,const char * = "UNNAMED");
01453 
01454   void allocate(const ad_integer& rmin, const ad_integer& rmax,
01455     const index_type& cmin,
01456                 const index_type& cmax, double _minb, double _maxb,
01457                 const char * = "UNNAMED");
01458   void allocate(const ad_integer& rmin, const ad_integer& rmax,
01459     const index_type& cmin,
01460                 const index_type& cmax, double _minb, double _maxb,
01461                 int phase_start = 1, const char * = "UNNAMED");
01462 };
01463 
01468 class data_int: public model_name_tag
01469 {
01470 protected:
01471   int val;
01472 
01473   data_int& operator=(const int);
01474 
01475   void allocate(int n,const char * s="UNNAMED");
01476   void allocate(const char * s="UNNAMED");
01477   void allocate(init_xml_doc&, const char * s="UNNAMED");
01478 
01479 public:
01480   virtual ~data_int(){;}
01481   operator int() const { return val; }
01482 
01483   friend class model_data;
01484   friend class model_parameters;
01485   friend int operator + (int n,data_int v);
01486   friend int operator + (data_int v,int n);
01487   friend int operator + (data_int v,data_int n);
01488 };
01489 
01494 class named_adstring : public adstring, public model_name_tag
01495 {
01496 protected:
01497   void allocate(const char * s1,const char * s="UNNAMED");
01498   void operator = (const adstring&);
01499   void operator = (const char *);
01500 };
01501 
01506 class named_line_adstring : public line_adstring, public model_name_tag
01507 {
01508 protected:
01509   void allocate(const char * s1,const char * s="UNNAMED");
01510   void operator = (const adstring&);
01511   void operator = (const char *);
01512 };
01513 
01518 class init_adstring: public named_adstring
01519 {
01520 public:
01521   void allocate(const char * s="UNNAMED");
01522 };
01523 
01528 class init_line_adstring: public named_line_adstring
01529 {
01530 public:
01531   void allocate(const char * s="UNNAMED");
01532 };
01533 
01538 class dll_named_adstring : public named_adstring
01539 {
01540   char ** d;
01541 public:
01542   void allocate(char ** ps1,const char * s="UNNAMED");
01543   void operator = (const adstring&);
01544   void operator = (const char *);
01545   ~dll_named_adstring();
01546   dll_named_adstring(void){d=NULL;}
01547 };
01548 
01553 class dll_data_int : public data_int
01554 {
01555 public:
01556   int *pi;
01557   void allocate(int *_pi,const char * s);
01558   virtual ~dll_data_int();
01559 };
01560 
01565 class data_matrix : public named_dmatrix
01566 {
01567 public:
01568   data_matrix(void) : named_dmatrix() {;}
01569   data_matrix& operator=(const dmatrix& m);
01570   data_matrix& operator=(const double& m);
01571 private:
01572   void allocate(int rmin,int rmax,int cmin,int cmax,
01573     const char * = "UNNAMED");
01574   void allocate(int rmin, int rmax, const ivector& cmin, const ivector& cmax,
01575     const char * = "UNNAMED");
01576   void allocate(int rmin, int rmax, const ivector& cmin, int cmax,
01577     const char * = "UNNAMED");
01578   void allocate(int rmin, int rmax, int cmin, const ivector& cmax,
01579     const char * = "UNNAMED");
01580   void allocate(init_xml_doc&, const char * = "UNNAMED");
01581   friend class model_data;
01582 };
01583 
01588 class dll_data_matrix : public data_matrix
01589 {
01590   double * d;
01591 public:
01592   void allocate(double * _d,int rmin,int rmax,int cmin,int cmax,
01593     const char * _s = "UNNAMED");
01594   virtual ~dll_data_matrix();
01595   dll_data_matrix& operator=(const dmatrix &m);
01596   dll_data_matrix& operator=(const double &m);
01597 };
01598 
01603 class data_3array : public named_d3_array
01604 {
01605 public:
01606   data_3array(void) : named_d3_array() {;}
01607 private:
01608   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
01609     const char * ="UNNAMED");
01610   void allocate(int hsl, int hsu, const ivector& rmin, int rmax,
01611     int cmin,int cmax,const char * ="UNNAMED");
01612   void allocate(int hsl, int hsu, int rmin, const ivector& rmax,
01613     int cmin,int cmax,const char * ="UNNAMED");
01614   void allocate(int hsl, int hsu, const ivector& rmin, const ivector& rmax,
01615     int cmin,int cmax,const char * ="UNNAMED");
01616   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
01617     int cmax, const char * ="UNNAMED");
01618   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
01619                 const ivector& cmax, const char * ="UNNAMED");
01620   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,
01621                 const ivector& cmax, const char * ="UNNAMED");
01622   void allocate(int hsl,int hsu, const index_type& rmin, const index_type& rmax,
01623     const index_type& cmin, const index_type& cmax, const char * ="UNNAMED");
01624   friend class model_data;
01625 };
01626 
01627 class data_3iarray : public named_i3_array
01628 {
01629   data_3iarray(void) : named_i3_array() {;}
01630   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
01631     const char * ="UNNAMED");
01632   void allocate(int hsl, int hsu,const index_type& rmin, const index_type& rmax,
01633     const index_type& cmin, const index_type& cmax, const char * ="UNNAMED");
01634   friend class model_data;
01635 };
01636 
01641 class dll_data_3array : public data_3array
01642 {
01643   double * d;
01644 public:
01645   void allocate(double * _d,int hmin,int hmax,int rmin,int rmax,
01646     int cmin,int cmax,const char * _s = "UNNAMED");
01647   dll_data_3array& operator=(const d3_array &);
01648   virtual ~dll_data_3array();
01649   friend class model_data;
01650 };
01651 
01652 
01657 class data_5array : public named_d5_array
01658 {
01659   data_5array(void) : named_d5_array() {;}
01660   void allocate(int hhsl,int hhsu,
01661     int hhhsl,int hhhsu,
01662     int hsl,int hsu,int rmin,int rmax,
01663     int cmin,int cmax,const char * ="UNNAMED");
01664   void allocate(ad_integer hhhsl, ad_integer hhhsu, const index_type& hhsl,
01665     const index_type& hhsu, const index_type& hsl, const index_type& hsu,
01666     const index_type& rmin, const index_type& rmax, const index_type& cmin,
01667     const index_type& cmax, const char * ="UNNAMED");
01668   friend class model_data;
01669 };
01670 
01675 class data_4array : public named_d4_array
01676 {
01677   data_4array(void) : named_d4_array() {;}
01678   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,int rmax,
01679     int cmin,int cmax,const char * ="UNNAMED");
01680   void allocate(ad_integer hhsl, ad_integer hhsu, const index_type& hsl,
01681     const index_type& hsu, const index_type& rmin, const index_type& rmax,
01682     const index_type& cmin, const index_type& cmax, const char * ="UNNAMED");
01683   friend class model_data;
01684 };
01685 
01686 class data_4iarray : public named_i4_array
01687 {
01688   data_4iarray(void) : named_i4_array() {;}
01689   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,int rmax,
01690     int cmin,int cmax,const char * ="UNNAMED");
01691   void allocate(ad_integer hhsl, ad_integer hhsu, const index_type& hsl,
01692     const index_type& hsu, const index_type& rmin, const index_type& rmax,
01693     const index_type& cmin, const index_type& cmax, const char * ="UNNAMED");
01694   friend class model_data;
01695 };
01696 
01701 class data_imatrix : public named_imatrix
01702 {
01703   data_imatrix(void) : named_imatrix() {;}
01704   void allocate(int rmin,int rmax,int cmin,int cmax,const char * ="UNNAMED");
01705   void allocate(int rmin, int rmax, const index_type&, const index_type& cmax,
01706     const char * ="UNNAMED");
01707   friend class model_data;
01708 };
01709 
01714 class data_vector : public named_dvector
01715 {
01716 public:
01717   data_vector& operator=(const dvector& m);
01718   data_vector& operator=(const double m);
01719   data_vector(void) : named_dvector() {;}
01720 private:
01721   void allocate(int imin,int imax,const char * ="UNNAMED");
01722   void allocate(int imin, const ivector& imax, const char * ="UNNAMED");
01723   void allocate(init_xml_doc&, const char * ="UNNAMED");
01724   friend class model_data;
01725 };
01726 
01731 class dll_data_vector : public data_vector
01732 {
01733 public:
01734   double * pd;
01735   void allocate(double * pd,int imin,int imax,const char * ="UNNAMED");
01736   void allocate(double *pd, int imin, const ivector& imax,
01737     const char * ="UNNAMED");
01738   virtual ~dll_data_vector();
01739   dll_data_vector& operator = (const dvector& x);
01740   dll_data_vector& operator = (const double& x);
01741 };
01742 
01747 class data_ivector : public named_ivector
01748 {
01749 public:
01750   data_ivector(void) : named_ivector() {;}
01751 private:
01752   void allocate(int imin,int imax,const char * ="UNNAMED");
01753   friend class model_data;
01754 };
01755 
01760 class data_number : public model_name_tag
01761 {
01762 protected:
01763   double val;
01764   void allocate(const char * ="UNNAMED");
01765 public:
01766   void report_value(void);
01767   operator double() {return val;}
01768   double& value(void) {return val;}
01769   void initialize(void) {val=0.0;}
01770   friend class model_data;
01771   data_number & operator=(const double& m);
01772 };
01773 
01778 class dll_data_number : public data_number
01779 {
01780 public:
01781   double * pd;
01782   void allocate(double *_pd,const char * s);
01783   virtual ~dll_data_number();
01784   dll_data_number & operator=(const double& m);
01785 };
01786 
01787 typedef dvariable (model_parameters::*PMF) (const dvariable&);
01788 typedef dvariable (model_parameters::*PMFI) (const dvariable&,int n);
01789 typedef dvariable (model_parameters::*PMFVI) (const dvar_vector&,int n);
01790 typedef void (model_parameters::*PMFVIV4) (const dvar_vector&,int n,
01791   dvariable& f1, const dvariable& f2, const dvariable& f3, const dvariable& f4);
01792 
01793   class init_df1b2vector;
01794   class df1b2vector;
01795   class df1b2variable;
01796 
01801 class function_minimizer
01802 {
01803 public:
01804   static int bad_step_flag;
01805   static int likeprof_flag;
01806   static int first_hessian_flag;
01807   static int test_trust_flag;
01808   static int random_effects_flag;
01809   dmatrix * negdirections;
01810   static int negative_eigenvalue_flag;
01811   static int inner_opt_flag;
01812   static int inner_opt(void);
01813   laplace_approximation_calculator * lapprox;
01814   dvector * multinomial_weights;
01815   void set_multinomial_weights(dvector&d);
01816   //init_df1b2vector* py;
01817   virtual void AD_uf_inner();
01818   virtual void AD_uf_outer();
01819   virtual void user_function();
01820   void pre_user_function(void);
01821   //void df1b2_pre_user_function(void);
01822   //virtual void user_function(const init_df1b2vector& x,df1b2variable& f);
01823   //static int hesstype;
01824   //int set_hessian_type(int);
01825   void hess_routine_noparallel_randomeffects(void);
01826   void other_separable_stuff_begin(void);
01827   void other_separable_stuff_end(void);
01828   void begin_df1b2_funnel(void){;}
01829   void setup_quadprior_calcs(void){;}
01830   void end_df1b2_funnel(void){;}
01831   void get_function_difference(void);
01832   void start_get_importance_sampling_comnponent(void);
01833   void end_get_importance_sampling_comnponent(void);
01834   int spminflag;
01835   int repeatminflag;
01836   int mcmc2_flag;
01837   int robust_hybrid_flag;
01838   long ifn;
01839   int maxfn;
01840   int iprint;
01841   double crit;
01842   int imax;
01843   double dfn;
01844   long iexit;
01845   long ihflag;
01846   long ihang;
01847   int scroll_flag;
01848   int maxfn_flag;
01849   int quit_flag;
01850   double min_improve;
01851   void pre_userfunction(void);
01852   virtual void userfunction(void)=0;
01853   virtual void allocate(void){;}
01854   static named_dvar_vector * ph;
01855   static named_dvar_vector * pg;
01856 protected:
01857   double ffbest;
01858 private:
01859   gradient_structure * pgs;
01860   adstring_array param_labels;
01861   ivector param_size;
01862 protected:
01863   void report_function_minimizer_stats(void){;}
01864   virtual void report(const dvector& gradients){;};
01865   static dvector convergence_criteria;
01866   static dvector maximum_function_evaluations;
01867   static int sd_flag;
01868   static adstring user_data_file;
01869   static adstring user_par_file;
01870   static int have_constraints;
01871 public:
01872   virtual dvariable user_randeff(const dvar_vector& x);
01873   virtual dvar_vector user_dfrandeff(const dvar_vector& x);
01874   virtual dvar_matrix user_d2frandeff(const dvar_vector& x);
01875   void limited_memory_quasi_newton(const independent_variables&,int);
01876   void limited_memory_quasi_newton(const independent_variables&,int,int);
01877   void limited_memory_quasi_newton(double& f, const independent_variables&,
01878     int, int, int,double);
01879   function_minimizer(long int sz=0L);
01880   void likeprof_routine(double global_min);
01881   virtual ~function_minimizer();
01882   virtual void other_calculations(void){;}
01883   virtual void final_calcs(void){;}
01884   virtual void minimize(void);
01885   virtual void constraints_minimize(void);
01886   virtual void between_phases_calculations(void){;}
01887   void computations(int argc,char * argv[]);
01888   void computations1(int argc,char * argv[]);
01889   void computations_np(int argc,char * argv[]);
01890   void computations(void);
01891   void hess_routine(void);
01892   void hess_routine_noparallel(void);
01893   void hess_routine_master(void);
01894   void hess_routine_slave(void);
01895   void constraint_report(void);
01896   void pvm_slave_likeprof_routine(void);
01897   void pvm_master_function_evaluation_profile(double& f,
01898     independent_variables& x,const dvector & g,int nvar,int iprof,double weight,
01899     double new_value);
01900   void pvm_slave_prof_minimize(int underflow_flag);
01901   void pvm_master_prof_minimize(int iprof, double sigma,
01902     double new_value, const double& _fprof, const int underflow_flag,
01903     double global_min, const double& _penalties,
01904     const double& _final_weight);
01905   //void pvm_master_function_evaluation_profile(double& f,
01906    // independent_variables& x,dvector & g,int nvar,int iprof);
01907   void pvm_slave_function_evaluation(void);
01908   void pvm_slave_function_evaluation_no_derivatives(void);
01909   void pvm_slave_function_evaluation_noder(void);
01910   void pvm_master_function_evaluation_no_derivatives(double& f,
01911     independent_variables& x,int nvar);
01912   void pvm_master_function_evaluation(double& f,
01913     independent_variables& x,const dvector & g,int nvar);
01914   dmatrix dep_hess_routine(const dvariable& dep);
01915   void top_mcmc_routine(int,int,double,int);
01916   void mcmc_routine(int,int,double,int);
01917   void sgibbs_mcmc_routine(int,int,double,int);
01918   void hybrid_mcmc_routine(int,int,double,int);
01919   double pvm_master_get_monte_carlo_value(int nvar,
01920     const dvector& x);
01921   void pvm_slave_get_monte_carlo_value(int nvar);
01922   void mcmc_eval(void);
01923   //void hess_routine_and_constraint(int);
01924   void hess_routine_and_constraint(int iprof, const dvector& g,
01925     dvector& fg);
01926   dmatrix diag_hess_routine(void);
01927   void hess_inv(void);
01928   void depvars_routine(void);
01929   void sd_routine(void);
01930   int ef_(double * f, double * x);
01931   int constrained_minimization2(int _n,int _nh, int _ng,dvector& __x);
01932   static int constraint_exit_number;
01933   void get_bigS(int ndvar,int nvar1,int nvar,
01934     dmatrix& S,dmatrix& BS,dvector& scale);
01935 
01936 #ifdef CURVE_CORRECT
01937   void constraint_hess_routine(int ip);
01938   void get_curvature_correction_factors(int ip,
01939     dvector& g, const int underflow_flag, const dvector& eigenvals,
01940     dvector& curvcor);
01941 #endif
01942   double projected_jacobian(const dvector& g, const dvector& xscale);
01943 
01944   dvariable random_effects_maximization(const dvar_vector& v);
01945   void prof_minimize_re(int iprof, double sigma,
01946     double new_value, const double& fprof, const int underflow_flag,
01947     double global_min, const double& penalties, const double& final_weight);
01948   void prof_minimize(int iprof, double sigma,
01949     double new_value, const double& fprof, const int underflow_flag,
01950     double global_min, const double& penalties, const double& final_weight);
01951   void begin_gauss_hermite_stuff(void);
01952   void begin_funnel_stuff(void);
01953   void end_gauss_hermite_stuff(void);
01954 
01955   void prof_hess_routine(int ip,double new_value,double final_weight);
01956 
01957   void quasi_newton_minimizer1(int nvar,double _crit,
01958     double& f, const independent_variables& x,const dvector& g);
01959 
01960   double hess_determinant(int underflow_flag);
01961 
01962 #ifndef CURVE_CORRECT
01963   void normalize_posterior_distribution(double udet,
01964     const dvector& siglevel, const ofstream& ofs2, int num_pp,
01965     const dvector& all_values, const dmatrix& actual_value, double global_min,
01966     int offset, const dmatrix& lprof, const dmatrix& ldet, const dmatrix& xdist,
01967     const dmatrix& penalties);
01968   //  dmatrix& penalties, const hdmatrix& lg_jacob);
01969 #else
01970   void normalize_posterior_distribution(double udet,
01971     const dvector& siglevel, const ofstream& ofs2, int num_pp,
01972     const dvector& all_values, const dmatrix& actual_value,
01973     double global_min,
01974     int offset, const dmatrix& lprof, const dmatrix& ldet,
01975     const dmatrix& xdist,
01976     const d3_array& eigenvals,d3_array& curvcor);
01977 #endif
01978   void get_particular_grad(int iprof, int nvar, const dvector& fg,
01979     const dvector& g);
01980   double projected_hess_determinant(const dvector& g, const int underflow_flag,
01981     const dvector& xscale, const double& ln_det_proj_jac);
01982 //  double projected_hess_determinant(const dvector& fg, const dvector& g,
01983   //const int underflow_flag, const dvector& xscale,
01984   // const double& ln_det_proj_jac);
01985   double projected_hess_determinant(const dvector& g,const int underflow_flag);
01986   double projected_hess_determinant(const dmatrix& hh, const dvector& g,
01987     const int underflow_flag);
01988  //double projected_hess_determinant(const dvector& g, const int underflow_flag,
01989     //dvector& xscale, const double& ln_det_proj_jac);
01990 
01991   double projected_hess_determinant(const dvector& fg, const dvector& g,
01992     const int underflow_flag, const dvector& curvscale, const dvector& xscale,
01993     double& ln_det_proj_jac, const double& tmp, const dmatrix& hesses);
01994   double diag_projected_hess_determinant(const dvector& g,
01995     const int underflow_flag, dmatrix& dh);
01996     double unrestricted_hess_determinant(void);
01997   void monte_carlo_routine(void);
01998   double get_monte_carlo_value(int nvar, const independent_variables& x);
01999   double get_monte_carlo_value(int nvar, const independent_variables& x,
02000     dvector& g);
02001   double get_hybrid_monte_carlo_value(int nvar,const independent_variables& x,
02002     dvector& g);
02003   void mcmc_computations(void);
02004 #if defined(USE_ADPVM)
02005   void pvm_slave_mcmc_computations(void);
02006   void pvm_master_mcmc_computations(void);
02007 #endif
02008   double get_monte_carlo_value(int nvar, const dvector& x);
02009   void sob_routine(int nmcmc,double dscale,int restart_flag);
02010   void sobol_importance_routine(int nmcmc,int iseed0,double dscale,
02011   int restart_flag);
02012   void pvm_master_mcmc_routine(int nmcmc,int iseed0,double dscale,
02013     int restart_flag);
02014 #if defined(USE_ADPVM)
02015   void pvm_slave_mcmc_routine(void);
02016 #else
02017   void pvm_slave_mcmc_routine(void) {}
02018 #endif
02019   void trust_region_update(int nvar,int _crit,
02020     independent_variables& x,const dvector& _g,const double& _f);
02021 
02022   void multint4(int n, const dvar_vector& a, const dvar_vector& b,
02023     const dvar_vector& h, double al, int m, const dvariable& e,
02024     const dvariable& aint1, const dvariable& aint2, dvariable& aint3,
02025     const dvariable& aint4, const int& key, PMFVIV4 f);
02026 
02027   void multint(int n, const dvar_vector& a, const dvar_vector& b,
02028     const dvar_vector& h, double al, int m, const dvariable& e,
02029     const dvariable& aint, const int& key, PMFVI f);
02030 
02031   virtual void set_runtime(void);
02032   virtual void set_runtime_maxfn(const char *);
02033   virtual void set_runtime_crit(const char *);
02034   dvariable adromb(PMF,double a,double b,int ns=9);
02035   dvariable adromb(PMF, const dvariable& a, double b, int ns = 9);
02036   dvariable adromb(PMF, const dvariable& a, const dvariable& b, int ns = 9);
02037   dvariable adromb(PMF, double a, const dvariable& b, int ns = 9);
02038 
02039   dvariable adrombo(PMF,double a,double b,int ns=9);
02040   dvariable adrombo(PMF, const dvariable& a, double b, int ns = 9);
02041   dvariable adrombo(PMF, const dvariable& a, const dvariable& b, int ns = 9);
02042   dvariable adrombo(PMF, double a, const dvariable& b,int ns = 9);
02043 
02044   dvariable trapzd(void*,double a,double b,int n);
02045   dvariable trapzd(PMF,double a,double b,int n);
02046   dvariable trapzd(PMF, double a, const dvariable& b, int n);
02047   dvariable trapzd(PMF, const dvariable& a, double b, int n);
02048   dvariable trapzd(PMF, const dvariable& a, const dvariable& b, int n);
02049 
02050   dvariable midpnt(PMF,double a,double b,int n);
02051   dvariable midpnt(PMF, double a, const dvariable& b, int n);
02052   dvariable midpnt(PMF, const dvariable& a, double b, int n);
02053   dvariable midpnt(PMF, const dvariable& a, const dvariable& b, int n);
02054 
02055   virtual void * mycast() { return (void*)this;}
02056 
02057   void neldmead(int n, dvector& _start, dvector& _xmin, double *ynewlo,
02058     double reqmin, double delta,int *icount, int *numres, int *ifault);
02059   void adamoeba(const dmatrix& p, const dvector& y, int ndim, double ftol,
02060      int maxfn);
02061   void set_initial_simplex(const dmatrix& p, const dvector& y, int nvar,
02062     const dvector& x, double delta);
02063   double amxxx(const dmatrix& p, const dvector& y, const dvector& psum,
02064     int ndim, int ihi, double fac);
02065   friend class equality_constraint_vector;
02066   friend class inequality_constraint_vector;
02067   void quasi_newton_block(int nvar,int crit,independent_variables& x,
02068     const dvector& g,const double& f);
02069   void limited_memory_quasi_newton_block(int nvar,int _crit,
02070     independent_variables& x,const dvector& _g,const double& _f,
02071     int nsteps);
02072 
02073   void function_evaluation_block_pvm_slave_random_effects(int nvar,int _crit,
02074     independent_variables& x,const dvector& g,const double& f);
02075   void quasi_newton_block_pvm_master_random_effects(int nvar,int _crit,
02076     independent_variables& x,const dvector& g,const double& f);
02077   void function_evaluation_block_pvm_slave_random_effects(void);
02078   void hess_routine_random_effects(void);
02079   void quasi_newton_block_pvm_master(int nvar,int _crit,
02080     independent_variables& x,const dvector& g,const double& f);
02081   void hess_routine_noparallel_random_effects(void);
02082 #if defined(USE_ADPVM)
02083   void function_evaluation_block_pvm_slave(void);
02084   void hess_routine_slave_random_effects(void);
02085 #endif
02086   dvariable do_gauss_hermite_integration(void);
02087   void end_df1b2_funnel_stuff(void);
02088 
02089 private:
02090   dvariable do_gauss_hermite_integration_multi(void);
02091 };
02092 
02093 cifstream& operator>>(const cifstream& s, const param_init_number& x);
02094 cifstream& operator>>(const cifstream& s, const param_init_vector& v);
02095 cifstream& operator>>(const cifstream& s, const param_init_matrix& m);
02096 ostream& operator<<(const ostream& s, const label_class& lc);
02097 
02102 class stddev_params
02103 {
02104 public:
02105   // this should be a resizeable array
02106   static stddev_params * stddevptr[150];
02107   // this should be a resizeable array
02108   static stddev_params * stddev_number_ptr[150];
02109   static void get_all_sd_values(const dvector& x, const int& ii);
02110   static int num_stddev_params;
02111   static int num_stddev_number_params;
02112   static ivector copy_all_number_offsets(void);
02113   void allocate(void){;};
02114   static int num_stddev_calc(void);
02115   static int num_stddev_number_calc(void);
02116   static void get_stddev_number_offset(void);
02117 public:
02118   stddev_params(void){;}
02119   virtual void setindex(int);
02120   virtual int getindex(void);
02121   virtual int size_count(void)=0; // get the number of active parameters
02122   virtual void set_dependent_variables(void)=0;
02123   virtual void copy_value_to_vector(const dvector&,const int&) = 0;
02124   virtual void get_sd_values(const dvector& x, const int& ii) = 0;
02125   //get the number of active parameters
02126   static void copy_all_values(const dvector& x, const int& ii);
02127   //get the number of active parameters
02128   static void copy_all_number_values(const dvector& x, const int& ii);
02129   virtual void add_to_list(void);
02130   virtual void add_to_gui_list(void);
02131   virtual const char * label()=0;
02132   friend class function_minimizer;
02133 };
02134 
02139 class likeprof_params
02140 {
02141   double stepsize;
02142   int    stepnumber;
02143 protected:
02144 public:
02145   static likeprof_params * likeprofptr[500]; // this should be a
02146                                                // resizeable array
02147   static int num_likeprof_params;
02148   void allocate(void){;};
02149   static int num_stddev_calc(void);
02150 public:
02151   likeprof_params(void);
02152   virtual void add_to_list(void);
02153   virtual const char * label()=0;
02154   virtual dvariable variable(void)=0;
02155   virtual double& get_sigma(void)=0;
02156   virtual double get_value(void)=0;
02157   double get_stepsize(void);
02158   int get_stepnumber(void);
02159   void set_stepsize(double);
02160   void set_stepnumber(int);
02161   friend class function_minimizer;
02162 };
02163 
02168 class param_stddev_vector: public named_dvar_vector , stddev_params
02169 {
02170   dvector sd;
02171     virtual int size_count(void); // get the number of active parameters
02172     virtual const char * label(void);
02173     param_stddev_vector();
02174     void allocate(int imin,int imax,const char * s="UNNAMED");
02175     virtual void set_dependent_variables(void);
02176     friend class model_parameters;
02177     virtual void copy_value_to_vector(const dvector& x, const int& ii);
02178     virtual void get_sd_values(const dvector& x, const int& ii);
02179   param_stddev_vector& operator=(const dvar_vector& m);
02180   param_stddev_vector& operator=(const dvector& m);
02181   param_stddev_vector& operator=(const double m);
02182 };
02183 
02188 class param_stddev_number: public named_dvariable , public stddev_params
02189 {
02190   double sd;
02191   int index;
02192   void allocate(const char *s="UNNAMED");
02193   virtual void setindex(int);
02194   virtual int getindex(void);
02195   virtual int size_count(void); // get the number of active parameters
02196   virtual const char * label(void);
02197   virtual void copy_value_to_vector(const dvector& x, const int& ii);
02198   virtual void get_sd_values(const dvector& x, const int& ii);
02199 protected:
02200   param_stddev_number();
02201   friend class model_parameters;
02202   virtual void set_dependent_variables(void);
02203   param_stddev_number& operator=(const prevariable&);
02204   param_stddev_number& operator=(const double);
02205 };
02206 
02211 class param_likeprof_number: public param_stddev_number ,
02212     public likeprof_params
02213 {
02214     double sigma;
02215     void allocate(const char *s="UNNAMED");
02216     virtual int size_count(void); // get the number of active parameters
02217     virtual const char * label(void);
02218     virtual double& get_sigma(void){return sigma;}
02219     virtual double get_value(void){return value(*this);}
02220     //void copy_value_to_vector(const dvector& x, const int& ii);
02221     virtual dvariable variable(void){ return dvariable(*this);}
02222     param_likeprof_number();
02223     friend class model_parameters;
02224 public:
02225     param_likeprof_number& operator=(const prevariable&);
02226     param_likeprof_number& operator=(const double);
02227 };
02228 
02233 class param_stddev_matrix: public named_dvar_matrix , stddev_params
02234 {
02235   dmatrix sd;
02236     virtual int size_count(void);
02237     //virtual void read_value(void);
02238     virtual const char * label(void);
02239     void allocate(int rmin,int rmax,int cmin,int cmax,
02240     const char * s="UNNAMED");
02241   param_stddev_matrix(void);
02242   friend class model_parameters;
02243   virtual void set_dependent_variables(void);
02244   virtual void get_sd_values(const dvector& x, const int& ii);
02245   virtual void copy_value_to_vector(const dvector& x, const int& ii);
02246   param_stddev_matrix& operator=(const double m);
02247   param_stddev_matrix& operator=(const dmatrix& m);
02248   param_stddev_matrix& operator=(const dvar_matrix& m);
02249 };
02250 
02255   class objective_function_value : public named_dvariable
02256   {
02257   public:
02258     static objective_function_value * pobjfun;
02259     static double fun_without_pen;
02260     static double gmax;
02261     objective_function_value();
02262     objective_function_value& operator=(const prevariable& v);
02263     objective_function_value& operator=(const double v);
02264   };
02265 
02266   int withinbound(int lb,int n,int ub);
02267 
02268 double cumd_cauchy(const double& x);
02269 double density_cauchy(const double& x);
02270 double log_density_cauchy(const double& x);
02271 double inv_cumd_cauchy(const double& x);
02272 
02273 double cumd_mixture(const double& x);
02274 double inv_cumd_mixture(const double& y);
02275 double cumd_mixture_02(const double& x);
02276 double inv_cumd_mixture_02(const double& y);
02277 
02278 #if defined _ADM_HIGHER_ARRAYS__
02279 
02284 class param_init_matrix: public named_dvar_matrix, public initial_params
02285 {
02286   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
02287 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
02288 //  const int& ii);
02289   virtual void add_value(const dvector&, const dvector&, const int&,
02290     const double&, const dvector&);
02291   virtual void add_value(const dvector&, const int&);
02292   virtual void get_jacobian(const dvector&, const dvector&, const int&);
02293 public:
02294   virtual void set_value(const dvar_vector& x, const int& ii,
02295     const dvariable& pen);
02296   virtual void copy_value_to_vector(const dvector& x, const int& ii);
02297   virtual void restore_value_from_vector(const dvector&, const int&);
02298   virtual void set_value_inv(const dvector& x, const int& ii);
02299   virtual int size_count(void);
02300   virtual void save_value(ofstream& ofs);
02301   virtual void bsave_value(uostream& uos);
02302   virtual void save_value(const ofstream& ofs, int prec);
02303   virtual void restore_value(const ifstream& ifs);
02304   void report_value(void);
02305   //virtual void read_value(void);
02306   virtual const char * label(void);
02307   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
02308   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
02309   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
02310   virtual void hess_scale(const dvector& d, const dvector& x, const int& ii){;};
02311 
02312 public:
02313   void allocate(int rmin,int rmax,int cmin,int cmax,
02314     int phase_start=1,const char * = "UNNAMED");
02315   void allocate(int rmin,int rmax,int cmin,int cmax,
02316     const char * = "UNNAMED");
02317   param_init_matrix(void);
02318   param_init_matrix& operator = (const dmatrix& m);
02319   param_init_matrix& operator = (const dvar_matrix& m);
02320   param_init_matrix& operator = (const dvariable& m);
02321   param_init_matrix& operator = (const double& m);
02322 };
02323 #endif // #if defined _ADM_HIGER_ARRAYS__
02324 
02329 class param_init_d3array: public named_dvar3_array, public initial_params
02330 {
02331 public:
02332 #if defined(USE_SHARE_FLAGS)
02333   virtual void setshare(const index_type& sf,const index_type& af);
02334   virtual void shared_set_value_inv(const dvector&,const int&);
02335   virtual void shared_set_value(const dvar_vector&,const int&,
02336     const dvariable& pen);
02337   virtual int shared_size_count(void); // get the number of active parameters
02338 #endif
02339 
02340   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
02341     const int& ii);
02342   virtual void dev_correction(const dmatrix&, const int&);
02343   virtual void curv_scale(const dvector& d, const dvector& x, const int& ii);
02344   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
02345 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
02346 //  const int& ii);
02347   virtual void add_value(const dvector&, const dvector&, const int&,
02348     const double&, const dvector&);
02349   virtual void add_value(const dvector&, const int&);
02350   virtual void get_jacobian(const dvector&, const dvector&, const int&);
02351   virtual void set_value(const dvar_vector& x, const int& ii,
02352     const dvariable& pen);
02353   virtual void copy_value_to_vector(const dvector& x, const int& ii);
02354   virtual void restore_value_from_vector(const dvector&,const int&);
02355   virtual void set_value_inv(const dvector& x, const int& ii);
02356   virtual int size_count(void);
02357   virtual void save_value(ofstream& ofs);
02358   virtual void bsave_value(uostream& uos);
02359   virtual void save_value(const ofstream& ofs, int prec);
02360   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
02361     int& offset){}
02362   virtual void restore_value(const ifstream& ifs);
02363   void report_value(void);
02364   //virtual void read_value(void);
02365   virtual const char * label(void);
02366   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
02367   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
02368   virtual void hess_scale(const dvector& d, const dvector& x, const  int& ii);
02369 
02370 public:
02371 #if defined(USE_ADPVM)
02372   virtual void pvm_pack(void) { cerr << "Error" << endl; ad_exit(1);}
02373   virtual void pvm_unpack(void) { cerr << "Error" << endl; ad_exit(1);}
02374 #endif
02375 
02376   void allocate(const ad_integer& sl,const ad_integer& sh,
02377     const index_type& nrl,const index_type& nrh,
02378     const index_type& ncl,const index_type& nch,const char * s="UNNAMED");
02379 
02380   void allocate(const ad_integer& sl,const ad_integer& sh,
02381     const index_type& nrl,const index_type& nrh,
02382     const index_type& ncl,const index_type& nch,int phase_start=1,
02383     const char * s="UNNAMED");
02384 
02385   void allocate(int smin,int smax,int rmin,int rmax,int cmin,int cmax,
02386     int phase_start=1,const char * = "UNNAMED");
02387   void allocate(int smin,int smax,int rmin,int rmax,int cmin,int cmax,
02388     const char * = "UNNAMED");
02389   param_init_d3array(void);
02390 };
02391 
02396 class dll_param_init_d3array: public param_init_d3array
02397 {
02398   double * d;
02399 public:
02400   void allocate(double* _d,int hmin,int hmax,
02401     int rmin,int rmax,int cmin,int cmax,
02402     int phase_start=1,const char * = "UNNAMED");
02403   void allocate(double * _d,int hmin,int hmax,int rmin,int rmax,
02404     int cmin,int cmax,const char * = "UNNAMED");
02405   virtual ~dll_param_init_d3array();
02406   dll_param_init_d3array(){d=NULL;}
02407   dll_param_init_d3array& operator=(const d3_array&);
02408   dll_param_init_d3array& operator=(const dvar3_array&);
02409 };
02410 
02415 class dll_param_d3array: public named_dvar3_array
02416 {
02417   double * d;
02418 public:
02419   void allocate(double* _d,int hmin,int hmax,
02420     int rmin,int rmax,int cmin,int cmax,
02421     int phase_start=1,const char * = "UNNAMED");
02422   void allocate(double * _d,int hmin,int hmax,int rmin,int rmax,
02423     int cmin,int cmax,const char * = "UNNAMED");
02424   virtual ~dll_param_d3array();
02425   dll_param_d3array(){d=NULL;}
02426   dll_param_d3array& operator=(const d3_array&);
02427   dll_param_d3array& operator=(const dvar3_array&);
02428 };
02429 
02430 
02431 //double set_value_mc(const double& x, const double fmin, const double fmax);
02432 
02433 void set_value_mc(const dvar_vector& x, const dvar_vector& v, const int& ii,
02434   const double fmin, const double fmax);
02435 
02436 double set_value_inv_mc(double v,double fmin,double fmax);
02437 
02438 double set_value_inv_mc(const prevariable& v, double fmin, double fmax);
02439 
02440 void set_value_inv_mc(const dvar_vector& x, const dvector& v, const int& ii,
02441   const double fmin, const double fmax);
02442 
02443 //double set_value_inv_mc(const dvector&, const dvector& x, int ii, double minb,
02444 //  double maxb);
02445 
02446 double set_value_mc(double z,double min,double max);
02447 double ndfboundp( double x, double fmin, double fmax,const double& fpen);
02448 double ndfboundp_mc( double x, double fmin, double fmax,const double& fpen);
02449 
02450 void copy_value_from_vector(const double& _sd,const dvector& x,const int & _ii);
02451 void copy_value_from_vector(const dvector& _sd,const dvector& x,
02452   const int & _ii);
02453 void copy_value_from_vector(const dmatrix& _sd,const dvector& x,
02454   const int & _ii);
02455 
02456 dvector bounded_multivariate_normal(int nvar, const dvector& a1,
02457   const dvector& b1, dmatrix& ch, const double& _wght,
02458   const random_number_generator & rng);
02459 
02460 void bounded_multivariate_normal_mcmc(int nvar, const dvector& a1,
02461   const dvector& b1, dmatrix& ch, const double& wght, const dvector& y,
02462   const random_number_generator & rng);
02463 
02464 void bounded_multivariate_uniform_mcmc(int nvar, const dvector& a1,
02465   const dvector& b1, dmatrix& ch, const double& wght, const dvector& y,
02466   const random_number_generator & rng);
02467 
02468 dvector bounded_multivariate_normal(int nvar, const dvector& a1,
02469   const dvector& b1, dmatrix& ch, const double& lprob,
02470   random_number_generator &rng);
02471 
02472 dvector bounded_multivariate_normal_sobol(int nvar, const dvector& a1,
02473   const dvector& b1, dmatrix& ch, const double& lprob,
02474   const random_number_generator &rng);
02475 
02476 dvector probing_bounded_multivariate_normal(int nvar, const dvector& a1,
02477   const dvector& b1, dmatrix& ch, const double& lprob, double pprobe,
02478   const random_number_generator &rng);
02479 
02480 dvector bounded_multivariate_uniform(int nvar, const dvector& a1,
02481   const dvector& b1, dmatrix& ch, const double& lprob,
02482   random_number_generator &rng);
02483 
02484 dvector bounded_robust_multivariate_normal(int nvar, const dvector& a1,
02485   const dvector& b1, dmatrix& ch, const dmatrix& ch3, double contaminant,
02486   const double& _wght, random_number_generator &rng);
02487 
02488 void probing_bounded_multivariate_normal_mcmc(int nvar, const dvector& a1,
02489   const dvector& b1, dmatrix& ch, const double& wght, const dvector& v,
02490   double pprobe, const random_number_generator &rng);
02491 
02492 /*
02493 int option_match(int argc,char * argv[], const char * string);
02494 int option_match(int argc,char * argv[], const char * string, const int& nopt);
02495 int option_match(char * s, const char * string, const int& _nopt);
02496 int option_match(char * s, const char * string);
02497 */
02498 
02499 double inv_cumd_exp(double x);
02500 double cumd_exp(double x);
02501 
02502 double ffmax(double a,double b);
02503 double ffmin(double a,double b);
02504 
02505 void check_datafile_pointer(void * p);
02506 
02507 adstring get_reportfile_name(void);
02508 
02509 void ad_make_code_reentrant(void);
02510 
02511 char** parse_dll_options(char *pname, const int& argc,
02512   char * dll_options);
02513 
02514 void parse_dll_options(char *pname, const int& argc, char *dll_options,
02515   char *** pargv);
02516 
02517 char** no_dll_options(char *pname, const int& argc);
02518 
02519 void cleanup_argv(int nopt,char *** pa);
02520 
02521 void get_sp_printf(void);
02522 
02523 void do_dll_housekeeping(int argc,char ** argv);
02524 
02525 void adwait(double);
02526 
02527 int ad_get_commandline_option(const char *option_label, const int &option_value,
02528   const char * error_message);
02529 
02534 class param_init_vector_vector
02535 {
02536   param_init_vector * v;
02537   int index_min;
02538   int index_max;
02539   double_index_type * it;
02540 
02541 public:
02542   param_init_vector_vector();
02543   ~param_init_vector_vector();
02544 
02545   void set_scalefactor(double s);
02546   void set_scalefactor(const dvector& s);
02547   dvector get_scalefactor(void);
02548 
02549 #if defined(OPT_LIB)
02550   param_init_vector& operator [] (int i) { return v[i];}
02551   param_init_vector& operator () (int i) { return v[i];}
02552   prevariable operator () (int i,int j) { return v[i][j];}
02553 #else
02554   param_init_vector& operator [] (int i);
02555   param_init_vector& operator () (int i);
02556   prevariable operator () (int i,int j);
02557 #endif
02558 
02559   void allocate(int min1,int max1,const index_type& min,
02560      const index_type& max,const index_type& phase_start,
02561      const char * s);
02562 
02563   void allocate(int min1,int max1,const index_type& min,
02564      const index_type& max,const char * s);
02565 
02566   bool allocated() const { return v != NULL; }
02567   int indexmin(void) {return (index_min);}
02568   int indexmax(void) {return (index_max);}
02569   void set_initial_value(const double_index_type& it);
02570   void deallocate(void);
02571 };
02572 
02577 class param_init_bounded_vector_vector
02578 {
02579   param_init_bounded_vector* v;
02580   int index_min;
02581   int index_max;
02582   double_index_type* it;
02583 public:
02584   param_init_bounded_vector_vector();
02585   ~param_init_bounded_vector_vector();
02586 
02587   void set_scalefactor(double s);
02588   void set_scalefactor(const dvector& s);
02589   dvector get_scalefactor(void);
02590 #if defined(OPT_LIB)
02591   param_init_bounded_vector& operator [] (int i) { return v[i];}
02592   param_init_bounded_vector& operator () (int i) { return v[i];}
02593   prevariable operator () (int i,int j) { return v[i][j];}
02594 #else
02595   param_init_bounded_vector& operator [] (int i);
02596   param_init_bounded_vector& operator () (int i);
02597   prevariable operator () (int i,int j);
02598 #endif
02599 
02600   void allocate(int min1,int max1,
02601     const index_type& min,
02602     const index_type& max,
02603     const double_index_type& dmin,
02604     const double_index_type& dmax,
02605     const index_type& phase_start,
02606     const char* s);
02607 
02608   void allocate(int min1,int max1,
02609     const index_type& min,
02610     const index_type& max,
02611     const double_index_type& dmin,
02612     const double_index_type& dmax,
02613     const char* s);
02614 
02615   bool allocated() const { return v != NULL; }
02616   int indexmin() const {return index_min;}
02617   int indexmax() const {return index_max;}
02618   void deallocate(void);
02619   void set_initial_value(const double_index_type& it);
02620 };
02621 
02626 class param_init_matrix_vector
02627 {
02628   param_init_matrix* v;
02629   int index_min;
02630   int index_max;
02631   double_index_type* it;
02632 
02633 public:
02634   param_init_matrix_vector();
02635    ~param_init_matrix_vector();
02636 
02637   void set_scalefactor(double s);
02638   void set_scalefactor(const dvector& s);
02639   dvector get_scalefactor(void);
02640 
02641 #if defined(OPT_LIB)
02642   param_init_matrix& operator [] (int i) { return v[i];}
02643   param_init_matrix& operator () (int i) { return v[i];}
02644   dvar_vector& operator () (int i,int j) { return v[i][j];}
02645   prevariable operator () (int i,int j,int k) { return v[i](j,k);}
02646 #else
02647   param_init_matrix& operator [] (int i);
02648   param_init_matrix& operator () (int i);
02649   dvar_vector& operator () (int i,int j);
02650   prevariable operator () (int i,int j,int k);
02651 #endif
02652   void allocate(int min0,int max0,const index_type& min,
02653      const index_type& max,const index_type& min1,
02654      const index_type& max1,const index_type& phase_start,
02655      const char * s);
02656 
02657   void allocate(int min0,int max0,const index_type& min,
02658      const index_type& max,const index_type& min1,
02659      const index_type& max1,const char * s);
02660 
02661   bool allocated() const { return v != NULL; }
02662   int indexmin() const { return index_min;}
02663   int indexmax() const { return index_max;}
02664   void set_initial_value(const double_index_type& it);
02665   void deallocate(void);
02666 };
02667 
02672 class param_init_bounded_matrix_vector
02673 {
02674   param_init_bounded_matrix * v;
02675   int index_min;
02676   int index_max;
02677   double_index_type * it;
02678 public:
02679   param_init_bounded_matrix_vector();
02680   ~param_init_bounded_matrix_vector();
02681 
02682   void set_scalefactor(double s);
02683   void set_scalefactor(const dvector& s);
02684   dvector get_scalefactor(void);
02685 
02686   void allocate(int min1,int max1,
02687     const index_type& min, const index_type& max, const index_type& min2,
02688     const index_type& max2, const double_index_type& dmin2,
02689     const double_index_type& dmax2, const index_type& phase_start,
02690     const char * s);
02691 
02692   void allocate(int min1,int max1,
02693     const index_type& min, const index_type& max, const index_type& min2,
02694     const index_type& max2, const double_index_type& dmin2,
02695     const double_index_type& dmax2,const char * s);
02696 
02697 #if defined(OPT_LIB)
02698   param_init_bounded_matrix& operator [] (int i) { return v[i];}
02699   param_init_bounded_matrix& operator () (int i) { return v[i];}
02700   dvar_vector& operator () (int i,int j) { return v[i][j];}
02701   prevariable operator () (int i,int j,int k) { return v[i](j,k);}
02702 #else
02703   param_init_bounded_matrix& operator [] (int i);
02704   param_init_bounded_matrix& operator () (int i);
02705   dvar_vector& operator () (int i,int j);
02706   prevariable operator () (int i,int j,int k);
02707 #endif
02708 
02709   bool allocated() const { return v != NULL; }
02710   int indexmin() const { return index_min; }
02711   int indexmax() const { return index_max; }
02712   void set_initial_value(const double_index_type& it);
02713   void deallocate(void);
02714 };
02715 
02720 class param_init_number_vector
02721 {
02722   param_init_number* v;
02723   int index_min;
02724   int index_max;
02725   double_index_type* it;
02726 public:
02727   param_init_number_vector();
02728   ~param_init_number_vector();
02729 
02730   void set_scalefactor(double s);
02731   void set_scalefactor(const dvector& s);
02732   dvector get_scalefactor(void);
02733 
02734 #if defined(OPT_LIB)
02735   param_init_number& operator[](int i) { return v[i];}
02736   param_init_number& operator()(int i) { return v[i];}
02737 #else
02738   param_init_number& operator[](int i);
02739   param_init_number& operator()(int i);
02740 #endif
02741 
02742   void allocate(int min1,int max1,const index_type& phase_start,
02743      const char * s);
02744 
02745   void allocate(int min1,int max1,const char * s);
02746 
02747   bool allocated() const { return v != NULL; }
02748   int indexmin() const { return index_min; }
02749   int indexmax() const { return index_max; }
02750   void set_initial_value(const double_index_type& it);
02751   void deallocate(void);
02752 };
02753 
02761  class data_matrix;
02762 class param_init_bounded_number_vector
02763 {
02764   param_init_bounded_number* v;
02765   int index_min;
02766   int index_max;
02767   double_index_type* it;
02768 public:
02769   param_init_bounded_number_vector();
02770   ~param_init_bounded_number_vector();
02771 
02772   void set_scalefactor(double s);
02773   void set_scalefactor(const dvector& s);
02774   dvector get_scalefactor(void);
02775 
02776 #if defined(OPT_LIB)
02777   param_init_bounded_number& operator [] (int i) { return v[i];}
02778   param_init_bounded_number& operator () (int i) { return v[i];}
02779 #else
02780   param_init_bounded_number& operator [] (int i);
02781   param_init_bounded_number& operator () (int i);
02782 #endif
02783 
02784   void allocate(int min1,int max1,const double_index_type & bmin,
02785     const double_index_type & bmax,const index_type& phase_start,
02786     const char * s);
02787 
02788   void allocate(int min1,int max1,const double_index_type & bmin,
02789      const double_index_type & bmax,const char * s);
02790 
02791   // Added by Steve Martell, Jan 18, 2014.
02792   void allocate(const data_matrix &m, const char *s);
02793 
02794   bool allocated() const { return v != NULL; }
02795   int indexmin() const { return (index_min); }
02796   int indexmax() const { return (index_max); }
02797   void set_initial_value(const double_index_type& it);
02798   void deallocate(void);
02799 };
02800   extern int traceflag;
02801   void tracing_message(int traceflag,const char *s,int *pn);
02802   void tracing_message(int traceflag,const char *s,double *pn);
02803   void set_gauss_covariance_matrix(const dll_data_matrix& m);
02804   void set_gauss_covariance_matrix(const dmatrix& m);
02805   void set_covariance_matrix(const dll_data_matrix& m);
02806   void set_covariance_matrix(const dmatrix& m);
02807 
02808   //ostream& operator<<(const ostream&, const param_init_number_vector);
02809   //ostream& operator<<(const ostream&, const param_init_bounded_number_vector);
02810   //ostream& operator<<(const ostream&, const param_init_vector_vector);
02811   //ostream& operator<<(const ostream&, const param_init_bounded_vector_vector);
02812   //ostream& operator<<(const ostream&, const param_init_matrix_vector);
02813   //ostream& operator<<(const ostream&, const param_init_bounded_matrix_vector);
02814 
02819   class vector_kludge : public dvar_vector
02820   {
02821     public:
02822      vector_kludge(const param_init_number_vector &);
02823      vector_kludge(const param_init_bounded_number_vector &);
02824   };
02829   class matrix_kludge : public dvar_matrix
02830   {
02831     public:
02832      matrix_kludge(const param_init_vector_vector &);
02833      matrix_kludge(const param_init_bounded_vector_vector &);
02834   };
02835 
02836 void read_covariance_matrix(const dmatrix& S, int nvar, int& hbf,
02837   dvector& sscale);
02838 
02839 dvector value(const param_init_number_vector& t);
02840 dvector value(const param_init_bounded_number_vector& t);
02841 //dvector value(const param_init_bounded_number_matrix& t);
02842 //dvector value(const param_init_vector_vector& t);
02843 //dvector value(const param_init_bounded_vector_vector& t);
02844 
02845 dvector read_old_scale(int & old_nvar);
02846 
02847 int withinbound(int lb,int n,int ub);
02848 
02849 #if defined(_MSC_VER)
02850 #  if defined(min)
02851 #    undef min
02852 #  endif
02853 #  if defined(max)
02854 #    undef max
02855 #  endif
02856 #endif
02857 
02858 #include "param_init_bounded_number_matrix.h"
02859 
02860 #endif