ADMB Documentation  11.1.2192
 All Classes Files Functions Variables Typedefs Friends Defines
admodel.h
Go to the documentation of this file.
00001 /*
00002  * $Id: admodel.h 2143 2014-06-25 02:35:44Z johnoel $
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 #if !defined(_MSC_VER)
00069   #include <unistd.h>
00070 #endif
00071 
00072 //#define _ADSTD_ std::
00073 #define param_matrix named_dvar_matrix
00074 #define param_vector named_dvar_vector
00075 #define param_number named_dvariable
00076 #define param_3array named_dvar3_array
00077 #define param_4array named_dvar4_array
00078 #define param_5array named_dvar5_array
00079 #define param_6array named_dvar6_array
00080 #define param_7array named_dvar7_array
00081 
00082 #define SPparam_matrix SPnamed_dvar_matrix
00083 #define SPparam_vector SPnamed_dvar_vector
00084 #define SPparam_number SPnamed_dvariable
00085 #define SPparam_3array SPnamed_dvar3_array
00086 #define SPparam_4array SPnamed_dvar4_array
00087   double mfexp(const double);
00088   dvariable mfexp(const prevariable& v1);
00089   dvar_vector mfexp(const dvar_vector&);
00090   dvector mfexp(const dvector&);
00091 
00092   class param_init_bounded_number_vector;
00093   class model_parameters;
00094   extern int AD_gaussflag;
00095 
00096   extern function_minimizer * pfm;
00097   extern int traceflag;
00098   extern int ADqd_flag;
00099 /*
00100  void set_value_inv(const dvariable& x, const dvector& v, const int& ii);
00101  void set_value_inv(const dvar_matrix& x, const dvector& v, const int& ii);
00102  void set_value_inv(const dvar_vector&, const dvector&, const int &ii);
00103  void set_value_inv(const dvariable& u, const dvector& x, const int& ii,
00104    const double fmin, const double fmax);
00105  void set_value_inv(const dvar_matrix& u, const dvector& x, const int& ii,
00106   const double fmin,const double fmax);
00107  void set_value_inv(const dvar3_array& u, const dvector& x, const int& ii,
00108   const double fmin,const double fmax);
00109  void set_value_inv(const dvar3_array& u, const dvector& x, const int& ii);
00110 */
00111 
00112   void copy_value_to_vector(const prevariable& x, const dvector& v,
00113     const int& ii);
00114   void copy_value_to_vector(const dvar_vector& x, const dvector& v,
00115     const int& ii);
00116   void copy_value_to_vector(const dvar_matrix& x, const dvector& v,
00117     const int& ii);
00118   void copy_value_to_vector(const dvar3_array& x, const dvector& v,
00119     const int& ii);
00120 
00121   void restore_value_from_vector(const prevariable& x, const dvector& v,
00122     const int& ii);
00123   void restore_value_from_vector(const dvar_vector& x, const dvector& v,
00124     const int& ii);
00125   void restore_value_from_vector(const dvar_matrix& x, const dvector& v,
00126     const int& ii);
00127   void restore_value_from_vector(dvar3_array& x, const dvector& v,
00128     const int& ii);
00129 
00134 class AD_matherror
00135 {
00136 public:
00137 #if defined (_BORLANDC_)
00138   exception * err;
00139   AD_matherror(exception * _err) : err(_err) {;}
00140 #endif
00141 #if defined (_MSC_VER)
00142   _exception * err;
00143   AD_matherror(_exception * _err) : err(_err) {;}
00144 #endif
00145 };
00146 
00147 class model_data;
00148 
00153 class label_class
00154 {
00155   const char * name;
00156   friend ostream& operator<<(const ostream& s, const label_class& lc);
00157   friend class model_name_tag;
00158 public:
00159   const char * mychar(void) { return name;}
00160   label_class(const char * s){name=s;}
00161 };
00162 
00167 class model_name_tag
00168 {
00169 protected:
00170   adstring name;
00171   //friend ostream& operator<<(const ostream& os, const model_name_tag& mnt);
00172 public:
00173   model_name_tag(void){;}
00174   void allocate(const char * s);
00175   label_class label(void){return (char*)(name);}
00176   const char * get_name(void) { return name;}
00177 };
00178 
00183 class named_dvar_vector : public dvar_vector, public model_name_tag
00184 {
00185 protected:
00186   named_dvar_vector(void) : dvar_vector(), model_name_tag() {;}
00187   named_dvar_vector& operator=(const dvar_vector& m);
00188   named_dvar_vector& operator=(const dvector& m);
00189   named_dvar_vector& operator=(const double m);
00190   named_dvar_vector& operator=(const prevariable& m);
00191   friend class model_parameters;
00192   void allocate(int mmin,int mmax,const char * s);
00193   void allocate(const char * s);
00194 };
00195 
00200 class equality_constraint_vector : public named_dvar_vector
00201 {
00202 protected:
00203   equality_constraint_vector(void);
00204   equality_constraint_vector& operator=(const dvar_vector& m);
00205   equality_constraint_vector& operator=(const dvector& m);
00206   equality_constraint_vector& operator=(const double m);
00207   equality_constraint_vector& operator=(const prevariable& m);
00208   friend class model_parameters;
00209   void allocate(int mmin,int mmax,const char * s);
00210 };
00211 
00216 class inequality_constraint_vector : public named_dvar_vector
00217 {
00218 protected:
00219   inequality_constraint_vector(void);
00220   inequality_constraint_vector& operator=(const dvar_vector& m);
00221   inequality_constraint_vector& operator=(const dvector& m);
00222   inequality_constraint_vector& operator=(const double m);
00223   inequality_constraint_vector& operator=(const prevariable& m);
00224   friend class model_parameters;
00225   void allocate(int mmin,int mmax,const char * s);
00226 };
00227 
00232 class dll_param_vector : public named_dvar_vector
00233 {
00234   double * pd;
00235 public:
00236   ~dll_param_vector();
00237   void allocate(double *_pd,int mmin,int mmax,const char * s);
00238   dll_param_vector& operator=(const dvar_vector& m);
00239   dll_param_vector& operator=(const dvector& m);
00240   dll_param_vector& operator=(const double m);
00241   dll_param_vector& operator=(const prevariable& m);
00242 };
00243 
00248 class named_dvariable : public dvariable, public model_name_tag
00249 {
00250   //named_dvariable& operator=(const dvariable& m);
00251 protected:
00252   named_dvariable(void) : dvariable(), model_name_tag() {;}
00253   void allocate(const char * s);
00254   named_dvariable& operator=(const prevariable& m);
00255   named_dvariable& operator=(const double m);
00256   friend class model_parameters;
00257 };
00258 
00263 class dll_param_number : public named_dvariable
00264 {
00265   double * pd;
00266 protected:
00267   //named_dvariable(void) : dvariable(), model_name_tag() {;}
00268   void allocate(double *_d,const char * s);
00269   dll_param_number& operator=(const prevariable& m);
00270   dll_param_number& operator=(const double m);
00271   virtual ~dll_param_number();
00272   friend class model_parameters;
00273 };
00274 
00279 class named_dvar_matrix : public dvar_matrix, public model_name_tag
00280 {
00281 protected:
00282   named_dvar_matrix(void) : dvar_matrix(), model_name_tag() {;}
00283   void allocate(int rmin,int rmax,int cmin,int cmax,const char * s);
00284   void allocate(int rmin,int rmax,const char * s);
00285   void allocate(const char * s);
00286   //void allocate(int rmin,int rmax,int,const ivector&,
00287     //const char * s);
00288   void allocate(int rmin,int rmax,const index_type&,const index_type&,
00289     const char * s);
00290 public:
00291   named_dvar_matrix& operator=(const double m);
00292   named_dvar_matrix& operator=(const dmatrix& m);
00293   named_dvar_matrix& operator=(const dvar_matrix& m);
00294   named_dvar_matrix& operator=(const dvariable& m);
00295   friend class model_parameters;
00296 };
00297 
00302 class dll_param_matrix : public named_dvar_matrix
00303 {
00304   double * pd;
00305 public:
00306   //named_dvar_matrix(void) : dvar_matrix(), model_name_tag() {;}
00307   void allocate(double * ,int rmin,int rmax,int cmin,int cmax,const char * s);
00308   void allocate(double * ,int rmin,int rmax,const index_type&
00309      ,const index_type& ,const char * s);
00310   dll_param_matrix& operator=(const dvariable& m);
00311   dll_param_matrix & operator = (const double m);
00312   dll_param_matrix& operator=(const dmatrix& m);
00313   dll_param_matrix& operator=(const dvar_matrix& m);
00314   virtual ~dll_param_matrix();
00315 };
00316 
00321 class named_dvar3_array : public dvar3_array, public model_name_tag
00322 {
00323 protected:
00324   named_dvar3_array(void) : dvar3_array(), model_name_tag() {;}
00325  void allocate(const ad_integer& sl,const ad_integer& sh,
00326   const index_type& nrl,const index_type& nrh,const index_type& ncl,
00327     const index_type& nch,const char * s="UNNAMED");
00328   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
00329     const char * s="UNNAMED");
00330 
00331   void allocate(int hsl,int hsu,int rmin,int rmax,const char * s="UNNAMED");
00332   void allocate(int hsl,int hsu,const index_type& rmin,const index_type& rmax,
00333     const char * s="UNNAMED");
00334   void allocate(int hsl,int hsu,const char * s="UNNAMED");
00335   void allocate(const char * s="UNNAMED");
00336 
00337 
00338   void allocate(int hsl, int hsu, const ivector& rmin, int rmax,
00339                 int cmin, int cmax, const char *s = "UNNAMED");
00340   void allocate(int hsl, int hsu, int rmin, const ivector& rmax,
00341     int cmin,int cmax,const char * s="UNNAMED");
00342   void allocate(int hsl, int hsu, const ivector& rmin, const ivector& rmax,
00343     int cmin,int cmax,const char * s="UNNAMED");
00344   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
00345     int cmax,const char * s="UNNAMED");
00346   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
00347                 const ivector& cmax, const char *s = "UNNAMED");
00348   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,
00349                 const ivector& cmax, const char *s = "UNNAMED");
00350   named_dvar3_array& operator=(const dvar3_array& m);
00351   named_dvar3_array& operator=(const d3_array& m);
00352   friend class model_parameters;
00353 };
00354 
00359 class named_dvar4_array : public dvar4_array, public model_name_tag
00360 {
00361 protected:
00362   named_dvar4_array(void) : dvar4_array(), model_name_tag() {;}
00363   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,int rmax,
00364   int cmin,int cmax,const char * s);
00365   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00366     const index_type&,const index_type&,const index_type&,const index_type&,
00367     const char *s);
00368 
00369   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00370     const index_type&,const index_type&,const char *s);
00371 
00372   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00373     const char *s);
00374   void allocate(ad_integer,ad_integer,const char *s);
00375   void allocate(const char *s);
00376 
00377   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,int rmax,
00378     const char * s);
00379   void allocate(int hhsl,int hhsu,int hsl,int hsu,const char * s);
00380   void allocate(int hhsl,int hhsu,const char * s);
00381 
00382   named_dvar4_array& operator=(const dvar4_array& m);
00383   named_dvar4_array& operator=(const d4_array& m);
00384   friend class model_parameters;
00385 };
00386 
00391 class named_dvar5_array : public dvar5_array, public model_name_tag
00392 {
00393 protected:
00394   named_dvar5_array(void) : dvar5_array(), model_name_tag() {;}
00395   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,int rmax,
00396   int cmin,int cmax,
00397   int l5,int u5,
00398   const char * s);
00399 
00400   void allocate(const ad_integer& hhsl,const ad_integer& hhsu,
00401     const index_type& hsl,const index_type& hsu,
00402     const index_type& sl,const index_type& sh,
00403     const index_type& nrl,const index_type& nrh,
00404     const index_type& ncl,const index_type& nch,
00405     const char * s);
00406 
00407   named_dvar5_array& operator=(const dvar5_array& m);
00408   named_dvar5_array& operator=(const d5_array& m);
00409   friend class model_parameters;
00410 };
00411 
00416 class named_dvar6_array : public dvar6_array, public model_name_tag
00417 {
00418 protected:
00419   named_dvar6_array(void) : dvar6_array(), model_name_tag() {;}
00420   void allocate(int hhsl,int hhsu,
00421     int hsl,int hsu,
00422     int rmin,int rmax,
00423     int cmin,int cmax,
00424     int l5,int u5,
00425     int l6,int u6,
00426     const char * s);
00427 
00428   void allocate(const ad_integer& hhsl,const ad_integer& hhsu,
00429     const index_type& hsl,const index_type& hsu,
00430     const index_type& sl,const index_type& sh,
00431     const index_type& nrl,const index_type& nrh,
00432     const index_type& ncl,const index_type& nch,
00433     const index_type& l5,const index_type& u5,
00434     const char * s);
00435 
00436   named_dvar6_array& operator=(const dvar6_array& m);
00437   named_dvar6_array& operator=(const d6_array& m);
00438   friend class model_parameters;
00439 };
00440 
00445 class named_dvar7_array : public dvar7_array, public model_name_tag
00446 {
00447 protected:
00448   named_dvar7_array(void) : dvar7_array(), model_name_tag() {;}
00449   void allocate(
00450     int hhsl,int hhsu,
00451     int hsl,int hsu,
00452     int rmin,int rmax,
00453     int cmin,int cmax,
00454     int l5,int u5,
00455     int l6,int u6,
00456     int l7,int u7,
00457     const char * s);
00458 
00459   void allocate(const ad_integer& hhsl,const ad_integer& hhsu,
00460     const index_type& hsl,const index_type& hsu,
00461     const index_type& sl,const index_type& sh,
00462     const index_type& nrl,const index_type& nrh,
00463     const index_type& ncl,const index_type& nch,
00464     const index_type& l5,const index_type& u5,
00465     const index_type& l6,const index_type& u6,
00466     const char * s);
00467 
00468   named_dvar7_array& operator=(const dvar7_array& m);
00469   named_dvar7_array& operator=(const d7_array& m);
00470   friend class model_parameters;
00471 };
00472 
00477 class named_dvector : public dvector, public model_name_tag
00478 {
00479 protected:
00480   named_dvector(void) : dvector(), model_name_tag() {;}
00481   void allocate(int mmin,int mmax,const char * s);
00482   void allocate(const char * s);
00483   void allocate(int mmin, const ivector& mmax, const char *s);
00484   named_dvector& operator=(const dvector& m);
00485   named_dvector& operator=(const double m);
00486   friend class model_data;
00487 };
00488 
00493 class named_ivector : public ivector, public model_name_tag
00494 {
00495 protected:
00496   named_ivector(void) : ivector(), model_name_tag() {;}
00497   void allocate(int mmin,int mmax,const char * s);
00498 };
00499 
00504 class named_dmatrix : public dmatrix, public model_name_tag
00505 {
00506 protected:
00507   named_dmatrix(void) : dmatrix(), model_name_tag() {;}
00508   void allocate(int rmin,int rmax,const char * s);
00509   void allocate(const char * s);
00510   void allocate(int rmin,int rmax,int cmin,int cmax,const char * s);
00511   void allocate(int rmin,int rmax,const ivector& cmin,int cmax,const char * s);
00512   void allocate(int rmin,int rmax,const ivector& cmin,const ivector& cmax,
00513     const char * s);
00514   void allocate(int rmin,int rmax,int cmin,const ivector& cmax,const char * s);
00515   named_dmatrix& operator=(const dmatrix& m);
00516   named_dmatrix& operator=(const double m);
00517 };
00518 
00523 class named_imatrix : public imatrix, public model_name_tag
00524 {
00525 protected:
00526   named_imatrix(void) : imatrix(), model_name_tag() {;}
00527   void allocate(int rmin,int rmax,int cmin,int cmax,const char * s);
00528   void allocate(int rmin, int rmax, const index_type& cmin,
00529                 const index_type& cmax, const char *s);
00530   named_imatrix& operator=(const imatrix& m);
00531   named_imatrix& operator=(const int& m);
00532 };
00533 
00538 class named_d3_array : public d3_array, public model_name_tag
00539 {
00540 protected:
00541   named_d3_array(void) : d3_array(), model_name_tag() {;}
00542   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
00543     const char * s);
00544   void allocate(int hsl, int hsu, const index_type& rmin,
00545                 const index_type& rmax, const index_type& cmin,
00546                 const index_type& cmax, const char *s);
00547   void allocate(int hsl, int hsu, const ivector& rmin, int rmax,
00548     int cmin,int cmax,const char * s);
00549   void allocate(int hsl, int hsu, int rmin, const ivector& rmax,
00550     int cmin,int cmax,const char * s);
00551   void allocate(int hsl, int hsu, const ivector& rmin, const ivector& rmax,
00552     int cmin,int cmax,const char * s);
00553   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
00554     int cmax, const char * s);
00555   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,
00556                 const ivector& cmax, const char *s);
00557   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
00558                 const ivector& cmax, const char *s);
00559   named_d3_array& operator=(const d3_array& m);
00560 
00561   void allocate(int hsl,int hsu,const index_type& rmin,const index_type& rmax,
00562     const char * s);
00563   void allocate(int hsl,int hsu,int rmin,int rmax,const char * s);
00564   void allocate(int hsl,int hsu,const char * s);
00565   void allocate(const char * s);
00566 };
00567 
00572 class named_i3_array : public i3_array, public model_name_tag
00573 {
00574 protected:
00575   named_i3_array(void) : i3_array(), model_name_tag() {;}
00576   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
00577     const char * s);
00578   void allocate(int hsl, int hsu, const index_type& rmin,
00579                 const index_type& rmax, const index_type& cmin,
00580                 const index_type& cmax, const char *s);
00581   named_i3_array& operator=(const i3_array& m);
00582 };
00583 
00588 class named_d4_array : public d4_array, public model_name_tag
00589 {
00590 protected:
00591   named_d4_array(void) : d4_array(), model_name_tag() {;}
00592   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00593     const index_type&,const index_type&,const index_type&,const index_type&,
00594     const char * s);
00595   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,
00596     int rmax,int cmin,int cmax,const char * s);
00597 
00598   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,
00599     int rmax,const char * s);
00600   void allocate(int hhsl,int hhsu,int hsl,int hsu,const char * s);
00601   void allocate(int hhsl,int hhsu,const char * s);
00602   void allocate(const char * s);
00603 
00604   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00605     const index_type&,const index_type&,const char * s);
00606   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00607     const char * s);
00608   void allocate(ad_integer,ad_integer,const char * s);
00609 
00610 
00611   named_d4_array& operator=(const d4_array& m);
00612 };
00613 
00618 class named_d5_array : public d5_array, public model_name_tag
00619 {
00620 protected:
00621   named_d5_array(void) : d5_array(), model_name_tag() {;}
00622   void allocate(int l5,int u5,int hhsl,int hhsu,int hsl,int hsu,int rmin,
00623     int rmax,int cmin,int cmax,const char * s);
00624   void allocate(const ad_integer& hhsl,const ad_integer& hhsu,
00625     const index_type& hsl,const index_type& hsu, const index_type& sl,
00626     const index_type& sh,const index_type& nrl,const index_type& nrh,
00627     const index_type& ncl,const index_type& nch,const char * s);
00628 
00629   named_d5_array& operator=(const d5_array& m);
00630 };
00631 
00636 class named_d6_array : public d6_array, public model_name_tag
00637 {
00638 protected:
00639   named_d6_array(void) : d6_array(), model_name_tag() {;}
00640   void allocate(int l6,int u6,int l5,int u5,int hhsl,int hhsu,int hsl,
00641     int hsu,int rmin,int rmax,int cmin,int cmax,const char * s);
00642   void allocate(const ad_integer& l6,const ad_integer& u6,
00643     const index_type& l5,const index_type& u5,
00644     const index_type& hhsl,const index_type& hhsu,
00645     const index_type& hsl,const index_type& hsu,
00646     const index_type& sl,const index_type& sh,
00647     const index_type& nrl,const index_type& nrh,
00648     const char * s);
00649 
00650   named_d6_array& operator=(const d6_array& m);
00651 };
00652 
00657 class named_d7_array : public d7_array, public model_name_tag
00658 {
00659 protected:
00660   named_d7_array(void) : d7_array(), model_name_tag() {;}
00661   void allocate(int l7,int u7,int l6,int u6,int l5,int u5,int hhsl,
00662     int hhsu,int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
00663     const char * s);
00664   void allocate(const ad_integer& l7,const ad_integer& u7,
00665     const index_type& l6,const index_type& u6,
00666     const index_type& l5,const index_type& u5,
00667     const index_type& hhsl,const index_type& hhsu,
00668     const index_type& hsl,const index_type& hsu,
00669     const index_type& sl,const index_type& sh,
00670     const index_type& nrl,const index_type& nrh,
00671     const char * s);
00672 
00673   named_d7_array& operator=(const d7_array& m);
00674 };
00675 
00676 
00677 class function_minimizer;
00678 
00679 #if defined(USE_ADPVM)
00680 
00684 class pvm_params
00685 {
00686   static pvm_params * varsptr[]; // this should be a resizeable array
00687   static int num_pvm_params;
00688   static const int maxnum_pvm_params;
00689   void add_to_list(void);
00690   virtual void send_to_slaves(void)=0;
00691   virtual void get_from_master(void)=0;
00692 public:
00693   static void pvm_params::send_all_to_slaves(void);
00694   static void pvm_params::get_all_from_master(void);
00695   void allocate(const char *);
00696   void allocate(void);
00697 };
00698 
00703 class pvm_number : public pvm_params
00704 {
00705 public:
00706   virtual void send_to_slaves(void);
00707   virtual void get_from_master(void);
00708   dvector v;
00709   double d;
00710   operator double();
00711   void assign(const dvector&);
00712   void assign(double);
00713 };
00714 
00719 class pvm_int : public pvm_params
00720 {
00721 public:
00722   virtual void send_to_slaves(void);
00723   virtual void get_from_master(void);
00724   ivector v;
00725   int d;
00726   operator int();
00727   void assign(const ivector&);
00728   void assign(int);
00729 };
00730 #endif // #if defined(USE_ADPVM)
00731 
00732 class initial_params;
00733 typedef initial_params* pinitial_params;
00734 typedef void* ptovoid;
00738 class adlist_ptr
00739 {
00740   ptovoid * ptr;
00741   int current_size;
00742   int current;
00743   void resize(void);
00744   void add_to_list(void* p);
00745 public:
00746   adlist_ptr(int init_size);
00747   ~adlist_ptr();
00748 
00749   void initialize();
00750 
00751   pinitial_params& operator[](int i);
00752 
00753   friend class initial_params;
00754 };
00755 
00756 #if defined(USE_SHARE_FLAGS)
00757 
00762   class shareinfo
00763   {
00764     index_type * sflags;
00765     index_type * original_sflags;
00766     index_type * aflags;
00767     index_type * invflags;
00768     i3_array * bmap;
00769     int dimension;
00770     int maxshare;
00771     int current_phase;
00772   public:
00773     void get_inv_matrix_shared( int _cf);
00774     void get_inv_vector_shared( int _cf);
00775     int & get_maxshare(void) { return maxshare; }
00776     i3_array & get_bmap(void);
00777     int & get_dimension(void){ return dimension;}
00778     index_type * get_original_shareflags(void);
00779     index_type * get_shareflags(void);
00780     int& get_current_phase(void);
00781     index_type * get_activeflags(void);
00782     index_type * get_invflags(void);
00783     void set_shareflags(const index_type& sf);
00784     void set_original_shareflags(const index_type& sf);
00785     void reset_shareflags(const index_type& sf);
00786     void set_activeflags(const index_type& af);
00787     void set_bmap(const i3_array& af);
00788     void reset_bmap(const i3_array& af);
00789     void set_invflags(const index_type& af);
00790     shareinfo(const index_type& sf,const index_type& af);
00791     ~shareinfo();
00792   };
00793 #endif
00794 
00799 class initial_params
00800 {
00801 protected:
00802 #if defined(USE_SHARE_FLAGS)
00803   shareinfo * share_flags;
00804 #endif
00805   virtual ~initial_params();
00806   int active_flag;
00807   int initial_value_flag;
00808   double initial_value;
00809   double scalefactor;
00810 public:
00811 #if defined(USE_SHARE_FLAGS)
00812   virtual void setshare(const index_type& sf,const index_type& af);
00813   virtual void shared_set_value_inv(const dvector&,const int&);
00814   virtual void shared_set_value(const dvar_vector&,const int&,
00815     const dvariable& pen);
00816   virtual int shared_size_count(void); // get the number of active parameters
00817 #endif
00818   double get_scalefactor();
00819   void set_scalefactor(const double);
00820 #if !defined(BIG_INIT_PARAMS)
00821   static initial_params  varsptr[]; // this should be a resizeable array
00822 #else
00823   static adlist_ptr varsptr; // this should be a resizeable array
00824 #endif
00825   static int num_initial_params;
00826   static const int max_num_initial_params;
00827   static int straight_through_flag;
00828   static int num_active_initial_params;
00829   static int max_number_phases;
00830   static int current_phase;
00831   static int restart_phase;
00832   static int sd_phase;
00833   static int mc_phase;
00834   static int mceval_phase;
00835   int phase_start;
00836   int phase_save;
00837   int phase_stop;
00838   virtual void set_random_effects_active();
00839   void restore_phase_start(void);
00840   virtual void set_random_effects_inactive();
00841   virtual void set_only_random_effects_active();
00842   virtual void set_only_random_effects_inactive();
00843   virtual void set_value(const dvar_vector&, const int&,
00844     const dvariable& pen) = 0;
00845   virtual void dev_correction(const dmatrix&, const int&) = 0;
00846   void set_initial_value(double x);
00847   double get_initial_value(void);
00848   void set_phase_start(int x);
00849   int get_phase_start(void);
00850   static void set_all_simulation_bounds(const dmatrix& symbds);
00851   static void set_all_simulation_bounds(const dmatrix& symbds, const dvector&);
00852   static void get_jacobian_value(const dvector& y, const dvector& jac);
00853   static int correct_for_dev_objects(const dmatrix &H);
00854   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii) = 0;
00855   //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
00856   //  const int& ii) = 0;
00857   virtual void set_value_inv(const dvector&, const int&) = 0;
00858   virtual void add_value(const dvector&, const int&) = 0;
00859   virtual void add_value(const dvector&, const dvector&, const int&,
00860     const double&, const dvector&) = 0;
00861   virtual void get_jacobian(const dvector&, const dvector&, const int&) = 0;
00862   //virtual void check_tightness(const dvector&, const int&) = 0;
00863   virtual void copy_value_to_vector(const dvector&, const int&) = 0;
00864   virtual void restore_value_from_vector(const dvector&, const int&) = 0;
00865   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii) = 0;
00866   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
00867     const int& ii)=0;
00868   virtual void mc_scale(const dvector& d, const dvector& x, const int& ii) = 0;
00869   virtual void curv_scale(const dvector& d, const dvector& x,
00870     const int& ii) = 0;
00871   virtual void hess_scale(const dvector& d, const dvector& x,
00872     const int& ii) = 0;
00873   virtual int size_count(void)=0; // get the number of active parameters
00874   virtual void save_value(void)=0; // save the objects value in a ascii file
00875   virtual void bsave_value(void)=0; // save the objects value in a binary file
00876   virtual void save_value(const ofstream& ofs, int prec) = 0;
00877   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
00878     int& offset){}
00879   //virtual void bsave_value(const uostream& ofs) = 0;
00880     virtual const char * label()=0;
00881   void allocate(int _phase_start);
00882   void set_active_flag(void);
00883   void set_inactive_flag(void);
00884   friend int active(const initial_params& ip);
00885   static adstring get_reportfile_name(void);
00886   initial_params(void);
00887   static void xinit(const dvector& x); // get the number of active parameters
00888   // get the number of active parameters
00889   static void xinit_all(const dvector& x);
00890   static void save_all(const ofstream& _ofs,int prec,const dvector&g);
00891   // get the number of active parameters
00892   static void set_active_random_effects(void);
00893   // get the number of active parameters
00894   static void set_active_only_random_effects(void);
00895   // get the number of active parameters
00896   static void set_inactive_only_random_effects(void);
00897   // get the number of active parameters
00898   static void set_inactive_random_effects(void);
00899   // get the number of active parameters
00900   static void restore_start_phase(void);
00901   static void xinit1(const dvector& x, const dvector& g);
00902   //save all initial parameter values in a vector
00903   static void copy_all_values(const dvector& x, const int& ii);
00904   //get ivalues for all active parameters from a vector
00905   static void restore_all_values(const dvector& x, const int& ii);
00906   // get the number of active parameters
00907   static dvariable reset(const dvar_vector& x);
00908   // get the number of active parameters
00909   static dvariable reset(const dvector& x);
00910   static dvariable reset1(const dvar_vector& x, const dvector& g);
00911   // get the number of active parameters
00912   static dvariable reset(const dvar_vector& x, const dvector& pen);
00913   // get the number of active parameters
00914   static dvariable reset_all(const dvar_vector& x,const dvector& pen);
00915   static int nvarcalc(void);
00916   static int nvarcalc_all(void);
00917   static int num_active_calc(void);
00918   static int stddev_scale(const dvector& d, const dvector& x);
00919   static int stddev_vscale(const dvar_vector& d,const dvar_vector& x);
00920   static int montecarlo_scale(const dvector& d, const dvector& x);
00921   static int stddev_curvscale(const dvector& d, const dvector& x);
00922   static void read(void);
00923   static void save(void);
00924   static void save(const ofstream& ofs, int prec);
00925   static void restore(const ifstream& ifs);
00926   static void add_random_vector(const dvector& x);
00927   static void add_random_vector(const dvector& y, const dvector& x,
00928     const double& ll, const dvector& diag);
00929   virtual void restore_value(const ifstream& ifs) = 0;
00930   virtual void add_to_list(void);
00931 #if defined(USE_ADPVM)
00932   virtual void pvm_pack(void)=0;
00933   virtual void pvm_unpack(void)=0;
00934 #endif
00935 
00936   friend class function_minimizer;
00937 };
00938 
00939 void pvm_pack(const dvar_vector&);
00940 void pvm_unpack(const dvar_vector&);
00941 void pvm_pack(const prevariable&);
00942 void pvm_unpack(const prevariable&);
00943 void pvm_pack(const dvar_matrix&);
00944 void pvm_unpack(const dvar_matrix&);
00945 
00950 class param_init_vector: public named_dvar_vector , public initial_params
00951 {
00952 public:
00953 #if defined(USE_SHARE_FLAGS)
00954   virtual void setshare(const index_type& sf,const index_type& af);
00955   virtual void shared_set_value_inv(const dvector&,const int&);
00956   virtual void shared_set_value(const dvar_vector&,const int&,
00957     const dvariable& pen);
00958   virtual int shared_size_count(void); // get the number of active parameters
00959 #endif
00960   virtual void dev_correction(const dmatrix& H, const int& ii);
00961   virtual const char * label(void);
00962   param_init_vector();
00963 #if defined(USE_ADPVM)
00964   void pvm_pack(void){::pvm_pack(*this);}
00965   void pvm_unpack(void){::pvm_unpack(*this);}
00966 #endif
00967 private:
00968   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
00969   //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
00970   //  const int& ii);
00971   virtual void add_value(const dvector&, const dvector&, const int&,
00972     const double&, const dvector&);
00973   virtual void add_value(const dvector&, const int&);
00974   virtual void get_jacobian(const dvector&, const dvector&, const int&);
00975   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
00976     const int& ii);
00977 
00978   virtual void set_value(const dvar_vector& x, const int& ii,
00979     const dvariable& pen);
00980   virtual void set_value_inv(const dvector& x, const int& ii);
00981   virtual void copy_value_to_vector(const dvector& x, const int& ii);
00982   virtual void restore_value_from_vector(const dvector&, const int&);
00983   virtual int size_count(void);
00984   virtual void sd_scale(const dvector& d, const dvector& x, const int&);
00985   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
00986   virtual void curv_scale(const dvector& d, const dvector& x,const int&);
00987   virtual void hess_scale(const dvector&, const dvector&, const int&){}
00988   virtual void save_value(void);
00989   virtual void bsave_value(void);
00990   virtual void save_value(const ofstream& ofs, int prec);
00991   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
00992     int& offset);
00993   virtual void restore_value(const ifstream& ifs);
00994   void report_value(void);
00995   //virtual void read_value(void);
00996   void allocate(int imin,int imax,int phasestart=1,const char * s="UNNAMED");
00997   void allocate(const ad_integer& imin,const ad_integer& imax,
00998     const ad_integer& phasestart=1,const char * s="UNNAMED");
00999   void allocate(int imin,int imax,const char * s="UNNAMED");
01000 #if defined(USE_SHARE_FLAGS)
01001   void allocate(int imin,int imax,const ivector& ishare,
01002     const char * s="UNNAMED");
01003 #endif
01004   friend class model_parameters;
01005   friend class param_init_vector_vector;
01006 public:
01007   param_init_vector& operator = (const dvector&);
01008   param_init_vector& operator = (const dvar_vector&);
01009   param_init_vector& operator = (const prevariable&);
01010   param_init_vector& operator = (const double&);
01011 };
01012 
01017 class dll_param_init_vector: public param_init_vector
01018 {
01019   double * pd;
01020 public:
01021   dll_param_init_vector& operator = (const dvector&);
01022   dll_param_init_vector& operator = (const dvar_vector&);
01023   dll_param_init_vector& operator = (const prevariable&);
01024   dll_param_init_vector& operator = (const double&);
01025   void allocate(double * _pd,int imin,int imax,
01026     int phasestart=1,const char * s="UNNAMED");
01027   void allocate(double * _pd,int imin,int imax,
01028     const char * s="UNNAMED");
01029 
01030   virtual ~dll_param_init_vector();
01031   friend class model_parameters;
01032 };
01033 
01038 class param_init_bounded_vector: public named_dvar_vector,public initial_params
01039 {
01040   virtual void* parent_this(void){return this;}
01041 public:
01042   double get_minb(void);
01043   void set_minb(double b);
01044   double get_maxb(void);
01045   void set_maxb(double b);
01046 protected:
01047   double minb;
01048   double maxb;
01049   param_init_bounded_vector();
01050 private:
01051   virtual void dev_correction(const dmatrix&, const int&);
01052   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01053   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01054     const int& ii);
01055   //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01056   //  const int& ii);
01057   virtual void add_value(const dvector&, const dvector&, const int&,
01058     const double&, const dvector&);
01059   virtual void add_value(const dvector&, const int&);
01060   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01061   virtual void set_value(const dvar_vector& x, const int& ii,
01062     const dvariable& pen);
01063   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01064   virtual void restore_value_from_vector(const dvector&, const int&);
01065   virtual void set_value_inv(const dvector& x, const int& ii);
01066   virtual int size_count(void);
01067   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01068   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01069   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01070   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01071   void allocate(int imin,int imax,double _minb,double _maxb,
01072     int phasestart=1, const char * name="UNNAMED");
01073   void allocate(int imin,int imax,double _minb,double _maxb,
01074     const char * name="UNNAMED");
01075  //void param_init_bounded_vector::allocate(const ad_integer& imin,
01076    //const ad_integer& imax,const ad_double& _minb,const ad_double& _maxb,
01077    //const ad_integer& phase_start,const char * s);
01078   friend class model_parameters;
01079   friend class param_init_bounded_vector_vector;
01080   virtual const char * label(void);
01081   virtual void save_value(const ofstream& ofs, int prec);
01082   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
01083     int& offset);
01084   virtual void restore_value(const ifstream& ifs);
01085   virtual void save_value(void);
01086   virtual void bsave_value(void);
01087   void report_value(void);
01088   //virtual void read_value(void);
01089 public:
01090   param_init_bounded_vector& operator = (const dvector&);
01091   param_init_bounded_vector& operator = (const dvar_vector&);
01092   param_init_bounded_vector& operator = (const prevariable&);
01093   param_init_bounded_vector& operator = (const double&);
01094 #if defined(USE_ADPVM)
01095   void pvm_pack(void){::pvm_pack(*this);}
01096   void pvm_unpack(void){::pvm_unpack(*this);}
01097 #endif
01098 };
01099 
01104 class dll_param_init_bounded_vector: public param_init_bounded_vector
01105 {
01106   double * pd;
01107 public:
01108   void allocate(double * _pd,int imin,int imax,double _minb,double _maxb,
01109     int phasestart=1, const char * name="UNNAMED");
01110   void allocate(double * _pd,int imin,int imax,double _minb,double _maxb,
01111     const char * name="UNNAMED");
01112   ~dll_param_init_bounded_vector();
01113 };
01114 
01119 class param_init_bounded_dev_vector: public param_init_bounded_vector
01120 {
01121   virtual void set_value(const dvar_vector& x, const int& ii,
01122     const dvariable& pen);
01123   virtual void dev_correction(const dmatrix& H, const int& ii);
01124 public:
01125   param_init_bounded_dev_vector& operator = (const dvar_vector& m);
01126   param_init_bounded_dev_vector& operator = (const dvector& m);
01127   param_init_bounded_dev_vector& operator = (const prevariable& m);
01128   param_init_bounded_dev_vector& operator = (const double& m);
01129 };
01130 
01135 class param_init_number: public named_dvariable , public initial_params
01136 {
01137   virtual void dev_correction(const dmatrix&, const int&);
01138   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01139 
01140 #if defined(USE_ADPVM)
01141   void pvm_pack(void){::pvm_pack(*this);}
01142   void pvm_unpack(void){::pvm_unpack(*this);}
01143 #endif
01144 
01145 //  virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01146 //    const int& ii);
01147   virtual void add_value(const dvector&, const dvector&, const int&,
01148     const double&, const dvector&);
01149   virtual void add_value(const dvector&, const int&);
01150   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01151   virtual void set_value(const dvar_vector& x, const int& ii,
01152     const dvariable& pen);
01153   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01154   virtual void restore_value_from_vector(const dvector&, const int&);
01155   virtual void set_value_inv(const dvector& x, const int& ii);
01156   virtual int size_count(void);
01157   virtual void save_value(const ofstream& ofs, int prec);
01158   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
01159     int& offset);
01160   virtual void restore_value(const ifstream& ifs);
01161   virtual void save_value(void);
01162   virtual void bsave_value(void);
01163   void report_value(void);
01164   virtual const char * label(void);
01165   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01166   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01167   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01168   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01169   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01170     const int& ii);
01171   //virtual void read_value(void);
01172 protected:
01173   void allocate(int phase_start=1,const char *s="UNNAMED");
01174   void allocate(const char *s="UNNAMED");
01175   void allocate(init_xml_doc&, const char *s="UNNAMED");
01176   friend class model_parameters;
01177   friend class param_init_number_vector;
01178   param_init_number();
01179   param_init_number& operator=(const double m);
01180   param_init_number& operator=(const prevariable& m);
01181 };
01182 
01187 class dll_param_init_number: public param_init_number
01188 {
01189   double * pd;
01190 public:
01191   void allocate(double * pd,int phase_start=1,const char *s="UNNAMED");
01192   void allocate(double *pd,const char *s="UNNAMED");
01193   virtual ~dll_param_init_number();
01194   dll_param_init_number& operator=(const double m);
01195   dll_param_init_number& operator=(const prevariable& m);
01196 };
01197 
01198 //forward declaration
01199 class data_vector;
01206 class param_init_bounded_number: public param_init_number
01207 {
01208 public:
01209   double get_minb(void);
01210   void set_minb(double b);
01211   double get_maxb(void);
01212   void set_maxb(double b);
01213 protected:
01214   double minb;
01215   double maxb;
01216   void allocate(double _minb,double _maxb,int phase_start=1,
01217     const char * s="UNNAMED");
01218   void allocate(double _minb,double _maxb,const char * s="UNNAMED");
01219   // Added by Steve Martell for using input data for allocation.
01220   void allocate(const data_vector& v,const char * s="UNNAMED");
01221   void allocate(init_xml_doc&, const char * s="UNNAMED");
01222 
01223 public:
01224 #if defined(USE_ADPVM)
01225   void pvm_pack(void){::pvm_pack(*this);}
01226   void pvm_unpack(void){::pvm_unpack(*this);}
01227 #endif
01228   param_init_bounded_number();
01229 private:
01230   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01231 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01232 //  const int& ii);
01233   virtual void add_value(const dvector&, const dvector&, const int&,
01234     const double&, const dvector&);
01235   virtual void add_value(const dvector&, const int&);
01236   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01237   virtual void set_value(const dvar_vector& x, const int& ii,
01238     const dvariable& pen);
01239   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01240   virtual void restore_value_from_vector(const dvector&, const int&);
01241   virtual void set_value_inv(const dvector& x, const int& ii);
01242   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01243   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01244   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01245   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01246   virtual const char * label(void);
01247   void report_value(void);
01248   param_init_bounded_number& operator=(const double m);
01249   param_init_bounded_number& operator=(const prevariable& m);
01250   friend class model_parameters;
01251   friend class param_init_bounded_number_vector;
01252   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01253     const int& ii);
01254 };
01255 
01260 class dll_param_init_bounded_number: public param_init_bounded_number
01261 {
01262   double * pd;
01263 public:
01264   void allocate(double * _pd,double _minb,double _maxb,int phase_start=1,
01265     const char * s="UNNAMED");
01266   void allocate(double* _pd,double _minb, double _maxb,const char* s="UNNAMED");
01267 public:
01268   virtual ~dll_param_init_bounded_number();
01269   void report_value(void);
01270 };
01271 
01276 class param_init_matrix: public named_dvar_matrix,public initial_params
01277 {
01278 #if defined(USE_SHARE_FLAGS)
01279   virtual void shared_set_value_inv(const dvector& x,const int& ii);
01280   virtual int shared_size_count(void); // get the number of active parameters
01281   virtual void shared_set_value(const dvar_vector&,const int&,
01282     const dvariable& pen);
01283 #endif
01284   virtual void dev_correction(const dmatrix&, const int&);
01285   virtual void set_simulation_bounds(const dmatrix& symbds,const int& ii);
01286 #if defined(USE_ADPVM)
01287   void pvm_pack(void){::pvm_pack(*this);}
01288   void pvm_unpack(void){::pvm_unpack(*this);}
01289 #endif
01290 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01291 //  const int& ii);
01292   virtual void add_value(const dvector&, const dvector&, const int&,
01293     const double&, const dvector&);
01294   virtual void add_value(const dvector&, const int&);
01295   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01296 public:
01297   virtual void set_value(const dvar_vector& x, const int& ii,
01298     const dvariable& pen);
01299   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01300   virtual void restore_value_from_vector(const dvector&, const int&);
01301   virtual void set_value_inv(const dvector& x, const int& ii);
01302   virtual int size_count(void);
01303   virtual void save_value(void);
01304   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
01305     int& offset);
01306   virtual void bsave_value(void);
01307   virtual void save_value(const ofstream& ofs, int prec);
01308   virtual void restore_value(const ifstream& ifs);
01309   void report_value(void);
01310   //virtual void read_value(void);
01311   virtual const char * label(void);
01312   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01313   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01314   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01315   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01316   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01317     const int& ii);
01318 
01319 public:
01320   void allocate(int rmin,int rmax,int cmin,int cmax,
01321     int phase_start=1,const char * = "UNNAMED");
01322 #if defined(USE_SHARE_FLAGS)
01323   //void allocate(int rmin,int rmax,int cmin,int cmax, const imatrix& jshare,
01324   // const char * = "UNNAMED");
01325   virtual void setshare(const index_type& sf,const index_type& af);
01326 #endif
01327   void allocate(int rmin,int rmax,int cmin,int cmax,
01328     const char * = "UNNAMED");
01329  void allocate(const ad_integer& imin,const ad_integer&imax,
01330    const index_type& imin2,const index_type& imax2,
01331    const ad_integer& phase_start, const char * s);
01332   void allocate(const ad_integer& rmin, const ad_integer& rmax,
01333                 const index_type& cmin, const index_type& cmax,
01334                 const char * = "UNNAMED");
01335   void allocate(const ad_integer& rmin, const ad_integer& rmax,
01336                 const index_type& cmin, const index_type& cmax,
01337                 int phase_start = 1, const char * = "UNNAMED");
01338   param_init_matrix(void);
01339   param_init_matrix& operator = (const dmatrix& m);
01340   param_init_matrix& operator = (const dvar_matrix& m);
01341   param_init_matrix& operator = (const dvariable& m);
01342   param_init_matrix& operator = (const double& m);
01343 };
01344 
01349 class dll_param_init_matrix: public param_init_matrix
01350 {
01351   double * d;
01352 public:
01353   void allocate(double* _d,int rmin,int rmax,int cmin,int cmax,
01354     int phase_start=1,const char * = "UNNAMED");
01355   void allocate(double * _d,int rmin,int rmax,int cmin,int cmax,
01356     const char * = "UNNAMED");
01357   virtual ~dll_param_init_matrix();
01358   dll_param_init_matrix(){d=NULL;}
01359   dll_param_init_matrix& operator = (const dmatrix& m);
01360   dll_param_init_matrix& operator = (const dvar_matrix& m);
01361   dll_param_init_matrix& operator = (const dvariable& m);
01362   dll_param_init_matrix& operator = (const double& m);
01363 };
01364 
01369 class param_init_bounded_matrix: public param_init_matrix
01370 {
01371 public:
01372   double get_minb(void);
01373   void set_minb(double b);
01374   double get_maxb(void);
01375   void set_maxb(double b);
01376 protected:
01377   double minb;
01378   double maxb;
01379   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01380 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01381 //  const int& ii);
01382   virtual void add_value(const dvector&, const dvector&, const int&,
01383     const double&, const dvector&);
01384   virtual void add_value(const dvector&, const int&);
01385   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01386 public:
01387 #if defined(USE_ADPVM)
01388   void pvm_pack(void){::pvm_pack(*this);}
01389   void pvm_unpack(void){::pvm_unpack(*this);}
01390 #endif
01391   virtual void set_value(const dvar_vector& x, const int& ii,
01392     const dvariable& pen);
01393 #if defined(USE_SHARE_FLAGS)
01394   virtual void shared_set_value_inv(const dvector&,const int&);
01395   virtual void shared_set_value(const dvar_vector&,const int&,
01396     const dvariable& pen);
01397 #endif
01398   virtual void set_value_inv(const dvector& x, const int& ii);
01399   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01400   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01401   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01402   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01403   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01404     const int& ii);
01405 
01406 public:
01407  void allocate(const ad_integer& imin,
01408    const ad_integer& imax, const ad_integer& imin2,
01409    const ad_integer& imax2, const ad_double& _bmin,
01410    const ad_double& _bmax, const ad_integer& phase_start,
01411    const char * s);
01412 
01413   param_init_bounded_matrix(void);
01414   void allocate(int rmin,int rmax,int cmin,int cmax,
01415     double _minb,double _maxb,
01416     int phase_start=1,const char * = "UNNAMED");
01417   void allocate(int rmin,int rmax,int cmin,int cmax,
01418     double _minb,double _maxb,const char * = "UNNAMED");
01419 
01420   void allocate(const ad_integer& rmin, const ad_integer& rmax,
01421     const index_type& cmin,
01422                 const index_type& cmax, double _minb, double _maxb,
01423                 const char * = "UNNAMED");
01424   void allocate(const ad_integer& rmin, const ad_integer& rmax,
01425     const index_type& cmin,
01426                 const index_type& cmax, double _minb, double _maxb,
01427                 int phase_start = 1, const char * = "UNNAMED");
01428 };
01429 
01434 class data_int : public model_name_tag
01435 {
01436 protected:
01437   int val;
01438   data_int& operator=(const int);
01439   void allocate(int n,const char * s="UNNAMED");
01440   void allocate(const char * s="UNNAMED");
01441   void allocate(init_xml_doc&, const char * s="UNNAMED");
01442   friend class model_data;
01443   friend class model_parameters;
01444   friend int operator + (int n,data_int v);
01445   friend int operator + (data_int v,int n);
01446   friend int operator + (data_int v,data_int n);
01447 public:
01448   operator int() {return val;}
01449 #if !defined (__BORLANDC__)
01450   //operator const int() const {return val;}
01451 #endif
01452   virtual ~data_int(){;}
01453 };
01454 
01459 class named_adstring : public adstring, public model_name_tag
01460 {
01461 protected:
01462   void allocate(const char * s1,const char * s="UNNAMED");
01463   void operator = (const adstring&);
01464   void operator = (const char *);
01465 };
01466 
01471 class named_line_adstring : public line_adstring, public model_name_tag
01472 {
01473 protected:
01474   void allocate(const char * s1,const char * s="UNNAMED");
01475   void operator = (const adstring&);
01476   void operator = (const char *);
01477 };
01478 
01483 class init_adstring: public named_adstring
01484 {
01485 public:
01486   void allocate(const char * s="UNNAMED");
01487 };
01488 
01493 class init_line_adstring: public named_line_adstring
01494 {
01495 public:
01496   void allocate(const char * s="UNNAMED");
01497 };
01498 
01503 class dll_named_adstring : public named_adstring
01504 {
01505   char ** d;
01506 public:
01507   void allocate(char ** ps1,const char * s="UNNAMED");
01508   void operator = (const adstring&);
01509   void operator = (const char *);
01510   ~dll_named_adstring();
01511   dll_named_adstring(void){d=NULL;}
01512 };
01513 
01518 class dll_data_int : public data_int
01519 {
01520 public:
01521   int *pi;
01522   void allocate(int *_pi,const char * s);
01523   virtual ~dll_data_int();
01524 };
01525 
01530 class data_matrix : public named_dmatrix
01531 {
01532 public:
01533   data_matrix(void) : named_dmatrix() {;}
01534   data_matrix& operator=(const dmatrix& m);
01535   data_matrix& operator=(const double& m);
01536 private:
01537   void allocate(int rmin,int rmax,int cmin,int cmax,
01538     const char * = "UNNAMED");
01539   void allocate(int rmin, int rmax, const ivector& cmin, const ivector& cmax,
01540     const char * = "UNNAMED");
01541   void allocate(int rmin, int rmax, const ivector& cmin, int cmax,
01542     const char * = "UNNAMED");
01543   void allocate(int rmin, int rmax, int cmin, const ivector& cmax,
01544     const char * = "UNNAMED");
01545   void allocate(init_xml_doc&, const char * = "UNNAMED");
01546   friend class model_data;
01547 };
01548 
01553 class dll_data_matrix : public data_matrix
01554 {
01555   double * d;
01556 public:
01557   void allocate(double * _d,int rmin,int rmax,int cmin,int cmax,
01558     const char * _s = "UNNAMED");
01559   virtual ~dll_data_matrix();
01560   dll_data_matrix& operator=(const dmatrix &m);
01561   dll_data_matrix& operator=(const double &m);
01562 };
01563 
01568 class data_3array : public named_d3_array
01569 {
01570 public:
01571   data_3array(void) : named_d3_array() {;}
01572 private:
01573   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
01574     const char * ="UNNAMED");
01575   void allocate(int hsl, int hsu, const ivector& rmin, int rmax,
01576     int cmin,int cmax,const char * ="UNNAMED");
01577   void allocate(int hsl, int hsu, int rmin, const ivector& rmax,
01578     int cmin,int cmax,const char * ="UNNAMED");
01579   void allocate(int hsl, int hsu, const ivector& rmin, const ivector& rmax,
01580     int cmin,int cmax,const char * ="UNNAMED");
01581   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
01582     int cmax, const char * ="UNNAMED");
01583   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
01584                 const ivector& cmax, const char * ="UNNAMED");
01585   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,
01586                 const ivector& cmax, const char * ="UNNAMED");
01587   void allocate(int hsl,int hsu, const index_type& rmin, const index_type& rmax,
01588     const index_type& cmin, const index_type& cmax, const char * ="UNNAMED");
01589   friend class model_data;
01590 };
01591 
01596 class dll_data_3array : public data_3array
01597 {
01598   double * d;
01599 public:
01600   void allocate(double * _d,int hmin,int hmax,int rmin,int rmax,
01601     int cmin,int cmax,const char * _s = "UNNAMED");
01602   dll_data_3array& operator=(const d3_array &);
01603   virtual ~dll_data_3array();
01604   friend class model_data;
01605 };
01606 
01611 class data_3iarray : public named_i3_array
01612 {
01613   data_3iarray(void) : named_i3_array() {;}
01614   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
01615     const char * ="UNNAMED");
01616   void allocate(int hsl, int hsu,const index_type& rmin, const index_type& rmax,
01617     const index_type& cmin, const index_type& cmax, const char * ="UNNAMED");
01618   friend class model_data;
01619 };
01620 
01625 class data_5array : public named_d5_array
01626 {
01627   data_5array(void) : named_d5_array() {;}
01628   void allocate(int hhsl,int hhsu,
01629     int hhhsl,int hhhsu,
01630     int hsl,int hsu,int rmin,int rmax,
01631     int cmin,int cmax,const char * ="UNNAMED");
01632   void allocate(ad_integer hhhsl, ad_integer hhhsu, const index_type& hhsl,
01633     const index_type& hhsu, const index_type& hsl, const index_type& hsu,
01634     const index_type& rmin, const index_type& rmax, const index_type& cmin,
01635     const index_type& cmax, const char * ="UNNAMED");
01636   friend class model_data;
01637 };
01638 
01643 class data_4array : public named_d4_array
01644 {
01645   data_4array(void) : named_d4_array() {;}
01646   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,int rmax,
01647     int cmin,int cmax,const char * ="UNNAMED");
01648   void allocate(ad_integer hhsl, ad_integer hhsu, const index_type& hsl,
01649     const index_type& hsu, const index_type& rmin, const index_type& rmax,
01650     const index_type& cmin, const index_type& cmax, const char * ="UNNAMED");
01651   friend class model_data;
01652 };
01653 
01658 class data_imatrix : public named_imatrix
01659 {
01660   data_imatrix(void) : named_imatrix() {;}
01661   void allocate(int rmin,int rmax,int cmin,int cmax,const char * ="UNNAMED");
01662   void allocate(int rmin, int rmax, const index_type&, const index_type& cmax,
01663     const char * ="UNNAMED");
01664   friend class model_data;
01665 };
01666 
01671 class data_vector : public named_dvector
01672 {
01673 public:
01674   data_vector& operator=(const dvector& m);
01675   data_vector& operator=(const double m);
01676   data_vector(void) : named_dvector() {;}
01677 private:
01678   void allocate(int imin,int imax,const char * ="UNNAMED");
01679   void allocate(int imin, const ivector& imax, const char * ="UNNAMED");
01680   void allocate(init_xml_doc&, const char * ="UNNAMED");
01681   friend class model_data;
01682 };
01683 
01688 class dll_data_vector : public data_vector
01689 {
01690 public:
01691   double * pd;
01692   void allocate(double * pd,int imin,int imax,const char * ="UNNAMED");
01693   void allocate(double *pd, int imin, const ivector& imax,
01694     const char * ="UNNAMED");
01695   virtual ~dll_data_vector();
01696   dll_data_vector& operator = (const dvector& x);
01697   dll_data_vector& operator = (const double& x);
01698 };
01699 
01704 class data_ivector : public named_ivector
01705 {
01706 public:
01707   data_ivector(void) : named_ivector() {;}
01708 private:
01709   void allocate(int imin,int imax,const char * ="UNNAMED");
01710   friend class model_data;
01711 };
01712 
01717 class data_number : public model_name_tag
01718 {
01719 protected:
01720   double val;
01721   void allocate(const char * ="UNNAMED");
01722 public:
01723   void report_value(void);
01724   operator double() {return val;}
01725   double& value(void) {return val;}
01726   void initialize(void) {val=0.0;}
01727   friend class model_data;
01728   data_number & operator=(const double& m);
01729 };
01730 
01735 class dll_data_number : public data_number
01736 {
01737 public:
01738   double * pd;
01739   void allocate(double *_pd,const char * s);
01740   virtual ~dll_data_number();
01741   dll_data_number & operator=(const double& m);
01742 };
01743 
01744 typedef dvariable (model_parameters::*PMF) (const dvariable&);
01745 typedef dvariable (model_parameters::*PMFI) (const dvariable&,int n);
01746 typedef dvariable (model_parameters::*PMFVI) (const dvar_vector&,int n);
01747 typedef void (model_parameters::*PMFVIV4) (const dvar_vector&,int n,
01748   dvariable& f1, const dvariable& f2, const dvariable& f3, const dvariable& f4);
01749 
01750   class init_df1b2vector;
01751   class df1b2vector;
01752   class df1b2variable;
01753 
01758 class function_minimizer
01759 {
01760 public:
01761   static int bad_step_flag;
01762   static int likeprof_flag;
01763   static int first_hessian_flag;
01764   static int test_trust_flag;
01765   static int random_effects_flag;
01766   dmatrix * negdirections;
01767   static int negative_eigenvalue_flag;
01768   static int inner_opt_flag;
01769   static int inner_opt(void);
01770   laplace_approximation_calculator * lapprox;
01771   dvector * multinomial_weights;
01772   void set_multinomial_weights(dvector&d);
01773   //init_df1b2vector* py;
01774   virtual void AD_uf_inner();
01775   virtual void AD_uf_outer();
01776   virtual void user_function();
01777   void pre_user_function(void);
01778   //void df1b2_pre_user_function(void);
01779   //virtual void user_function(const init_df1b2vector& x,df1b2variable& f);
01780   //static int hesstype;
01781   //int set_hessian_type(int);
01782   void hess_routine_noparallel_randomeffects(void);
01783   void other_separable_stuff_begin(void);
01784   void other_separable_stuff_end(void);
01785   void begin_df1b2_funnel(void){;}
01786   void setup_quadprior_calcs(void){;}
01787   void end_df1b2_funnel(void){;}
01788   void get_function_difference(void);
01789   void start_get_importance_sampling_comnponent(void);
01790   void end_get_importance_sampling_comnponent(void);
01791   int spminflag;
01792   int repeatminflag;
01793   int mcmc2_flag;
01794   int robust_hybrid_flag;
01795   int ifn;
01796   int maxfn;
01797   int iprint;
01798   double crit;
01799   int imax;
01800   double dfn;
01801   int iexit;
01802   int ihflag;
01803   int ihang;
01804   int scroll_flag;
01805   int maxfn_flag;
01806   int quit_flag;
01807   double min_improve;
01808   void pre_userfunction(void);
01809   virtual void userfunction(void)=0;
01810   virtual void allocate(void){;}
01811   static named_dvar_vector * ph;
01812   static named_dvar_vector * pg;
01813 protected:
01814   double ffbest;
01815 private:
01816   gradient_structure * pgs;
01817   adstring_array param_labels;
01818   ivector param_size;
01819 protected:
01820   void report_function_minimizer_stats(void){;}
01821   virtual void report(const dvector& gradients){;};
01822   static dvector convergence_criteria;
01823   static dvector maximum_function_evaluations;
01824   static int sd_flag;
01825   static adstring user_data_file;
01826   static adstring user_par_file;
01827   static int have_constraints;
01828 public:
01829   virtual dvariable user_randeff(const dvar_vector& x);
01830   virtual dvar_vector user_dfrandeff(const dvar_vector& x);
01831   virtual dvar_matrix user_d2frandeff(const dvar_vector& x);
01832   void limited_memory_quasi_newton(const independent_variables&,int);
01833   void limited_memory_quasi_newton(const independent_variables&,int,int);
01834   void limited_memory_quasi_newton(double& f, const independent_variables&,
01835     int, int, int,double);
01836   function_minimizer(long int sz=0L);
01837   void likeprof_routine(double global_min);
01838   virtual ~function_minimizer();
01839   virtual void other_calculations(void){;}
01840   virtual void final_calcs(void){;}
01841   virtual void minimize(void);
01842   virtual void constraints_minimize(void);
01843   virtual void between_phases_calculations(void){;}
01844   void computations(int argc,char * argv[]);
01845   void computations1(int argc,char * argv[]);
01846   void computations_np(int argc,char * argv[]);
01847   void computations(void);
01848   void hess_routine(void);
01849   void hess_routine_noparallel(void);
01850   void hess_routine_master(void);
01851   void hess_routine_slave(void);
01852   void constraint_report(void);
01853   void pvm_slave_likeprof_routine(void);
01854   void pvm_master_function_evaluation_profile(double& f,
01855     independent_variables& x,const dvector & g,int nvar,int iprof,double weight,
01856     double new_value);
01857   void pvm_slave_prof_minimize(int underflow_flag);
01858   void pvm_master_prof_minimize(int iprof, double sigma,
01859     double new_value, const double& _fprof, const int underflow_flag,
01860     double global_min, const double& _penalties,
01861     const double& _final_weight);
01862   //void pvm_master_function_evaluation_profile(double& f,
01863    // independent_variables& x,dvector & g,int nvar,int iprof);
01864   void pvm_slave_function_evaluation(void);
01865   void pvm_slave_function_evaluation_no_derivatives(void);
01866   void pvm_slave_function_evaluation_noder(void);
01867   void pvm_master_function_evaluation_no_derivatives(double& f,
01868     independent_variables& x,int nvar);
01869   void pvm_master_function_evaluation(double& f,
01870     independent_variables& x,const dvector & g,int nvar);
01871   dmatrix dep_hess_routine(const dvariable& dep);
01872   void top_mcmc_routine(int,int,double,int);
01873   void mcmc_routine(int,int,double,int);
01874   void sgibbs_mcmc_routine(int,int,double,int);
01875   void hybrid_mcmc_routine(int,int,double,int);
01876   double pvm_master_get_monte_carlo_value(int nvar,
01877     const dvector& x);
01878   void pvm_slave_get_monte_carlo_value(int nvar);
01879   void mcmc_eval(void);
01880   //void hess_routine_and_constraint(int);
01881   void hess_routine_and_constraint(int iprof, const dvector& g,
01882     dvector& fg);
01883   dmatrix diag_hess_routine(void);
01884   void hess_inv(void);
01885   void depvars_routine(void);
01886   void sd_routine(void);
01887   int ef_(double * f, double * x);
01888   int constrained_minimization2(int _n,int _nh, int _ng,dvector& __x);
01889   static int constraint_exit_number;
01890   void get_bigS(int ndvar,int nvar1,int nvar,
01891     dmatrix& S,dmatrix& BS,dvector& scale);
01892 
01893 #ifdef CURVE_CORRECT
01894   void constraint_hess_routine(int ip);
01895   void get_curvature_correction_factors(int ip,
01896     dvector& g, const int underflow_flag, const dvector& eigenvals,
01897     dvector& curvcor);
01898 #endif
01899   double projected_jacobian(const dvector& g, const dvector& xscale);
01900 
01901   dvariable random_effects_maximization(const dvar_vector& v);
01902   void prof_minimize_re(int iprof, double sigma,
01903     double new_value, const double& fprof, const int underflow_flag,
01904     double global_min, const double& penalties, const double& final_weight);
01905   void prof_minimize(int iprof, double sigma,
01906     double new_value, const double& fprof, const int underflow_flag,
01907     double global_min, const double& penalties, const double& final_weight);
01908   void begin_gauss_hermite_stuff(void);
01909   void begin_funnel_stuff(void);
01910   void end_gauss_hermite_stuff(void);
01911 
01912   void prof_hess_routine(int ip,double new_value,double final_weight);
01913 
01914   void quasi_newton_minimizer1(int nvar,double _crit,
01915     double& f, const independent_variables& x,const dvector& g);
01916 
01917   double hess_determinant(int underflow_flag);
01918 
01919 #ifndef CURVE_CORRECT
01920   void normalize_posterior_distribution(double udet,
01921     const dvector& siglevel, const ofstream& ofs2, int num_pp,
01922     const dvector& all_values, const dmatrix& actual_value, double global_min,
01923     int offset, const dmatrix& lprof, const dmatrix& ldet, const dmatrix& xdist,
01924     const dmatrix& penalties);
01925   //  dmatrix& penalties, const hdmatrix& lg_jacob);
01926 #else
01927   void normalize_posterior_distribution(double udet,
01928     const dvector& siglevel, const ofstream& ofs2, int num_pp,
01929     const dvector& all_values, const dmatrix& actual_value,
01930     double global_min,
01931     int offset, const dmatrix& lprof, const dmatrix& ldet,
01932     const dmatrix& xdist,
01933     const d3_array& eigenvals,d3_array& curvcor);
01934 #endif
01935   void get_particular_grad(int iprof, int nvar, const dvector& fg,
01936     const dvector& g);
01937   double projected_hess_determinant(const dvector& g, const int underflow_flag,
01938     const dvector& xscale, const double& ln_det_proj_jac);
01939 //  double projected_hess_determinant(const dvector& fg, const dvector& g,
01940   //const int underflow_flag, const dvector& xscale,
01941   // const double& ln_det_proj_jac);
01942   double projected_hess_determinant(const dvector& g,const int underflow_flag);
01943   double projected_hess_determinant(const dmatrix& hh, const dvector& g,
01944     const int underflow_flag);
01945  //double projected_hess_determinant(const dvector& g, const int underflow_flag,
01946     //dvector& xscale, const double& ln_det_proj_jac);
01947 
01948   double projected_hess_determinant(const dvector& fg, const dvector& g,
01949     const int underflow_flag, const dvector& curvscale, const dvector& xscale,
01950     double& ln_det_proj_jac, const double& tmp, const dmatrix& hesses);
01951   double diag_projected_hess_determinant(const dvector& g,
01952     const int underflow_flag, dmatrix& dh);
01953     double unrestricted_hess_determinant(void);
01954   void monte_carlo_routine(void);
01955   double get_monte_carlo_value(int nvar, const independent_variables& x);
01956   double get_monte_carlo_value(int nvar, const independent_variables& x,
01957     dvector& g);
01958   double get_hybrid_monte_carlo_value(int nvar,const independent_variables& x,
01959     dvector& g);
01960   void mcmc_computations(void);
01961   void pvm_slave_mcmc_computations(void);
01962   void pvm_master_mcmc_computations(void);
01963   double get_monte_carlo_value(int nvar, const dvector& x);
01964   void sob_routine(int nmcmc,double dscale,int restart_flag);
01965   void sobol_importance_routine(int nmcmc,int iseed0,double dscale,
01966   int restart_flag);
01967   void pvm_master_mcmc_routine(int nmcmc,int iseed0,double dscale,
01968     int restart_flag);
01969 #if defined(USE_ADPVM)
01970   void pvm_slave_mcmc_routine(void);
01971 #else
01972   void pvm_slave_mcmc_routine(void) {}
01973 #endif
01974   void trust_region_update(int nvar,int _crit,
01975     independent_variables& x,const dvector& _g,const double& _f);
01976 
01977   void multint4(int n, const dvar_vector& a, const dvar_vector& b,
01978     const dvar_vector& h, double al, int m, const dvariable& e,
01979     const dvariable& aint1, const dvariable& aint2, dvariable& aint3,
01980     const dvariable& aint4, const int& key, PMFVIV4 f);
01981 
01982   void multint(int n, const dvar_vector& a, const dvar_vector& b,
01983     const dvar_vector& h, double al, int m, const dvariable& e,
01984     const dvariable& aint, const int& key, PMFVI f);
01985 
01986   virtual void set_runtime(void);
01987   virtual void set_runtime_maxfn(const char *);
01988   virtual void set_runtime_crit(const char *);
01989   dvariable adromb(PMF,double a,double b,int ns=9);
01990   dvariable adromb(PMF, const dvariable& a, double b, int ns = 9);
01991   dvariable adromb(PMF, const dvariable& a, const dvariable& b, int ns = 9);
01992   dvariable adromb(PMF, double a, const dvariable& b, int ns = 9);
01993 
01994   dvariable adrombo(PMF,double a,double b,int ns=9);
01995   dvariable adrombo(PMF, const dvariable& a, double b, int ns = 9);
01996   dvariable adrombo(PMF, const dvariable& a, const dvariable& b, int ns = 9);
01997   dvariable adrombo(PMF, double a, const dvariable& b,int ns = 9);
01998 
01999   dvariable trapzd(void*,double a,double b,int n);
02000   dvariable trapzd(PMF,double a,double b,int n);
02001   dvariable trapzd(PMF, double a, const dvariable& b, int n);
02002   dvariable trapzd(PMF, const dvariable& a, double b, int n);
02003   dvariable trapzd(PMF, const dvariable& a, const dvariable& b, int n);
02004 
02005   dvariable midpnt(PMF,double a,double b,int n);
02006   dvariable midpnt(PMF, double a, const dvariable& b, int n);
02007   dvariable midpnt(PMF, const dvariable& a, double b, int n);
02008   dvariable midpnt(PMF, const dvariable& a, const dvariable& b, int n);
02009 
02010   virtual void * mycast() { return (void*)this;}
02011 
02012   void neldmead(int n, dvector& _start, dvector& _xmin, double *ynewlo,
02013     double reqmin, double delta,int *icount, int *numres, int *ifault);
02014   void adamoeba(const dmatrix& p, const dvector& y, int ndim, double ftol,
02015      int maxfn);
02016   void set_initial_simplex(const dmatrix& p, const dvector& y, int nvar,
02017     const dvector& x, double delta);
02018   double amxxx(const dmatrix& p, const dvector& y, const dvector& psum,
02019     int ndim, int ihi, double fac);
02020   friend class equality_constraint_vector;
02021   friend class inequality_constraint_vector;
02022   void quasi_newton_block(int nvar,int crit,independent_variables& x,
02023     const dvector& g,const double& f);
02024   void limited_memory_quasi_newton_block(int nvar,int _crit,
02025     independent_variables& x,const dvector& _g,const double& _f,
02026     int nsteps);
02027 
02028   void function_evaluation_block_pvm_slave_random_effects(int nvar,int _crit,
02029     independent_variables& x,const dvector& g,const double& f);
02030   void quasi_newton_block_pvm_master_random_effects(int nvar,int _crit,
02031     independent_variables& x,const dvector& g,const double& f);
02032   void function_evaluation_block_pvm_slave_random_effects(void);
02033   void hess_routine_random_effects(void);
02034   void quasi_newton_block_pvm_master(int nvar,int _crit,
02035     independent_variables& x,const dvector& g,const double& f);
02036   void hess_routine_noparallel_random_effects(void);
02037 #if defined(USE_ADPVM)
02038   void function_evaluation_block_pvm_slave(void);
02039   void hess_routine_slave_random_effects(void);
02040 #endif
02041   dvariable do_gauss_hermite_integration(void);
02042   void end_df1b2_funnel_stuff(void);
02043 
02044 private:
02045   dvariable do_gauss_hermite_integration_multi(void);
02046 };
02047 
02048 cifstream& operator>>(const cifstream& s, const param_init_number& x);
02049 cifstream& operator>>(const cifstream& s, const param_init_vector& v);
02050 cifstream& operator>>(const cifstream& s, const param_init_matrix& m);
02051 ostream& operator<<(const ostream& s, const label_class& lc);
02052 
02057 class stddev_params
02058 {
02059 public:
02060   // this should be a resizeable array
02061   static stddev_params * stddevptr[150];
02062   // this should be a resizeable array
02063   static stddev_params * stddev_number_ptr[150];
02064   static void get_all_sd_values(const dvector& x, const int& ii);
02065   static int num_stddev_params;
02066   static int num_stddev_number_params;
02067   static ivector copy_all_number_offsets(void);
02068   void allocate(void){;};
02069   static int num_stddev_calc(void);
02070   static int num_stddev_number_calc(void);
02071   static void get_stddev_number_offset(void);
02072 public:
02073   stddev_params(void){;}
02074   virtual void setindex(int);
02075   virtual int getindex(void);
02076   virtual int size_count(void)=0; // get the number of active parameters
02077   virtual void set_dependent_variables(void)=0;
02078   virtual void copy_value_to_vector(const dvector&,const int&) = 0;
02079   virtual void get_sd_values(const dvector& x, const int& ii) = 0;
02080   //get the number of active parameters
02081   static void copy_all_values(const dvector& x, const int& ii);
02082   //get the number of active parameters
02083   static void copy_all_number_values(const dvector& x, const int& ii);
02084   virtual void add_to_list(void);
02085   virtual void add_to_gui_list(void);
02086   virtual const char * label()=0;
02087   friend class function_minimizer;
02088 };
02089 
02094 class likeprof_params
02095 {
02096   double stepsize;
02097   int    stepnumber;
02098 protected:
02099 public:
02100   static likeprof_params * likeprofptr[500]; // this should be a
02101                                                // resizeable array
02102   static int num_likeprof_params;
02103   void allocate(void){;};
02104   static int num_stddev_calc(void);
02105 public:
02106   likeprof_params(void);
02107   virtual void add_to_list(void);
02108   virtual const char * label()=0;
02109   virtual dvariable variable(void)=0;
02110   virtual double& get_sigma(void)=0;
02111   virtual double get_value(void)=0;
02112   double get_stepsize(void);
02113   int get_stepnumber(void);
02114   void set_stepsize(double);
02115   void set_stepnumber(int);
02116   friend class function_minimizer;
02117 };
02118 
02123 class param_stddev_vector: public named_dvar_vector , stddev_params
02124 {
02125   dvector sd;
02126     virtual int size_count(void); // get the number of active parameters
02127     virtual const char * label(void);
02128     param_stddev_vector();
02129     void allocate(int imin,int imax,const char * s="UNNAMED");
02130     virtual void set_dependent_variables(void);
02131     friend class model_parameters;
02132     virtual void copy_value_to_vector(const dvector& x, const int& ii);
02133     virtual void get_sd_values(const dvector& x, const int& ii);
02134   param_stddev_vector& operator=(const dvar_vector& m);
02135   param_stddev_vector& operator=(const dvector& m);
02136   param_stddev_vector& operator=(const double m);
02137 };
02138 
02143 class param_stddev_number: public named_dvariable , public stddev_params
02144 {
02145   double sd;
02146   int index;
02147   void allocate(const char *s="UNNAMED");
02148   virtual void setindex(int);
02149   virtual int getindex(void);
02150   virtual int size_count(void); // get the number of active parameters
02151   virtual const char * label(void);
02152   virtual void copy_value_to_vector(const dvector& x, const int& ii);
02153   virtual void get_sd_values(const dvector& x, const int& ii);
02154 protected:
02155   param_stddev_number();
02156   friend class model_parameters;
02157   virtual void set_dependent_variables(void);
02158   param_stddev_number& operator=(const prevariable&);
02159   param_stddev_number& operator=(const double);
02160 };
02161 
02166 class param_likeprof_number: public param_stddev_number ,
02167     public likeprof_params
02168 {
02169     double sigma;
02170     void allocate(const char *s="UNNAMED");
02171     virtual int size_count(void); // get the number of active parameters
02172     virtual const char * label(void);
02173     virtual double& get_sigma(void){return sigma;}
02174     virtual double get_value(void){return value(*this);}
02175     //void copy_value_to_vector(const dvector& x, const int& ii);
02176     virtual dvariable variable(void){ return dvariable(*this);}
02177     param_likeprof_number();
02178     friend class model_parameters;
02179 public:
02180     param_likeprof_number& operator=(const prevariable&);
02181     param_likeprof_number& operator=(const double);
02182 };
02183 
02188 class param_stddev_matrix: public named_dvar_matrix , stddev_params
02189 {
02190   dmatrix sd;
02191     virtual int size_count(void);
02192     //virtual void read_value(void);
02193     virtual const char * label(void);
02194     void allocate(int rmin,int rmax,int cmin,int cmax,
02195     const char * s="UNNAMED");
02196   param_stddev_matrix(void);
02197   friend class model_parameters;
02198   virtual void set_dependent_variables(void);
02199   virtual void get_sd_values(const dvector& x, const int& ii);
02200   virtual void copy_value_to_vector(const dvector& x, const int& ii);
02201   param_stddev_matrix& operator=(const double m);
02202   param_stddev_matrix& operator=(const dmatrix& m);
02203   param_stddev_matrix& operator=(const dvar_matrix& m);
02204 };
02205 
02210   class objective_function_value : public named_dvariable
02211   {
02212   public:
02213     static objective_function_value * pobjfun;
02214     static double fun_without_pen;
02215     static double gmax;
02216     objective_function_value();
02217     objective_function_value& operator=(const prevariable& v);
02218     objective_function_value& operator=(const double v);
02219   };
02220 
02221   int withinbound(int lb,int n,int ub);
02222 
02223 double cumd_cauchy(const double& x);
02224 double density_cauchy(const double& x);
02225 double log_density_cauchy(const double& x);
02226 double inv_cumd_cauchy(const double& x);
02227 
02228 double cumd_mixture(const double& x);
02229 double inv_cumd_mixture(const double& y);
02230 double cumd_mixture_02(const double& x);
02231 double inv_cumd_mixture_02(const double& y);
02232 
02233 #if defined _ADM_HIGHER_ARRAYS__
02234 
02239 class param_init_matrix: public named_dvar_matrix,public initial_params
02240 {
02241   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
02242 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
02243 //  const int& ii);
02244   virtual void add_value(const dvector&, const dvector&, const int&,
02245     const double&, const dvector&);
02246   virtual void add_value(const dvector&, const int&);
02247   virtual void get_jacobian(const dvector&, const dvector&, const int&);
02248 public:
02249   virtual void set_value(const dvar_vector& x, const int& ii,
02250     const dvariable& pen);
02251   virtual void copy_value_to_vector(const dvector& x, const int& ii);
02252   virtual void restore_value_from_vector(const dvector&, const int&);
02253   virtual void set_value_inv(const dvector& x, const int& ii);
02254   virtual int size_count(void);
02255   virtual void save_value(void);
02256   virtual void bsave_value(void);
02257   virtual void save_value(const ofstream& ofs, int prec);
02258   virtual void restore_value(const ifstream& ifs);
02259   void report_value(void);
02260   //virtual void read_value(void);
02261   virtual const char * label(void);
02262   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
02263   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
02264   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
02265   virtual void hess_scale(const dvector& d, const dvector& x, const int& ii){;};
02266 
02267 public:
02268   void allocate(int rmin,int rmax,int cmin,int cmax,
02269     int phase_start=1,const char * = "UNNAMED");
02270   void allocate(int rmin,int rmax,int cmin,int cmax,
02271     const char * = "UNNAMED");
02272   param_init_matrix(void);
02273   param_init_matrix& operator = (const dmatrix& m);
02274   param_init_matrix& operator = (const dvar_matrix& m);
02275   param_init_matrix& operator = (const dvariable& m);
02276   param_init_matrix& operator = (const double& m);
02277 };
02278 #endif // #if defined _ADM_HIGER_ARRAYS__
02279 
02284 class param_init_d3array: public named_dvar3_array,public initial_params
02285 {
02286 public:
02287 #if defined(USE_SHARE_FLAGS)
02288   virtual void setshare(const index_type& sf,const index_type& af);
02289   virtual void shared_set_value_inv(const dvector&,const int&);
02290   virtual void shared_set_value(const dvar_vector&,const int&,
02291     const dvariable& pen);
02292   virtual int shared_size_count(void); // get the number of active parameters
02293 #endif
02294 
02295   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
02296     const int& ii);
02297   virtual void dev_correction(const dmatrix&, const int&);
02298   virtual void curv_scale(const dvector& d, const dvector& x, const int& ii);
02299   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
02300 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
02301 //  const int& ii);
02302   virtual void add_value(const dvector&, const dvector&, const int&,
02303     const double&, const dvector&);
02304   virtual void add_value(const dvector&, const int&);
02305   virtual void get_jacobian(const dvector&, const dvector&, const int&);
02306   virtual void set_value(const dvar_vector& x, const int& ii,
02307     const dvariable& pen);
02308   virtual void copy_value_to_vector(const dvector& x, const int& ii);
02309   virtual void restore_value_from_vector(const dvector&,const int&);
02310   virtual void set_value_inv(const dvector& x, const int& ii);
02311   virtual int size_count(void);
02312   virtual void save_value(void);
02313   virtual void bsave_value(void);
02314   virtual void save_value(const ofstream& ofs, int prec);
02315   virtual void restore_value(const ifstream& ifs);
02316   void report_value(void);
02317   //virtual void read_value(void);
02318   virtual const char * label(void);
02319   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
02320   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
02321   virtual void hess_scale(const dvector& d, const dvector& x, const  int& ii);
02322 
02323 public:
02324 #if defined(USE_ADPVM)
02325   virtual void pvm_pack(void) { cerr << "Error" << endl; ad_exit(1);}
02326   virtual void pvm_unpack(void) { cerr << "Error" << endl; ad_exit(1);}
02327 #endif
02328 
02329   void allocate(const ad_integer& sl,const ad_integer& sh,
02330     const index_type& nrl,const index_type& nrh,
02331     const index_type& ncl,const index_type& nch,const char * s="UNNAMED");
02332 
02333   void allocate(const ad_integer& sl,const ad_integer& sh,
02334     const index_type& nrl,const index_type& nrh,
02335     const index_type& ncl,const index_type& nch,int phase_start=1,
02336     const char * s="UNNAMED");
02337 
02338   void allocate(int smin,int smax,int rmin,int rmax,int cmin,int cmax,
02339     int phase_start=1,const char * = "UNNAMED");
02340   void allocate(int smin,int smax,int rmin,int rmax,int cmin,int cmax,
02341     const char * = "UNNAMED");
02342   param_init_d3array(void);
02343 };
02344 
02349 class dll_param_init_d3array: public param_init_d3array
02350 {
02351   double * d;
02352 public:
02353   void allocate(double* _d,int hmin,int hmax,
02354     int rmin,int rmax,int cmin,int cmax,
02355     int phase_start=1,const char * = "UNNAMED");
02356   void allocate(double * _d,int hmin,int hmax,int rmin,int rmax,
02357     int cmin,int cmax,const char * = "UNNAMED");
02358   virtual ~dll_param_init_d3array();
02359   dll_param_init_d3array(){d=NULL;}
02360   dll_param_init_d3array& operator=(const d3_array&);
02361   dll_param_init_d3array& operator=(const dvar3_array&);
02362 };
02363 
02368 class dll_param_d3array: public named_dvar3_array
02369 {
02370   double * d;
02371 public:
02372   void allocate(double* _d,int hmin,int hmax,
02373     int rmin,int rmax,int cmin,int cmax,
02374     int phase_start=1,const char * = "UNNAMED");
02375   void allocate(double * _d,int hmin,int hmax,int rmin,int rmax,
02376     int cmin,int cmax,const char * = "UNNAMED");
02377   virtual ~dll_param_d3array();
02378   dll_param_d3array(){d=NULL;}
02379   dll_param_d3array& operator=(const d3_array&);
02380   dll_param_d3array& operator=(const dvar3_array&);
02381 };
02382 
02383 
02384 //double set_value_mc(const double& x, const double fmin, const double fmax);
02385 
02386 void set_value_mc(const dvar_vector& x, const dvar_vector& v, const int& ii,
02387   const double fmin, const double fmax);
02388 
02389 double set_value_inv_mc(double v,double fmin,double fmax);
02390 
02391 double set_value_inv_mc(const prevariable& v, double fmin, double fmax);
02392 
02393 void set_value_inv_mc(const dvar_vector& x, const dvector& v, const int& ii,
02394   const double fmin, const double fmax);
02395 
02396 //double set_value_inv_mc(const dvector&, const dvector& x, int ii, double minb,
02397 //  double maxb);
02398 
02399 double set_value_mc(double z,double min,double max);
02400 double ndfboundp( double x, double fmin, double fmax,const double& fpen);
02401 double ndfboundp_mc( double x, double fmin, double fmax,const double& fpen);
02402 
02403 void copy_value_from_vector(const double& _sd,const dvector& x,const int & _ii);
02404 void copy_value_from_vector(const dvector& _sd,const dvector& x,
02405   const int & _ii);
02406 void copy_value_from_vector(const dmatrix& _sd,const dvector& x,
02407   const int & _ii);
02408 
02409 dvector bounded_multivariate_normal(int nvar, const dvector& a1,
02410   const dvector& b1, dmatrix& ch, const double& _wght,
02411   const random_number_generator & rng);
02412 
02413 void bounded_multivariate_normal_mcmc(int nvar, const dvector& a1,
02414   const dvector& b1, dmatrix& ch, const double& wght, const dvector& y,
02415   const random_number_generator & rng);
02416 
02417 void bounded_multivariate_uniform_mcmc(int nvar, const dvector& a1,
02418   const dvector& b1, dmatrix& ch, const double& wght, const dvector& y,
02419   const random_number_generator & rng);
02420 
02421 dvector bounded_multivariate_normal(int nvar, const dvector& a1,
02422   const dvector& b1, dmatrix& ch, const double& lprob,
02423   random_number_generator &rng);
02424 
02425 dvector bounded_multivariate_normal_sobol(int nvar, const dvector& a1,
02426   const dvector& b1, dmatrix& ch, const double& lprob,
02427   const random_number_generator &rng);
02428 
02429 dvector probing_bounded_multivariate_normal(int nvar, const dvector& a1,
02430   const dvector& b1, dmatrix& ch, const double& lprob, double pprobe,
02431   const random_number_generator &rng);
02432 
02433 dvector bounded_multivariate_uniform(int nvar, const dvector& a1,
02434   const dvector& b1, dmatrix& ch, const double& lprob,
02435   random_number_generator &rng);
02436 
02437 dvector bounded_robust_multivariate_normal(int nvar, const dvector& a1,
02438   const dvector& b1, dmatrix& ch, const dmatrix& ch3, double contaminant,
02439   const double& _wght, random_number_generator &rng);
02440 
02441 void probing_bounded_multivariate_normal_mcmc(int nvar, const dvector& a1,
02442   const dvector& b1, dmatrix& ch, const double& wght, const dvector& v,
02443   double pprobe, const random_number_generator &rng);
02444 
02445 /*
02446 int option_match(int argc,char * argv[], const char * string);
02447 int option_match(int argc,char * argv[], const char * string, const int& nopt);
02448 int option_match(char * s, const char * string, const int& _nopt);
02449 int option_match(char * s, const char * string);
02450 */
02451 
02452 double inv_cumd_exp(double x);
02453 double cumd_exp(double x);
02454 
02455 double ffmax(double a,double b);
02456 double ffmin(double a,double b);
02457 
02458 void check_datafile_pointer(void * p);
02459 
02460 adstring get_reportfile_name(void);
02461 
02462 void ad_make_code_reentrant(void);
02463 
02464 char** parse_dll_options(char *pname, const int& argc,
02465   char * dll_options);
02466 
02467 void parse_dll_options(char *pname, const int& argc, char *dll_options,
02468   char *** pargv);
02469 
02470 char** no_dll_options(char *pname, const int& argc);
02471 
02472 void cleanup_argv(int nopt,char *** pa);
02473 
02474 void get_sp_printf(void);
02475 
02476 void do_dll_housekeeping(int argc,char ** argv);
02477 
02478 void adwait(double);
02479 
02480 int ad_get_commandline_option(const char *option_label, const int &option_value,
02481   const char * error_message);
02482 
02487  class param_init_vector_vector
02488  {
02489    param_init_vector * v;
02490    int index_min;
02491    int index_max;
02492    double_index_type * it;
02493 
02494  public:
02495   void set_scalefactor(double s);
02496   void set_scalefactor(const dvector& s);
02497   dvector get_scalefactor(void);
02498 
02499 #if defined(OPT_LIB)
02500    param_init_vector& operator [] (int i) { return v[i];}
02501    param_init_vector& operator () (int i) { return v[i];}
02502    prevariable operator () (int i,int j) { return v[i][j];}
02503 #else
02504    param_init_vector& operator [] (int i);
02505    param_init_vector& operator () (int i);
02506    prevariable operator () (int i,int j);
02507 #endif
02508 
02509    void allocate(int min1,int max1,const index_type& min,
02510      const index_type& max,const index_type& phase_start,
02511      const char * s);
02512 
02513    void allocate(int min1,int max1,const index_type& min,
02514      const index_type& max,const char * s);
02515 
02516    param_init_vector_vector();
02517    int allocated(void) { return (v!=NULL); }
02518    int indexmin(void) {return (index_min);}
02519    int indexmax(void) {return (index_max);}
02520    ~param_init_vector_vector();
02521    void set_initial_value(const double_index_type& it);
02522    void deallocate(void);
02523  };
02524 
02529  class param_init_bounded_vector_vector
02530  {
02531    param_init_bounded_vector * v;
02532    int index_min;
02533    int index_max;
02534    double_index_type * it;
02535  public:
02536   void set_scalefactor(double s);
02537   void set_scalefactor(const dvector& s);
02538   dvector get_scalefactor(void);
02539   param_init_bounded_vector_vector();
02540 #if defined(OPT_LIB)
02541    param_init_bounded_vector& operator [] (int i) { return v[i];}
02542    param_init_bounded_vector& operator () (int i) { return v[i];}
02543    prevariable operator () (int i,int j) { return v[i][j];}
02544 #else
02545    param_init_bounded_vector& operator [] (int i);
02546    param_init_bounded_vector& operator () (int i);
02547    prevariable operator () (int i,int j);
02548 #endif
02549 
02550    void allocate(int min1,int max1,
02551      const index_type& min,
02552      const index_type& max,
02553      const double_index_type& dmin,
02554      const double_index_type& dmax,
02555      const index_type& phase_start,
02556      const char * s);
02557 
02558    void allocate(int min1,int max1,
02559      const index_type& min,
02560      const index_type& max,
02561      const double_index_type& dmin,
02562      const double_index_type& dmax,
02563      const char * s);
02564 
02565    int allocated(void) { return (v!=NULL); }
02566    int indexmin(void) {return (index_min);}
02567    int indexmax(void) {return (index_max);}
02568    ~param_init_bounded_vector_vector();
02569    void deallocate(void);
02570    void set_initial_value(const double_index_type& it);
02571  };
02572 
02577  class param_init_matrix_vector
02578  {
02579    param_init_matrix * v;
02580    int index_min;
02581    int index_max;
02582    double_index_type * it;
02583  public:
02584   param_init_matrix_vector();
02585   void set_scalefactor(double s);
02586   void set_scalefactor(const dvector& s);
02587   dvector get_scalefactor(void);
02588 
02589 #if defined(OPT_LIB)
02590    param_init_matrix& operator [] (int i) { return v[i];}
02591    param_init_matrix& operator () (int i) { return v[i];}
02592    dvar_vector& operator () (int i,int j) { return v[i][j];}
02593    prevariable operator () (int i,int j,int k) { return v[i](j,k);}
02594 #else
02595    param_init_matrix& operator [] (int i);
02596    param_init_matrix& operator () (int i);
02597    dvar_vector& operator () (int i,int j);
02598    prevariable operator () (int i,int j,int k);
02599 #endif
02600    void allocate(int min0,int max0,const index_type& min,
02601      const index_type& max,const index_type& min1,
02602      const index_type& max1,const index_type& phase_start,
02603      const char * s);
02604 
02605    void allocate(int min0,int max0,const index_type& min,
02606      const index_type& max,const index_type& min1,
02607      const index_type& max1,const char * s);
02608 
02609    int allocated(void) { return (v!=NULL); }
02610    int indexmin(void) {return (index_min);}
02611    int indexmax(void) {return (index_max);}
02612    ~param_init_matrix_vector();
02613    void set_initial_value(const double_index_type& it);
02614    void deallocate(void);
02615  };
02616 
02621  class param_init_bounded_matrix_vector
02622  {
02623    param_init_bounded_matrix * v;
02624    int index_min;
02625    int index_max;
02626    double_index_type * it;
02627  public:
02628   void set_scalefactor(double s);
02629   void set_scalefactor(const dvector& s);
02630   dvector get_scalefactor(void);
02631   param_init_bounded_matrix_vector();
02632 
02633    void allocate(int min1,int max1,
02634      const index_type& min, const index_type& max, const index_type& min2,
02635      const index_type& max2, const double_index_type& dmin2,
02636      const double_index_type& dmax2, const index_type& phase_start,
02637      const char * s);
02638 
02639    void allocate(int min1,int max1,
02640      const index_type& min, const index_type& max, const index_type& min2,
02641      const index_type& max2, const double_index_type& dmin2,
02642      const double_index_type& dmax2,const char * s);
02643 
02644 #if defined(OPT_LIB)
02645    param_init_bounded_matrix& operator [] (int i) { return v[i];}
02646    param_init_bounded_matrix& operator () (int i) { return v[i];}
02647    dvar_vector& operator () (int i,int j) { return v[i][j];}
02648    prevariable operator () (int i,int j,int k) { return v[i](j,k);}
02649 #else
02650    param_init_bounded_matrix& operator [] (int i);
02651    param_init_bounded_matrix& operator () (int i);
02652    dvar_vector& operator () (int i,int j);
02653    prevariable operator () (int i,int j,int k);
02654 #endif
02655 
02656    int allocated(void) { return (v!=NULL); }
02657    int indexmin(void) {return (index_min);}
02658    int indexmax(void) {return (index_max);}
02659    ~param_init_bounded_matrix_vector();
02660    void set_initial_value(const double_index_type& it);
02661    void deallocate(void);
02662  };
02663 
02668  class param_init_number_vector
02669  {
02670    param_init_number * v;
02671    int index_min;
02672    int index_max;
02673    double_index_type * it;
02674  public:
02675   void set_scalefactor(double s);
02676   void set_scalefactor(const dvector& s);
02677   dvector get_scalefactor(void);
02678   param_init_number_vector();
02679 
02680 #if defined(OPT_LIB)
02681    param_init_number& operator [] (int i) { return v[i];}
02682    param_init_number& operator () (int i) { return v[i];}
02683 #else
02684    param_init_number& operator [] (int i);
02685    param_init_number& operator () (int i);
02686 #endif
02687 
02688    void allocate(int min1,int max1,const index_type& phase_start,
02689      const char * s);
02690 
02691    void allocate(int min1,int max1,const char * s);
02692 
02693    int allocated(void) { return (v!=NULL); }
02694    int indexmin(void) {return (index_min);}
02695    int indexmax(void) {return (index_max);}
02696    ~param_init_number_vector();
02697    void set_initial_value(const double_index_type& it);
02698    void deallocate(void);
02699  };
02700 
02708  class data_matrix;
02709  class param_init_bounded_number_vector
02710  {
02711    param_init_bounded_number * v;
02712    int index_min;
02713    int index_max;
02714    double_index_type * it;
02715  public:
02716   void set_scalefactor(double s);
02717   void set_scalefactor(const dvector& s);
02718   dvector get_scalefactor(void);
02719   param_init_bounded_number_vector();
02720 
02721 #if defined(OPT_LIB)
02722    param_init_bounded_number& operator [] (int i) { return v[i];}
02723    param_init_bounded_number& operator () (int i) { return v[i];}
02724 #else
02725    param_init_bounded_number& operator [] (int i);
02726    param_init_bounded_number& operator () (int i);
02727 #endif
02728 
02729    void allocate(int min1,int max1,const double_index_type & bmin,
02730      const double_index_type & bmax,const index_type& phase_start,
02731      const char * s);
02732 
02733    void allocate(int min1,int max1,const double_index_type & bmin,
02734      const double_index_type & bmax,const char * s);
02735 
02736    // Added by Steve Martell, Jan 18, 2014.
02737    void allocate(const data_matrix &m, const char *s);
02738 
02739    int allocated(void) { return (v!=NULL); }
02740    int indexmin(void) {return (index_min);}
02741    int indexmax(void) {return (index_max);}
02742    ~param_init_bounded_number_vector();
02743    void set_initial_value(const double_index_type& it);
02744    void deallocate(void);
02745  };
02746   extern int traceflag;
02747   void tracing_message(int traceflag,const char *s,int *pn);
02748   void tracing_message(int traceflag,const char *s,double *pn);
02749   void set_gauss_covariance_matrix(const dll_data_matrix& m);
02750   void set_gauss_covariance_matrix(const dmatrix& m);
02751   void set_covariance_matrix(const dll_data_matrix& m);
02752   void set_covariance_matrix(const dmatrix& m);
02753 
02754   //ostream& operator<<(const ostream&, const param_init_number_vector);
02755   //ostream& operator<<(const ostream&, const param_init_bounded_number_vector);
02756   //ostream& operator<<(const ostream&, const param_init_vector_vector);
02757   //ostream& operator<<(const ostream&, const param_init_bounded_vector_vector);
02758   //ostream& operator<<(const ostream&, const param_init_matrix_vector);
02759   //ostream& operator<<(const ostream&, const param_init_bounded_matrix_vector);
02760 
02765   class vector_kludge : public dvar_vector
02766   {
02767     public:
02768      vector_kludge(const param_init_number_vector &);
02769      vector_kludge(const param_init_bounded_number_vector &);
02770   };
02775   class matrix_kludge : public dvar_matrix
02776   {
02777     public:
02778      matrix_kludge(const param_init_vector_vector &);
02779      matrix_kludge(const param_init_bounded_vector_vector &);
02780   };
02781 
02782 void read_covariance_matrix(const dmatrix& S, int nvar, int& hbf,
02783   dvector& sscale);
02784 
02785 dvector value(const param_init_number_vector& t);
02786 dvector value(const param_init_bounded_number_vector& t);
02787 //dvector value(const param_init_bounded_number_matrix& t);
02788 //dvector value(const param_init_vector_vector& t);
02789 //dvector value(const param_init_bounded_vector_vector& t);
02790 
02791 dvector read_old_scale(int & old_nvar);
02792 
02793 int withinbound(int lb,int n,int ub);
02794 
02795 #if defined(_MSC_VER)
02796 #  if defined(min)
02797 #    undef min
02798 #  endif
02799 #  if defined(max)
02800 #    undef max
02801 #  endif
02802 #endif
02803 
02804 #include "param_init_bounded_number_matrix.h"
02805 
02806 #endif