ADMB Documentation  11.1.2436
 All Classes Files Functions Variables Typedefs Friends Defines
admodel.h
Go to the documentation of this file.
00001 /*
00002  * $Id: admodel.h 2424 2014-09-29 20:53:15Z 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 class named_i4_array : public i4_array, public model_name_tag
00614 {
00615 protected:
00616   named_i4_array(void) : i4_array(), model_name_tag() {;}
00617   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00618     const index_type&,const index_type&,const index_type&,const index_type&,
00619     const char * s);
00620   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,
00621     int rmax,int cmin,int cmax,const char * s);
00622 
00623   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,
00624     int rmax,const char * s);
00625   void allocate(int hhsl,int hhsu,int hsl,int hsu,const char * s);
00626   void allocate(int hhsl,int hhsu,const char * s);
00627   void allocate(const char * s);
00628 
00629   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00630     const index_type&,const index_type&,const char * s);
00631   void allocate(ad_integer,ad_integer,const index_type&,const index_type&,
00632     const char * s);
00633   void allocate(ad_integer,ad_integer,const char * s);
00634 
00635 
00636   named_i4_array& operator=(const i4_array& m);
00637 };
00638 
00643 class named_d5_array : public d5_array, public model_name_tag
00644 {
00645 protected:
00646   named_d5_array(void) : d5_array(), model_name_tag() {;}
00647   void allocate(int l5,int u5,int hhsl,int hhsu,int hsl,int hsu,int rmin,
00648     int rmax,int cmin,int cmax,const char * s);
00649   void allocate(const ad_integer& hhsl,const ad_integer& hhsu,
00650     const index_type& hsl,const index_type& hsu, const index_type& sl,
00651     const index_type& sh,const index_type& nrl,const index_type& nrh,
00652     const index_type& ncl,const index_type& nch,const char * s);
00653 
00654   named_d5_array& operator=(const d5_array& m);
00655 };
00656 
00661 class named_d6_array : public d6_array, public model_name_tag
00662 {
00663 protected:
00664   named_d6_array(void) : d6_array(), model_name_tag() {;}
00665   void allocate(int l6,int u6,int l5,int u5,int hhsl,int hhsu,int hsl,
00666     int hsu,int rmin,int rmax,int cmin,int cmax,const char * s);
00667   void allocate(const ad_integer& l6,const ad_integer& u6,
00668     const index_type& l5,const index_type& u5,
00669     const index_type& hhsl,const index_type& hhsu,
00670     const index_type& hsl,const index_type& hsu,
00671     const index_type& sl,const index_type& sh,
00672     const index_type& nrl,const index_type& nrh,
00673     const char * s);
00674 
00675   named_d6_array& operator=(const d6_array& m);
00676 };
00677 
00682 class named_d7_array : public d7_array, public model_name_tag
00683 {
00684 protected:
00685   named_d7_array(void) : d7_array(), model_name_tag() {;}
00686   void allocate(int l7,int u7,int l6,int u6,int l5,int u5,int hhsl,
00687     int hhsu,int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
00688     const char * s);
00689   void allocate(const ad_integer& l7,const ad_integer& u7,
00690     const index_type& l6,const index_type& u6,
00691     const index_type& l5,const index_type& u5,
00692     const index_type& hhsl,const index_type& hhsu,
00693     const index_type& hsl,const index_type& hsu,
00694     const index_type& sl,const index_type& sh,
00695     const index_type& nrl,const index_type& nrh,
00696     const char * s);
00697 
00698   named_d7_array& operator=(const d7_array& m);
00699 };
00700 
00701 
00702 class function_minimizer;
00703 
00704 #if defined(USE_ADPVM)
00705 
00709 class pvm_params
00710 {
00711   static pvm_params * varsptr[]; // this should be a resizeable array
00712   static int num_pvm_params;
00713   static const int maxnum_pvm_params;
00714   void add_to_list(void);
00715   virtual void send_to_slaves(void)=0;
00716   virtual void get_from_master(void)=0;
00717 public:
00718   static void pvm_params::send_all_to_slaves(void);
00719   static void pvm_params::get_all_from_master(void);
00720   void allocate(const char *);
00721   void allocate(void);
00722 };
00723 
00728 class pvm_number : public pvm_params
00729 {
00730 public:
00731   virtual void send_to_slaves(void);
00732   virtual void get_from_master(void);
00733   dvector v;
00734   double d;
00735   operator double();
00736   void assign(const dvector&);
00737   void assign(double);
00738 };
00739 
00744 class pvm_int : public pvm_params
00745 {
00746 public:
00747   virtual void send_to_slaves(void);
00748   virtual void get_from_master(void);
00749   ivector v;
00750   int d;
00751   operator int();
00752   void assign(const ivector&);
00753   void assign(int);
00754 };
00755 #endif // #if defined(USE_ADPVM)
00756 
00757 class initial_params;
00758 typedef initial_params* pinitial_params;
00759 typedef void* ptovoid;
00763 class adlist_ptr
00764 {
00765   ptovoid * ptr;
00766   int current_size;
00767   int current;
00768   void resize(void);
00769   void add_to_list(void* p);
00770 public:
00771   adlist_ptr(int init_size);
00772   ~adlist_ptr();
00773 
00774   void initialize();
00775 
00776   pinitial_params& operator[](int i);
00777 
00778   friend class initial_params;
00779 };
00780 
00781 #if defined(USE_SHARE_FLAGS)
00782 
00787   class shareinfo
00788   {
00789     index_type * sflags;
00790     index_type * original_sflags;
00791     index_type * aflags;
00792     index_type * invflags;
00793     i3_array * bmap;
00794     int dimension;
00795     int maxshare;
00796     int current_phase;
00797   public:
00798     void get_inv_matrix_shared( int _cf);
00799     void get_inv_vector_shared( int _cf);
00800     int & get_maxshare(void) { return maxshare; }
00801     i3_array & get_bmap(void);
00802     int & get_dimension(void){ return dimension;}
00803     index_type * get_original_shareflags(void);
00804     index_type * get_shareflags(void);
00805     int& get_current_phase(void);
00806     index_type * get_activeflags(void);
00807     index_type * get_invflags(void);
00808     void set_shareflags(const index_type& sf);
00809     void set_original_shareflags(const index_type& sf);
00810     void reset_shareflags(const index_type& sf);
00811     void set_activeflags(const index_type& af);
00812     void set_bmap(const i3_array& af);
00813     void reset_bmap(const i3_array& af);
00814     void set_invflags(const index_type& af);
00815     shareinfo(const index_type& sf,const index_type& af);
00816     ~shareinfo();
00817   };
00818 #endif
00819 
00824 class initial_params
00825 {
00826 protected:
00827 #if defined(USE_SHARE_FLAGS)
00828   shareinfo * share_flags;
00829 #endif
00830   virtual ~initial_params();
00831   int active_flag;
00832   int initial_value_flag;
00833   double initial_value;
00834   double scalefactor;
00835 public:
00836 #if defined(USE_SHARE_FLAGS)
00837   virtual void setshare(const index_type& sf,const index_type& af);
00838   virtual void shared_set_value_inv(const dvector&,const int&);
00839   virtual void shared_set_value(const dvar_vector&,const int&,
00840     const dvariable& pen);
00841   virtual int shared_size_count(void); // get the number of active parameters
00842 #endif
00843   double get_scalefactor();
00844   void set_scalefactor(const double);
00845 #if !defined(BIG_INIT_PARAMS)
00846   static initial_params  varsptr[]; // this should be a resizeable array
00847 #else
00848   static adlist_ptr varsptr; // this should be a resizeable array
00849 #endif
00850   static int num_initial_params;
00851   static const int max_num_initial_params;
00852   static int straight_through_flag;
00853   static int num_active_initial_params;
00854   static int max_number_phases;
00855   static int current_phase;
00856   static int restart_phase;
00857   static int sd_phase;
00858   static int mc_phase;
00859   static int mceval_phase;
00860   int phase_start;
00861   int phase_save;
00862   int phase_stop;
00863   virtual void set_random_effects_active();
00864   void restore_phase_start(void);
00865   virtual void set_random_effects_inactive();
00866   virtual void set_only_random_effects_active();
00867   virtual void set_only_random_effects_inactive();
00868   virtual void set_value(const dvar_vector&, const int&,
00869     const dvariable& pen) = 0;
00870   virtual void dev_correction(const dmatrix&, const int&) = 0;
00871   void set_initial_value(double x);
00872   double get_initial_value(void);
00873   void set_phase_start(int x);
00874   int get_phase_start(void);
00875   static void set_all_simulation_bounds(const dmatrix& symbds);
00876   static void set_all_simulation_bounds(const dmatrix& symbds, const dvector&);
00877   static void get_jacobian_value(const dvector& y, const dvector& jac);
00878   static int correct_for_dev_objects(const dmatrix &H);
00879   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii) = 0;
00880   //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
00881   //  const int& ii) = 0;
00882   virtual void set_value_inv(const dvector&, const int&) = 0;
00883   virtual void add_value(const dvector&, const int&) = 0;
00884   virtual void add_value(const dvector&, const dvector&, const int&,
00885     const double&, const dvector&) = 0;
00886   virtual void get_jacobian(const dvector&, const dvector&, const int&) = 0;
00887   //virtual void check_tightness(const dvector&, const int&) = 0;
00888   virtual void copy_value_to_vector(const dvector&, const int&) = 0;
00889   virtual void restore_value_from_vector(const dvector&, const int&) = 0;
00890   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii) = 0;
00891   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
00892     const int& ii)=0;
00893   virtual void mc_scale(const dvector& d, const dvector& x, const int& ii) = 0;
00894   virtual void curv_scale(const dvector& d, const dvector& x,
00895     const int& ii) = 0;
00896   virtual void hess_scale(const dvector& d, const dvector& x,
00897     const int& ii) = 0;
00898   virtual int size_count(void)=0; // get the number of active parameters
00899   virtual void save_value(void)=0; // save the objects value in a ascii file
00900   virtual void bsave_value(void)=0; // save the objects value in a binary file
00901   virtual void save_value(const ofstream& ofs, int prec) = 0;
00902   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
00903     int& offset){}
00904   //virtual void bsave_value(const uostream& ofs) = 0;
00905     virtual const char * label()=0;
00906   void allocate(int _phase_start);
00907   void set_active_flag(void);
00908   void set_inactive_flag(void);
00909   friend int active(const initial_params& ip);
00910   static adstring get_reportfile_name(void);
00911   initial_params(void);
00912   static void xinit(const dvector& x); // get the number of active parameters
00913   // get the number of active parameters
00914   static void xinit_all(const dvector& x);
00915   static void save_all(const ofstream& _ofs,int prec,const dvector&g);
00916   // get the number of active parameters
00917   static void set_active_random_effects(void);
00918   // get the number of active parameters
00919   static void set_active_only_random_effects(void);
00920   // get the number of active parameters
00921   static void set_inactive_only_random_effects(void);
00922   // get the number of active parameters
00923   static void set_inactive_random_effects(void);
00924   // get the number of active parameters
00925   static void restore_start_phase(void);
00926   static void xinit1(const dvector& x, const dvector& g);
00927   //save all initial parameter values in a vector
00928   static void copy_all_values(const dvector& x, const int& ii);
00929   //get ivalues for all active parameters from a vector
00930   static void restore_all_values(const dvector& x, const int& ii);
00931   // get the number of active parameters
00932   static dvariable reset(const dvar_vector& x);
00933   // get the number of active parameters
00934   static dvariable reset(const dvector& x);
00935   static dvariable reset1(const dvar_vector& x, const dvector& g);
00936   // get the number of active parameters
00937   static dvariable reset(const dvar_vector& x, const dvector& pen);
00938   // get the number of active parameters
00939   static dvariable reset_all(const dvar_vector& x,const dvector& pen);
00940   static size_t nvarcalc();
00941   static int nvarcalc_all(void);
00942   static int num_active_calc(void);
00943   static int stddev_scale(const dvector& d, const dvector& x);
00944   static int stddev_vscale(const dvar_vector& d,const dvar_vector& x);
00945   static int montecarlo_scale(const dvector& d, const dvector& x);
00946   static int stddev_curvscale(const dvector& d, const dvector& x);
00947   static void read(void);
00948   static void save(void);
00949   static void save(const ofstream& ofs, int prec);
00950   static void restore(const ifstream& ifs);
00951   static void add_random_vector(const dvector& x);
00952   static void add_random_vector(const dvector& y, const dvector& x,
00953     const double& ll, const dvector& diag);
00954   virtual void restore_value(const ifstream& ifs) = 0;
00955   virtual void add_to_list(void);
00956 #if defined(USE_ADPVM)
00957   virtual void pvm_pack(void)=0;
00958   virtual void pvm_unpack(void)=0;
00959 #endif
00960 
00961   friend class function_minimizer;
00962 };
00963 
00964 void pvm_pack(const dvar_vector&);
00965 void pvm_unpack(const dvar_vector&);
00966 void pvm_pack(const prevariable&);
00967 void pvm_unpack(const prevariable&);
00968 void pvm_pack(const dvar_matrix&);
00969 void pvm_unpack(const dvar_matrix&);
00970 
00975 class param_init_vector: public named_dvar_vector , public initial_params
00976 {
00977 public:
00978 #if defined(USE_SHARE_FLAGS)
00979   virtual void setshare(const index_type& sf,const index_type& af);
00980   virtual void shared_set_value_inv(const dvector&,const int&);
00981   virtual void shared_set_value(const dvar_vector&,const int&,
00982     const dvariable& pen);
00983   virtual int shared_size_count(void); // get the number of active parameters
00984 #endif
00985   virtual void dev_correction(const dmatrix& H, const int& ii);
00986   virtual const char * label(void);
00987   param_init_vector();
00988 #if defined(USE_ADPVM)
00989   void pvm_pack(void){::pvm_pack(*this);}
00990   void pvm_unpack(void){::pvm_unpack(*this);}
00991 #endif
00992 private:
00993   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
00994   //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
00995   //  const int& ii);
00996   virtual void add_value(const dvector&, const dvector&, const int&,
00997     const double&, const dvector&);
00998   virtual void add_value(const dvector&, const int&);
00999   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01000   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01001     const int& ii);
01002 
01003   virtual void set_value(const dvar_vector& x, const int& ii,
01004     const dvariable& pen);
01005   virtual void set_value_inv(const dvector& x, const int& ii);
01006   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01007   virtual void restore_value_from_vector(const dvector&, const int&);
01008   virtual int size_count(void);
01009   virtual void sd_scale(const dvector& d, const dvector& x, const int&);
01010   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01011   virtual void curv_scale(const dvector& d, const dvector& x,const int&);
01012   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01013   virtual void save_value(void);
01014   virtual void bsave_value(void);
01015   virtual void save_value(const ofstream& ofs, int prec);
01016   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
01017     int& offset);
01018   virtual void restore_value(const ifstream& ifs);
01019   void report_value(void);
01020   //virtual void read_value(void);
01021   void allocate(int imin,int imax,int phasestart=1,const char * s="UNNAMED");
01022   void allocate(const ad_integer& imin,const ad_integer& imax,
01023     const ad_integer& phasestart=1,const char * s="UNNAMED");
01024   void allocate(int imin,int imax,const char * s="UNNAMED");
01025 #if defined(USE_SHARE_FLAGS)
01026   void allocate(int imin,int imax,const ivector& ishare,
01027     const char * s="UNNAMED");
01028 #endif
01029   friend class model_parameters;
01030   friend class param_init_vector_vector;
01031 public:
01032   param_init_vector& operator = (const dvector&);
01033   param_init_vector& operator = (const dvar_vector&);
01034   param_init_vector& operator = (const prevariable&);
01035   param_init_vector& operator = (const double&);
01036 };
01037 
01042 class dll_param_init_vector: public param_init_vector
01043 {
01044   double * pd;
01045 public:
01046   dll_param_init_vector& operator = (const dvector&);
01047   dll_param_init_vector& operator = (const dvar_vector&);
01048   dll_param_init_vector& operator = (const prevariable&);
01049   dll_param_init_vector& operator = (const double&);
01050   void allocate(double * _pd,int imin,int imax,
01051     int phasestart=1,const char * s="UNNAMED");
01052   void allocate(double * _pd,int imin,int imax,
01053     const char * s="UNNAMED");
01054 
01055   virtual ~dll_param_init_vector();
01056   friend class model_parameters;
01057 };
01058 
01063 class param_init_bounded_vector: public named_dvar_vector,public initial_params
01064 {
01065   virtual void* parent_this(void){return this;}
01066 public:
01067   double get_minb(void);
01068   void set_minb(double b);
01069   double get_maxb(void);
01070   void set_maxb(double b);
01071 protected:
01072   double minb;
01073   double maxb;
01074   param_init_bounded_vector();
01075 private:
01076   virtual void dev_correction(const dmatrix&, const int&);
01077   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01078   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01079     const int& ii);
01080   //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01081   //  const int& ii);
01082   virtual void add_value(const dvector&, const dvector&, const int&,
01083     const double&, const dvector&);
01084   virtual void add_value(const dvector&, const int&);
01085   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01086   virtual void set_value(const dvar_vector& x, const int& ii,
01087     const dvariable& pen);
01088   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01089   virtual void restore_value_from_vector(const dvector&, const int&);
01090   virtual void set_value_inv(const dvector& x, const int& ii);
01091   virtual int size_count(void);
01092   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01093   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01094   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01095   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01096   void allocate(int imin,int imax,double _minb,double _maxb,
01097     int phasestart=1, const char * name="UNNAMED");
01098   void allocate(int imin,int imax,double _minb,double _maxb,
01099     const char * name="UNNAMED");
01100  //void param_init_bounded_vector::allocate(const ad_integer& imin,
01101    //const ad_integer& imax,const ad_double& _minb,const ad_double& _maxb,
01102    //const ad_integer& phase_start,const char * s);
01103   friend class model_parameters;
01104   friend class param_init_bounded_vector_vector;
01105   virtual const char * label(void);
01106   virtual void save_value(const ofstream& ofs, int prec);
01107   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
01108     int& offset);
01109   virtual void restore_value(const ifstream& ifs);
01110   virtual void save_value(void);
01111   virtual void bsave_value(void);
01112   void report_value(void);
01113   //virtual void read_value(void);
01114 public:
01115   param_init_bounded_vector& operator = (const dvector&);
01116   param_init_bounded_vector& operator = (const dvar_vector&);
01117   param_init_bounded_vector& operator = (const prevariable&);
01118   param_init_bounded_vector& operator = (const double&);
01119 #if defined(USE_ADPVM)
01120   void pvm_pack(void){::pvm_pack(*this);}
01121   void pvm_unpack(void){::pvm_unpack(*this);}
01122 #endif
01123 };
01124 
01129 class dll_param_init_bounded_vector: public param_init_bounded_vector
01130 {
01131   double * pd;
01132 public:
01133   void allocate(double * _pd,int imin,int imax,double _minb,double _maxb,
01134     int phasestart=1, const char * name="UNNAMED");
01135   void allocate(double * _pd,int imin,int imax,double _minb,double _maxb,
01136     const char * name="UNNAMED");
01137   ~dll_param_init_bounded_vector();
01138 };
01139 
01144 class param_init_bounded_dev_vector: public param_init_bounded_vector
01145 {
01146   virtual void set_value(const dvar_vector& x, const int& ii,
01147     const dvariable& pen);
01148   virtual void dev_correction(const dmatrix& H, const int& ii);
01149 public:
01150   param_init_bounded_dev_vector& operator = (const dvar_vector& m);
01151   param_init_bounded_dev_vector& operator = (const dvector& m);
01152   param_init_bounded_dev_vector& operator = (const prevariable& m);
01153   param_init_bounded_dev_vector& operator = (const double& m);
01154 };
01155 
01160 class param_init_number: public named_dvariable , public initial_params
01161 {
01162   virtual void dev_correction(const dmatrix&, const int&);
01163   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01164 
01165 #if defined(USE_ADPVM)
01166   void pvm_pack(void){::pvm_pack(*this);}
01167   void pvm_unpack(void){::pvm_unpack(*this);}
01168 #endif
01169 
01170 //  virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01171 //    const int& ii);
01172   virtual void add_value(const dvector&, const dvector&, const int&,
01173     const double&, const dvector&);
01174   virtual void add_value(const dvector&, const int&);
01175   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01176   virtual void set_value(const dvar_vector& x, const int& ii,
01177     const dvariable& pen);
01178   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01179   virtual void restore_value_from_vector(const dvector&, const int&);
01180   virtual void set_value_inv(const dvector& x, const int& ii);
01181   virtual int size_count(void);
01182   virtual void save_value(const ofstream& ofs, int prec);
01183   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
01184     int& offset);
01185   virtual void restore_value(const ifstream& ifs);
01186   virtual void save_value(void);
01187   virtual void bsave_value(void);
01188   void report_value(void);
01189   virtual const char * label(void);
01190   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01191   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01192   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01193   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01194   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01195     const int& ii);
01196   //virtual void read_value(void);
01197 protected:
01198   void allocate(int phase_start=1,const char *s="UNNAMED");
01199   void allocate(const char *s="UNNAMED");
01200   void allocate(init_xml_doc&, const char *s="UNNAMED");
01201   friend class model_parameters;
01202   friend class param_init_number_vector;
01203   param_init_number();
01204   param_init_number& operator=(const double m);
01205   param_init_number& operator=(const prevariable& m);
01206 };
01207 
01212 class dll_param_init_number: public param_init_number
01213 {
01214   double * pd;
01215 public:
01216   void allocate(double * pd,int phase_start=1,const char *s="UNNAMED");
01217   void allocate(double *pd,const char *s="UNNAMED");
01218   virtual ~dll_param_init_number();
01219   dll_param_init_number& operator=(const double m);
01220   dll_param_init_number& operator=(const prevariable& m);
01221 };
01222 
01223 //forward declaration
01224 class data_vector;
01231 class param_init_bounded_number: public param_init_number
01232 {
01233 public:
01234   double get_minb(void);
01235   void set_minb(double b);
01236   double get_maxb(void);
01237   void set_maxb(double b);
01238 protected:
01239   double minb;
01240   double maxb;
01241   void allocate(double _minb,double _maxb,int phase_start=1,
01242     const char * s="UNNAMED");
01243   void allocate(double _minb,double _maxb,const char * s="UNNAMED");
01244   // Added by Steve Martell for using input data for allocation.
01245   void allocate(const data_vector& v,const char * s="UNNAMED");
01246   void allocate(init_xml_doc&, const char * s="UNNAMED");
01247 
01248 public:
01249 #if defined(USE_ADPVM)
01250   void pvm_pack(void){::pvm_pack(*this);}
01251   void pvm_unpack(void){::pvm_unpack(*this);}
01252 #endif
01253   param_init_bounded_number();
01254 private:
01255   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01256 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01257 //  const int& ii);
01258   virtual void add_value(const dvector&, const dvector&, const int&,
01259     const double&, const dvector&);
01260   virtual void add_value(const dvector&, const int&);
01261   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01262   virtual void set_value(const dvar_vector& x, const int& ii,
01263     const dvariable& pen);
01264   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01265   virtual void restore_value_from_vector(const dvector&, const int&);
01266   virtual void set_value_inv(const dvector& x, const int& ii);
01267   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01268   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01269   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01270   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01271   virtual const char * label(void);
01272   void report_value(void);
01273   param_init_bounded_number& operator=(const double m);
01274   param_init_bounded_number& operator=(const prevariable& m);
01275   friend class model_parameters;
01276   friend class param_init_bounded_number_vector;
01277   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01278     const int& ii);
01279 };
01280 
01285 class dll_param_init_bounded_number: public param_init_bounded_number
01286 {
01287   double * pd;
01288 public:
01289   void allocate(double * _pd,double _minb,double _maxb,int phase_start=1,
01290     const char * s="UNNAMED");
01291   void allocate(double* _pd,double _minb, double _maxb,const char* s="UNNAMED");
01292 public:
01293   virtual ~dll_param_init_bounded_number();
01294   void report_value(void);
01295 };
01296 
01301 class param_init_matrix: public named_dvar_matrix,public initial_params
01302 {
01303 #if defined(USE_SHARE_FLAGS)
01304   virtual void shared_set_value_inv(const dvector& x,const int& ii);
01305   virtual int shared_size_count(void); // get the number of active parameters
01306   virtual void shared_set_value(const dvar_vector&,const int&,
01307     const dvariable& pen);
01308 #endif
01309   virtual void dev_correction(const dmatrix&, const int&);
01310   virtual void set_simulation_bounds(const dmatrix& symbds,const int& ii);
01311 #if defined(USE_ADPVM)
01312   void pvm_pack(void){::pvm_pack(*this);}
01313   void pvm_unpack(void){::pvm_unpack(*this);}
01314 #endif
01315 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01316 //  const int& ii);
01317   virtual void add_value(const dvector&, const dvector&, const int&,
01318     const double&, const dvector&);
01319   virtual void add_value(const dvector&, const int&);
01320   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01321 public:
01322   virtual void set_value(const dvar_vector& x, const int& ii,
01323     const dvariable& pen);
01324   virtual void copy_value_to_vector(const dvector& x, const int& ii);
01325   virtual void restore_value_from_vector(const dvector&, const int&);
01326   virtual void set_value_inv(const dvector& x, const int& ii);
01327   virtual int size_count(void);
01328   virtual void save_value(void);
01329   virtual void save_value(const ofstream& ofs, int prec,const dvector&,
01330     int& offset);
01331   virtual void bsave_value(void);
01332   virtual void save_value(const ofstream& ofs, int prec);
01333   virtual void restore_value(const ifstream& ifs);
01334   void report_value(void);
01335   //virtual void read_value(void);
01336   virtual const char * label(void);
01337   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01338   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01339   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01340   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01341   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01342     const int& ii);
01343 
01344 public:
01345   void allocate(int rmin,int rmax,int cmin,int cmax,
01346     int phase_start=1,const char * = "UNNAMED");
01347 #if defined(USE_SHARE_FLAGS)
01348   //void allocate(int rmin,int rmax,int cmin,int cmax, const imatrix& jshare,
01349   // const char * = "UNNAMED");
01350   virtual void setshare(const index_type& sf,const index_type& af);
01351 #endif
01352   void allocate(int rmin,int rmax,int cmin,int cmax,
01353     const char * = "UNNAMED");
01354  void allocate(const ad_integer& imin,const ad_integer&imax,
01355    const index_type& imin2,const index_type& imax2,
01356    const ad_integer& phase_start, const char * s);
01357   void allocate(const ad_integer& rmin, const ad_integer& rmax,
01358                 const index_type& cmin, const index_type& cmax,
01359                 const char * = "UNNAMED");
01360   void allocate(const ad_integer& rmin, const ad_integer& rmax,
01361                 const index_type& cmin, const index_type& cmax,
01362                 int phase_start = 1, const char * = "UNNAMED");
01363   param_init_matrix(void);
01364   param_init_matrix& operator = (const dmatrix& m);
01365   param_init_matrix& operator = (const dvar_matrix& m);
01366   param_init_matrix& operator = (const dvariable& m);
01367   param_init_matrix& operator = (const double& m);
01368 };
01369 
01374 class dll_param_init_matrix: public param_init_matrix
01375 {
01376   double * d;
01377 public:
01378   void allocate(double* _d,int rmin,int rmax,int cmin,int cmax,
01379     int phase_start=1,const char * = "UNNAMED");
01380   void allocate(double * _d,int rmin,int rmax,int cmin,int cmax,
01381     const char * = "UNNAMED");
01382   virtual ~dll_param_init_matrix();
01383   dll_param_init_matrix(){d=NULL;}
01384   dll_param_init_matrix& operator = (const dmatrix& m);
01385   dll_param_init_matrix& operator = (const dvar_matrix& m);
01386   dll_param_init_matrix& operator = (const dvariable& m);
01387   dll_param_init_matrix& operator = (const double& m);
01388 };
01389 
01394 class param_init_bounded_matrix: public param_init_matrix
01395 {
01396 public:
01397   double get_minb(void);
01398   void set_minb(double b);
01399   double get_maxb(void);
01400   void set_maxb(double b);
01401 protected:
01402   double minb;
01403   double maxb;
01404   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
01405 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
01406 //  const int& ii);
01407   virtual void add_value(const dvector&, const dvector&, const int&,
01408     const double&, const dvector&);
01409   virtual void add_value(const dvector&, const int&);
01410   virtual void get_jacobian(const dvector&, const dvector&, const int&);
01411 public:
01412 #if defined(USE_ADPVM)
01413   void pvm_pack(void){::pvm_pack(*this);}
01414   void pvm_unpack(void){::pvm_unpack(*this);}
01415 #endif
01416   virtual void set_value(const dvar_vector& x, const int& ii,
01417     const dvariable& pen);
01418 #if defined(USE_SHARE_FLAGS)
01419   virtual void shared_set_value_inv(const dvector&,const int&);
01420   virtual void shared_set_value(const dvar_vector&,const int&,
01421     const dvariable& pen);
01422 #endif
01423   virtual void set_value_inv(const dvector& x, const int& ii);
01424   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
01425   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
01426   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
01427   virtual void hess_scale(const dvector&, const dvector&, const int&){}
01428   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
01429     const int& ii);
01430 
01431 public:
01432  void allocate(const ad_integer& imin,
01433    const ad_integer& imax, const ad_integer& imin2,
01434    const ad_integer& imax2, const ad_double& _bmin,
01435    const ad_double& _bmax, const ad_integer& phase_start,
01436    const char * s);
01437 
01438   param_init_bounded_matrix(void);
01439   void allocate(int rmin,int rmax,int cmin,int cmax,
01440     double _minb,double _maxb,
01441     int phase_start=1,const char * = "UNNAMED");
01442   void allocate(int rmin,int rmax,int cmin,int cmax,
01443     double _minb,double _maxb,const char * = "UNNAMED");
01444 
01445   void allocate(const ad_integer& rmin, const ad_integer& rmax,
01446     const index_type& cmin,
01447                 const index_type& cmax, double _minb, double _maxb,
01448                 const char * = "UNNAMED");
01449   void allocate(const ad_integer& rmin, const ad_integer& rmax,
01450     const index_type& cmin,
01451                 const index_type& cmax, double _minb, double _maxb,
01452                 int phase_start = 1, const char * = "UNNAMED");
01453 };
01454 
01459 class data_int : public model_name_tag
01460 {
01461 protected:
01462   int val;
01463   data_int& operator=(const int);
01464   void allocate(int n,const char * s="UNNAMED");
01465   void allocate(const char * s="UNNAMED");
01466   void allocate(init_xml_doc&, const char * s="UNNAMED");
01467   friend class model_data;
01468   friend class model_parameters;
01469   friend int operator + (int n,data_int v);
01470   friend int operator + (data_int v,int n);
01471   friend int operator + (data_int v,data_int n);
01472 public:
01473   operator int() {return val;}
01474 #if !defined (__BORLANDC__)
01475   //operator const int() const {return val;}
01476 #endif
01477   virtual ~data_int(){;}
01478 };
01479 
01484 class named_adstring : public adstring, public model_name_tag
01485 {
01486 protected:
01487   void allocate(const char * s1,const char * s="UNNAMED");
01488   void operator = (const adstring&);
01489   void operator = (const char *);
01490 };
01491 
01496 class named_line_adstring : public line_adstring, public model_name_tag
01497 {
01498 protected:
01499   void allocate(const char * s1,const char * s="UNNAMED");
01500   void operator = (const adstring&);
01501   void operator = (const char *);
01502 };
01503 
01508 class init_adstring: public named_adstring
01509 {
01510 public:
01511   void allocate(const char * s="UNNAMED");
01512 };
01513 
01518 class init_line_adstring: public named_line_adstring
01519 {
01520 public:
01521   void allocate(const char * s="UNNAMED");
01522 };
01523 
01528 class dll_named_adstring : public named_adstring
01529 {
01530   char ** d;
01531 public:
01532   void allocate(char ** ps1,const char * s="UNNAMED");
01533   void operator = (const adstring&);
01534   void operator = (const char *);
01535   ~dll_named_adstring();
01536   dll_named_adstring(void){d=NULL;}
01537 };
01538 
01543 class dll_data_int : public data_int
01544 {
01545 public:
01546   int *pi;
01547   void allocate(int *_pi,const char * s);
01548   virtual ~dll_data_int();
01549 };
01550 
01555 class data_matrix : public named_dmatrix
01556 {
01557 public:
01558   data_matrix(void) : named_dmatrix() {;}
01559   data_matrix& operator=(const dmatrix& m);
01560   data_matrix& operator=(const double& m);
01561 private:
01562   void allocate(int rmin,int rmax,int cmin,int cmax,
01563     const char * = "UNNAMED");
01564   void allocate(int rmin, int rmax, const ivector& cmin, const ivector& cmax,
01565     const char * = "UNNAMED");
01566   void allocate(int rmin, int rmax, const ivector& cmin, int cmax,
01567     const char * = "UNNAMED");
01568   void allocate(int rmin, int rmax, int cmin, const ivector& cmax,
01569     const char * = "UNNAMED");
01570   void allocate(init_xml_doc&, const char * = "UNNAMED");
01571   friend class model_data;
01572 };
01573 
01578 class dll_data_matrix : public data_matrix
01579 {
01580   double * d;
01581 public:
01582   void allocate(double * _d,int rmin,int rmax,int cmin,int cmax,
01583     const char * _s = "UNNAMED");
01584   virtual ~dll_data_matrix();
01585   dll_data_matrix& operator=(const dmatrix &m);
01586   dll_data_matrix& operator=(const double &m);
01587 };
01588 
01593 class data_3array : public named_d3_array
01594 {
01595 public:
01596   data_3array(void) : named_d3_array() {;}
01597 private:
01598   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
01599     const char * ="UNNAMED");
01600   void allocate(int hsl, int hsu, const ivector& rmin, int rmax,
01601     int cmin,int cmax,const char * ="UNNAMED");
01602   void allocate(int hsl, int hsu, int rmin, const ivector& rmax,
01603     int cmin,int cmax,const char * ="UNNAMED");
01604   void allocate(int hsl, int hsu, const ivector& rmin, const ivector& rmax,
01605     int cmin,int cmax,const char * ="UNNAMED");
01606   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
01607     int cmax, const char * ="UNNAMED");
01608   void allocate(int hsl, int hsu, int rmin, int rmax, const ivector& cmin,
01609                 const ivector& cmax, const char * ="UNNAMED");
01610   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,
01611                 const ivector& cmax, const char * ="UNNAMED");
01612   void allocate(int hsl,int hsu, const index_type& rmin, const index_type& rmax,
01613     const index_type& cmin, const index_type& cmax, const char * ="UNNAMED");
01614   friend class model_data;
01615 };
01616 
01617 class data_3iarray : public named_i3_array
01618 {
01619   data_3iarray(void) : named_i3_array() {;}
01620   void allocate(int hsl,int hsu,int rmin,int rmax,int cmin,int cmax,
01621     const char * ="UNNAMED");
01622   void allocate(int hsl, int hsu,const index_type& rmin, const index_type& rmax,
01623     const index_type& cmin, const index_type& cmax, const char * ="UNNAMED");
01624   friend class model_data;
01625 };
01626 
01631 class dll_data_3array : public data_3array
01632 {
01633   double * d;
01634 public:
01635   void allocate(double * _d,int hmin,int hmax,int rmin,int rmax,
01636     int cmin,int cmax,const char * _s = "UNNAMED");
01637   dll_data_3array& operator=(const d3_array &);
01638   virtual ~dll_data_3array();
01639   friend class model_data;
01640 };
01641 
01642 
01647 class data_5array : public named_d5_array
01648 {
01649   data_5array(void) : named_d5_array() {;}
01650   void allocate(int hhsl,int hhsu,
01651     int hhhsl,int hhhsu,
01652     int hsl,int hsu,int rmin,int rmax,
01653     int cmin,int cmax,const char * ="UNNAMED");
01654   void allocate(ad_integer hhhsl, ad_integer hhhsu, const index_type& hhsl,
01655     const index_type& hhsu, const index_type& hsl, const index_type& hsu,
01656     const index_type& rmin, const index_type& rmax, const index_type& cmin,
01657     const index_type& cmax, const char * ="UNNAMED");
01658   friend class model_data;
01659 };
01660 
01665 class data_4array : public named_d4_array
01666 {
01667   data_4array(void) : named_d4_array() {;}
01668   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,int rmax,
01669     int cmin,int cmax,const char * ="UNNAMED");
01670   void allocate(ad_integer hhsl, ad_integer hhsu, const index_type& hsl,
01671     const index_type& hsu, const index_type& rmin, const index_type& rmax,
01672     const index_type& cmin, const index_type& cmax, const char * ="UNNAMED");
01673   friend class model_data;
01674 };
01675 
01676 class data_4iarray : public named_i4_array
01677 {
01678   data_4iarray(void) : named_i4_array() {;}
01679   void allocate(int hhsl,int hhsu,int hsl,int hsu,int rmin,int rmax,
01680     int cmin,int cmax,const char * ="UNNAMED");
01681   void allocate(ad_integer hhsl, ad_integer hhsu, const index_type& hsl,
01682     const index_type& hsu, const index_type& rmin, const index_type& rmax,
01683     const index_type& cmin, const index_type& cmax, const char * ="UNNAMED");
01684   friend class model_data;
01685 };
01686 
01691 class data_imatrix : public named_imatrix
01692 {
01693   data_imatrix(void) : named_imatrix() {;}
01694   void allocate(int rmin,int rmax,int cmin,int cmax,const char * ="UNNAMED");
01695   void allocate(int rmin, int rmax, const index_type&, const index_type& cmax,
01696     const char * ="UNNAMED");
01697   friend class model_data;
01698 };
01699 
01704 class data_vector : public named_dvector
01705 {
01706 public:
01707   data_vector& operator=(const dvector& m);
01708   data_vector& operator=(const double m);
01709   data_vector(void) : named_dvector() {;}
01710 private:
01711   void allocate(int imin,int imax,const char * ="UNNAMED");
01712   void allocate(int imin, const ivector& imax, const char * ="UNNAMED");
01713   void allocate(init_xml_doc&, const char * ="UNNAMED");
01714   friend class model_data;
01715 };
01716 
01721 class dll_data_vector : public data_vector
01722 {
01723 public:
01724   double * pd;
01725   void allocate(double * pd,int imin,int imax,const char * ="UNNAMED");
01726   void allocate(double *pd, int imin, const ivector& imax,
01727     const char * ="UNNAMED");
01728   virtual ~dll_data_vector();
01729   dll_data_vector& operator = (const dvector& x);
01730   dll_data_vector& operator = (const double& x);
01731 };
01732 
01737 class data_ivector : public named_ivector
01738 {
01739 public:
01740   data_ivector(void) : named_ivector() {;}
01741 private:
01742   void allocate(int imin,int imax,const char * ="UNNAMED");
01743   friend class model_data;
01744 };
01745 
01750 class data_number : public model_name_tag
01751 {
01752 protected:
01753   double val;
01754   void allocate(const char * ="UNNAMED");
01755 public:
01756   void report_value(void);
01757   operator double() {return val;}
01758   double& value(void) {return val;}
01759   void initialize(void) {val=0.0;}
01760   friend class model_data;
01761   data_number & operator=(const double& m);
01762 };
01763 
01768 class dll_data_number : public data_number
01769 {
01770 public:
01771   double * pd;
01772   void allocate(double *_pd,const char * s);
01773   virtual ~dll_data_number();
01774   dll_data_number & operator=(const double& m);
01775 };
01776 
01777 typedef dvariable (model_parameters::*PMF) (const dvariable&);
01778 typedef dvariable (model_parameters::*PMFI) (const dvariable&,int n);
01779 typedef dvariable (model_parameters::*PMFVI) (const dvar_vector&,int n);
01780 typedef void (model_parameters::*PMFVIV4) (const dvar_vector&,int n,
01781   dvariable& f1, const dvariable& f2, const dvariable& f3, const dvariable& f4);
01782 
01783   class init_df1b2vector;
01784   class df1b2vector;
01785   class df1b2variable;
01786 
01791 class function_minimizer
01792 {
01793 public:
01794   static int bad_step_flag;
01795   static int likeprof_flag;
01796   static int first_hessian_flag;
01797   static int test_trust_flag;
01798   static int random_effects_flag;
01799   dmatrix * negdirections;
01800   static int negative_eigenvalue_flag;
01801   static int inner_opt_flag;
01802   static int inner_opt(void);
01803   laplace_approximation_calculator * lapprox;
01804   dvector * multinomial_weights;
01805   void set_multinomial_weights(dvector&d);
01806   //init_df1b2vector* py;
01807   virtual void AD_uf_inner();
01808   virtual void AD_uf_outer();
01809   virtual void user_function();
01810   void pre_user_function(void);
01811   //void df1b2_pre_user_function(void);
01812   //virtual void user_function(const init_df1b2vector& x,df1b2variable& f);
01813   //static int hesstype;
01814   //int set_hessian_type(int);
01815   void hess_routine_noparallel_randomeffects(void);
01816   void other_separable_stuff_begin(void);
01817   void other_separable_stuff_end(void);
01818   void begin_df1b2_funnel(void){;}
01819   void setup_quadprior_calcs(void){;}
01820   void end_df1b2_funnel(void){;}
01821   void get_function_difference(void);
01822   void start_get_importance_sampling_comnponent(void);
01823   void end_get_importance_sampling_comnponent(void);
01824   int spminflag;
01825   int repeatminflag;
01826   int mcmc2_flag;
01827   int robust_hybrid_flag;
01828   long ifn;
01829   int maxfn;
01830   int iprint;
01831   double crit;
01832   int imax;
01833   double dfn;
01834   long iexit;
01835   long ihflag;
01836   long ihang;
01837   long scroll_flag;
01838   int maxfn_flag;
01839   int quit_flag;
01840   double min_improve;
01841   void pre_userfunction(void);
01842   virtual void userfunction(void)=0;
01843   virtual void allocate(void){;}
01844   static named_dvar_vector * ph;
01845   static named_dvar_vector * pg;
01846 protected:
01847   double ffbest;
01848 private:
01849   gradient_structure * pgs;
01850   adstring_array param_labels;
01851   ivector param_size;
01852 protected:
01853   void report_function_minimizer_stats(void){;}
01854   virtual void report(const dvector& gradients){;};
01855   static dvector convergence_criteria;
01856   static dvector maximum_function_evaluations;
01857   static int sd_flag;
01858   static adstring user_data_file;
01859   static adstring user_par_file;
01860   static int have_constraints;
01861 public:
01862   virtual dvariable user_randeff(const dvar_vector& x);
01863   virtual dvar_vector user_dfrandeff(const dvar_vector& x);
01864   virtual dvar_matrix user_d2frandeff(const dvar_vector& x);
01865   void limited_memory_quasi_newton(const independent_variables&,int);
01866   void limited_memory_quasi_newton(const independent_variables&,int,int);
01867   void limited_memory_quasi_newton(double& f, const independent_variables&,
01868     int, int, int,double);
01869   function_minimizer(long int sz=0L);
01870   void likeprof_routine(double global_min);
01871   virtual ~function_minimizer();
01872   virtual void other_calculations(void){;}
01873   virtual void final_calcs(void){;}
01874   virtual void minimize(void);
01875   virtual void constraints_minimize(void);
01876   virtual void between_phases_calculations(void){;}
01877   void computations(int argc,char * argv[]);
01878   void computations1(int argc,char * argv[]);
01879   void computations_np(int argc,char * argv[]);
01880   void computations(void);
01881   void hess_routine(void);
01882   void hess_routine_noparallel(void);
01883   void hess_routine_master(void);
01884   void hess_routine_slave(void);
01885   void constraint_report(void);
01886   void pvm_slave_likeprof_routine(void);
01887   void pvm_master_function_evaluation_profile(double& f,
01888     independent_variables& x,const dvector & g,int nvar,int iprof,double weight,
01889     double new_value);
01890   void pvm_slave_prof_minimize(int underflow_flag);
01891   void pvm_master_prof_minimize(int iprof, double sigma,
01892     double new_value, const double& _fprof, const int underflow_flag,
01893     double global_min, const double& _penalties,
01894     const double& _final_weight);
01895   //void pvm_master_function_evaluation_profile(double& f,
01896    // independent_variables& x,dvector & g,int nvar,int iprof);
01897   void pvm_slave_function_evaluation(void);
01898   void pvm_slave_function_evaluation_no_derivatives(void);
01899   void pvm_slave_function_evaluation_noder(void);
01900   void pvm_master_function_evaluation_no_derivatives(double& f,
01901     independent_variables& x,int nvar);
01902   void pvm_master_function_evaluation(double& f,
01903     independent_variables& x,const dvector & g,int nvar);
01904   dmatrix dep_hess_routine(const dvariable& dep);
01905   void top_mcmc_routine(int,int,double,int);
01906   void mcmc_routine(int,int,double,int);
01907   void sgibbs_mcmc_routine(int,int,double,int);
01908   void hybrid_mcmc_routine(int,int,double,int);
01909   double pvm_master_get_monte_carlo_value(int nvar,
01910     const dvector& x);
01911   void pvm_slave_get_monte_carlo_value(int nvar);
01912   void mcmc_eval(void);
01913   //void hess_routine_and_constraint(int);
01914   void hess_routine_and_constraint(int iprof, const dvector& g,
01915     dvector& fg);
01916   dmatrix diag_hess_routine(void);
01917   void hess_inv(void);
01918   void depvars_routine(void);
01919   void sd_routine(void);
01920   int ef_(double * f, double * x);
01921   int constrained_minimization2(int _n,int _nh, int _ng,dvector& __x);
01922   static int constraint_exit_number;
01923   void get_bigS(int ndvar,int nvar1,int nvar,
01924     dmatrix& S,dmatrix& BS,dvector& scale);
01925 
01926 #ifdef CURVE_CORRECT
01927   void constraint_hess_routine(int ip);
01928   void get_curvature_correction_factors(int ip,
01929     dvector& g, const int underflow_flag, const dvector& eigenvals,
01930     dvector& curvcor);
01931 #endif
01932   double projected_jacobian(const dvector& g, const dvector& xscale);
01933 
01934   dvariable random_effects_maximization(const dvar_vector& v);
01935   void prof_minimize_re(int iprof, double sigma,
01936     double new_value, const double& fprof, const int underflow_flag,
01937     double global_min, const double& penalties, const double& final_weight);
01938   void prof_minimize(int iprof, double sigma,
01939     double new_value, const double& fprof, const int underflow_flag,
01940     double global_min, const double& penalties, const double& final_weight);
01941   void begin_gauss_hermite_stuff(void);
01942   void begin_funnel_stuff(void);
01943   void end_gauss_hermite_stuff(void);
01944 
01945   void prof_hess_routine(int ip,double new_value,double final_weight);
01946 
01947   void quasi_newton_minimizer1(int nvar,double _crit,
01948     double& f, const independent_variables& x,const dvector& g);
01949 
01950   double hess_determinant(int underflow_flag);
01951 
01952 #ifndef CURVE_CORRECT
01953   void normalize_posterior_distribution(double udet,
01954     const dvector& siglevel, const ofstream& ofs2, int num_pp,
01955     const dvector& all_values, const dmatrix& actual_value, double global_min,
01956     int offset, const dmatrix& lprof, const dmatrix& ldet, const dmatrix& xdist,
01957     const dmatrix& penalties);
01958   //  dmatrix& penalties, const hdmatrix& lg_jacob);
01959 #else
01960   void normalize_posterior_distribution(double udet,
01961     const dvector& siglevel, const ofstream& ofs2, int num_pp,
01962     const dvector& all_values, const dmatrix& actual_value,
01963     double global_min,
01964     int offset, const dmatrix& lprof, const dmatrix& ldet,
01965     const dmatrix& xdist,
01966     const d3_array& eigenvals,d3_array& curvcor);
01967 #endif
01968   void get_particular_grad(int iprof, int nvar, const dvector& fg,
01969     const dvector& g);
01970   double projected_hess_determinant(const dvector& g, const int underflow_flag,
01971     const dvector& xscale, const double& ln_det_proj_jac);
01972 //  double projected_hess_determinant(const dvector& fg, const dvector& g,
01973   //const int underflow_flag, const dvector& xscale,
01974   // const double& ln_det_proj_jac);
01975   double projected_hess_determinant(const dvector& g,const int underflow_flag);
01976   double projected_hess_determinant(const dmatrix& hh, const dvector& g,
01977     const int underflow_flag);
01978  //double projected_hess_determinant(const dvector& g, const int underflow_flag,
01979     //dvector& xscale, const double& ln_det_proj_jac);
01980 
01981   double projected_hess_determinant(const dvector& fg, const dvector& g,
01982     const int underflow_flag, const dvector& curvscale, const dvector& xscale,
01983     double& ln_det_proj_jac, const double& tmp, const dmatrix& hesses);
01984   double diag_projected_hess_determinant(const dvector& g,
01985     const int underflow_flag, dmatrix& dh);
01986     double unrestricted_hess_determinant(void);
01987   void monte_carlo_routine(void);
01988   double get_monte_carlo_value(int nvar, const independent_variables& x);
01989   double get_monte_carlo_value(int nvar, const independent_variables& x,
01990     dvector& g);
01991   double get_hybrid_monte_carlo_value(int nvar,const independent_variables& x,
01992     dvector& g);
01993   void mcmc_computations(void);
01994   void pvm_slave_mcmc_computations(void);
01995   void pvm_master_mcmc_computations(void);
01996   double get_monte_carlo_value(int nvar, const dvector& x);
01997   void sob_routine(int nmcmc,double dscale,int restart_flag);
01998   void sobol_importance_routine(int nmcmc,int iseed0,double dscale,
01999   int restart_flag);
02000   void pvm_master_mcmc_routine(int nmcmc,int iseed0,double dscale,
02001     int restart_flag);
02002 #if defined(USE_ADPVM)
02003   void pvm_slave_mcmc_routine(void);
02004 #else
02005   void pvm_slave_mcmc_routine(void) {}
02006 #endif
02007   void trust_region_update(int nvar,int _crit,
02008     independent_variables& x,const dvector& _g,const double& _f);
02009 
02010   void multint4(int n, const dvar_vector& a, const dvar_vector& b,
02011     const dvar_vector& h, double al, int m, const dvariable& e,
02012     const dvariable& aint1, const dvariable& aint2, dvariable& aint3,
02013     const dvariable& aint4, const int& key, PMFVIV4 f);
02014 
02015   void multint(int n, const dvar_vector& a, const dvar_vector& b,
02016     const dvar_vector& h, double al, int m, const dvariable& e,
02017     const dvariable& aint, const int& key, PMFVI f);
02018 
02019   virtual void set_runtime(void);
02020   virtual void set_runtime_maxfn(const char *);
02021   virtual void set_runtime_crit(const char *);
02022   dvariable adromb(PMF,double a,double b,int ns=9);
02023   dvariable adromb(PMF, const dvariable& a, double b, int ns = 9);
02024   dvariable adromb(PMF, const dvariable& a, const dvariable& b, int ns = 9);
02025   dvariable adromb(PMF, double a, const dvariable& b, int ns = 9);
02026 
02027   dvariable adrombo(PMF,double a,double b,int ns=9);
02028   dvariable adrombo(PMF, const dvariable& a, double b, int ns = 9);
02029   dvariable adrombo(PMF, const dvariable& a, const dvariable& b, int ns = 9);
02030   dvariable adrombo(PMF, double a, const dvariable& b,int ns = 9);
02031 
02032   dvariable trapzd(void*,double a,double b,int n);
02033   dvariable trapzd(PMF,double a,double b,int n);
02034   dvariable trapzd(PMF, double a, const dvariable& b, int n);
02035   dvariable trapzd(PMF, const dvariable& a, double b, int n);
02036   dvariable trapzd(PMF, const dvariable& a, const dvariable& b, int n);
02037 
02038   dvariable midpnt(PMF,double a,double b,int n);
02039   dvariable midpnt(PMF, double a, const dvariable& b, int n);
02040   dvariable midpnt(PMF, const dvariable& a, double b, int n);
02041   dvariable midpnt(PMF, const dvariable& a, const dvariable& b, int n);
02042 
02043   virtual void * mycast() { return (void*)this;}
02044 
02045   void neldmead(int n, dvector& _start, dvector& _xmin, double *ynewlo,
02046     double reqmin, double delta,int *icount, int *numres, int *ifault);
02047   void adamoeba(const dmatrix& p, const dvector& y, int ndim, double ftol,
02048      int maxfn);
02049   void set_initial_simplex(const dmatrix& p, const dvector& y, int nvar,
02050     const dvector& x, double delta);
02051   double amxxx(const dmatrix& p, const dvector& y, const dvector& psum,
02052     int ndim, int ihi, double fac);
02053   friend class equality_constraint_vector;
02054   friend class inequality_constraint_vector;
02055   void quasi_newton_block(int nvar,int crit,independent_variables& x,
02056     const dvector& g,const double& f);
02057   void limited_memory_quasi_newton_block(int nvar,int _crit,
02058     independent_variables& x,const dvector& _g,const double& _f,
02059     int nsteps);
02060 
02061   void function_evaluation_block_pvm_slave_random_effects(int nvar,int _crit,
02062     independent_variables& x,const dvector& g,const double& f);
02063   void quasi_newton_block_pvm_master_random_effects(int nvar,int _crit,
02064     independent_variables& x,const dvector& g,const double& f);
02065   void function_evaluation_block_pvm_slave_random_effects(void);
02066   void hess_routine_random_effects(void);
02067   void quasi_newton_block_pvm_master(int nvar,int _crit,
02068     independent_variables& x,const dvector& g,const double& f);
02069   void hess_routine_noparallel_random_effects(void);
02070 #if defined(USE_ADPVM)
02071   void function_evaluation_block_pvm_slave(void);
02072   void hess_routine_slave_random_effects(void);
02073 #endif
02074   dvariable do_gauss_hermite_integration(void);
02075   void end_df1b2_funnel_stuff(void);
02076 
02077 private:
02078   dvariable do_gauss_hermite_integration_multi(void);
02079 };
02080 
02081 cifstream& operator>>(const cifstream& s, const param_init_number& x);
02082 cifstream& operator>>(const cifstream& s, const param_init_vector& v);
02083 cifstream& operator>>(const cifstream& s, const param_init_matrix& m);
02084 ostream& operator<<(const ostream& s, const label_class& lc);
02085 
02090 class stddev_params
02091 {
02092 public:
02093   // this should be a resizeable array
02094   static stddev_params * stddevptr[150];
02095   // this should be a resizeable array
02096   static stddev_params * stddev_number_ptr[150];
02097   static void get_all_sd_values(const dvector& x, const int& ii);
02098   static int num_stddev_params;
02099   static int num_stddev_number_params;
02100   static ivector copy_all_number_offsets(void);
02101   void allocate(void){;};
02102   static int num_stddev_calc(void);
02103   static int num_stddev_number_calc(void);
02104   static void get_stddev_number_offset(void);
02105 public:
02106   stddev_params(void){;}
02107   virtual void setindex(int);
02108   virtual int getindex(void);
02109   virtual int size_count(void)=0; // get the number of active parameters
02110   virtual void set_dependent_variables(void)=0;
02111   virtual void copy_value_to_vector(const dvector&,const int&) = 0;
02112   virtual void get_sd_values(const dvector& x, const int& ii) = 0;
02113   //get the number of active parameters
02114   static void copy_all_values(const dvector& x, const int& ii);
02115   //get the number of active parameters
02116   static void copy_all_number_values(const dvector& x, const int& ii);
02117   virtual void add_to_list(void);
02118   virtual void add_to_gui_list(void);
02119   virtual const char * label()=0;
02120   friend class function_minimizer;
02121 };
02122 
02127 class likeprof_params
02128 {
02129   double stepsize;
02130   int    stepnumber;
02131 protected:
02132 public:
02133   static likeprof_params * likeprofptr[500]; // this should be a
02134                                                // resizeable array
02135   static int num_likeprof_params;
02136   void allocate(void){;};
02137   static int num_stddev_calc(void);
02138 public:
02139   likeprof_params(void);
02140   virtual void add_to_list(void);
02141   virtual const char * label()=0;
02142   virtual dvariable variable(void)=0;
02143   virtual double& get_sigma(void)=0;
02144   virtual double get_value(void)=0;
02145   double get_stepsize(void);
02146   int get_stepnumber(void);
02147   void set_stepsize(double);
02148   void set_stepnumber(int);
02149   friend class function_minimizer;
02150 };
02151 
02156 class param_stddev_vector: public named_dvar_vector , stddev_params
02157 {
02158   dvector sd;
02159     virtual int size_count(void); // get the number of active parameters
02160     virtual const char * label(void);
02161     param_stddev_vector();
02162     void allocate(int imin,int imax,const char * s="UNNAMED");
02163     virtual void set_dependent_variables(void);
02164     friend class model_parameters;
02165     virtual void copy_value_to_vector(const dvector& x, const int& ii);
02166     virtual void get_sd_values(const dvector& x, const int& ii);
02167   param_stddev_vector& operator=(const dvar_vector& m);
02168   param_stddev_vector& operator=(const dvector& m);
02169   param_stddev_vector& operator=(const double m);
02170 };
02171 
02176 class param_stddev_number: public named_dvariable , public stddev_params
02177 {
02178   double sd;
02179   int index;
02180   void allocate(const char *s="UNNAMED");
02181   virtual void setindex(int);
02182   virtual int getindex(void);
02183   virtual int size_count(void); // get the number of active parameters
02184   virtual const char * label(void);
02185   virtual void copy_value_to_vector(const dvector& x, const int& ii);
02186   virtual void get_sd_values(const dvector& x, const int& ii);
02187 protected:
02188   param_stddev_number();
02189   friend class model_parameters;
02190   virtual void set_dependent_variables(void);
02191   param_stddev_number& operator=(const prevariable&);
02192   param_stddev_number& operator=(const double);
02193 };
02194 
02199 class param_likeprof_number: public param_stddev_number ,
02200     public likeprof_params
02201 {
02202     double sigma;
02203     void allocate(const char *s="UNNAMED");
02204     virtual int size_count(void); // get the number of active parameters
02205     virtual const char * label(void);
02206     virtual double& get_sigma(void){return sigma;}
02207     virtual double get_value(void){return value(*this);}
02208     //void copy_value_to_vector(const dvector& x, const int& ii);
02209     virtual dvariable variable(void){ return dvariable(*this);}
02210     param_likeprof_number();
02211     friend class model_parameters;
02212 public:
02213     param_likeprof_number& operator=(const prevariable&);
02214     param_likeprof_number& operator=(const double);
02215 };
02216 
02221 class param_stddev_matrix: public named_dvar_matrix , stddev_params
02222 {
02223   dmatrix sd;
02224     virtual int size_count(void);
02225     //virtual void read_value(void);
02226     virtual const char * label(void);
02227     void allocate(int rmin,int rmax,int cmin,int cmax,
02228     const char * s="UNNAMED");
02229   param_stddev_matrix(void);
02230   friend class model_parameters;
02231   virtual void set_dependent_variables(void);
02232   virtual void get_sd_values(const dvector& x, const int& ii);
02233   virtual void copy_value_to_vector(const dvector& x, const int& ii);
02234   param_stddev_matrix& operator=(const double m);
02235   param_stddev_matrix& operator=(const dmatrix& m);
02236   param_stddev_matrix& operator=(const dvar_matrix& m);
02237 };
02238 
02243   class objective_function_value : public named_dvariable
02244   {
02245   public:
02246     static objective_function_value * pobjfun;
02247     static double fun_without_pen;
02248     static double gmax;
02249     objective_function_value();
02250     objective_function_value& operator=(const prevariable& v);
02251     objective_function_value& operator=(const double v);
02252   };
02253 
02254   int withinbound(int lb,int n,int ub);
02255 
02256 double cumd_cauchy(const double& x);
02257 double density_cauchy(const double& x);
02258 double log_density_cauchy(const double& x);
02259 double inv_cumd_cauchy(const double& x);
02260 
02261 double cumd_mixture(const double& x);
02262 double inv_cumd_mixture(const double& y);
02263 double cumd_mixture_02(const double& x);
02264 double inv_cumd_mixture_02(const double& y);
02265 
02266 #if defined _ADM_HIGHER_ARRAYS__
02267 
02272 class param_init_matrix: public named_dvar_matrix,public initial_params
02273 {
02274   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
02275 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
02276 //  const int& ii);
02277   virtual void add_value(const dvector&, const dvector&, const int&,
02278     const double&, const dvector&);
02279   virtual void add_value(const dvector&, const int&);
02280   virtual void get_jacobian(const dvector&, const dvector&, const int&);
02281 public:
02282   virtual void set_value(const dvar_vector& x, const int& ii,
02283     const dvariable& pen);
02284   virtual void copy_value_to_vector(const dvector& x, const int& ii);
02285   virtual void restore_value_from_vector(const dvector&, const int&);
02286   virtual void set_value_inv(const dvector& x, const int& ii);
02287   virtual int size_count(void);
02288   virtual void save_value(void);
02289   virtual void bsave_value(void);
02290   virtual void save_value(const ofstream& ofs, int prec);
02291   virtual void restore_value(const ifstream& ifs);
02292   void report_value(void);
02293   //virtual void read_value(void);
02294   virtual const char * label(void);
02295   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
02296   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
02297   virtual void curv_scale(const dvector& d, const dvector& x, const int&);
02298   virtual void hess_scale(const dvector& d, const dvector& x, const int& ii){;};
02299 
02300 public:
02301   void allocate(int rmin,int rmax,int cmin,int cmax,
02302     int phase_start=1,const char * = "UNNAMED");
02303   void allocate(int rmin,int rmax,int cmin,int cmax,
02304     const char * = "UNNAMED");
02305   param_init_matrix(void);
02306   param_init_matrix& operator = (const dmatrix& m);
02307   param_init_matrix& operator = (const dvar_matrix& m);
02308   param_init_matrix& operator = (const dvariable& m);
02309   param_init_matrix& operator = (const double& m);
02310 };
02311 #endif // #if defined _ADM_HIGER_ARRAYS__
02312 
02317 class param_init_d3array: public named_dvar3_array,public initial_params
02318 {
02319 public:
02320 #if defined(USE_SHARE_FLAGS)
02321   virtual void setshare(const index_type& sf,const index_type& af);
02322   virtual void shared_set_value_inv(const dvector&,const int&);
02323   virtual void shared_set_value(const dvar_vector&,const int&,
02324     const dvariable& pen);
02325   virtual int shared_size_count(void); // get the number of active parameters
02326 #endif
02327 
02328   virtual void sd_vscale(const dvar_vector& d,const dvar_vector& x,
02329     const int& ii);
02330   virtual void dev_correction(const dmatrix&, const int&);
02331   virtual void curv_scale(const dvector& d, const dvector& x, const int& ii);
02332   virtual void set_simulation_bounds(const dmatrix& symbds, const int& ii);
02333 //virtual void set_simulation_bounds(const dmatrix&, const dvector& symbds,
02334 //  const int& ii);
02335   virtual void add_value(const dvector&, const dvector&, const int&,
02336     const double&, const dvector&);
02337   virtual void add_value(const dvector&, const int&);
02338   virtual void get_jacobian(const dvector&, const dvector&, const int&);
02339   virtual void set_value(const dvar_vector& x, const int& ii,
02340     const dvariable& pen);
02341   virtual void copy_value_to_vector(const dvector& x, const int& ii);
02342   virtual void restore_value_from_vector(const dvector&,const int&);
02343   virtual void set_value_inv(const dvector& x, const int& ii);
02344   virtual int size_count(void);
02345   virtual void save_value(void);
02346   virtual void bsave_value(void);
02347   virtual void save_value(const ofstream& ofs, int prec);
02348   virtual void restore_value(const ifstream& ifs);
02349   void report_value(void);
02350   //virtual void read_value(void);
02351   virtual const char * label(void);
02352   virtual void sd_scale(const dvector& d, const dvector& x, const int& ii);
02353   virtual void mc_scale(const dvector& d, const dvector& x, const int&);
02354   virtual void hess_scale(const dvector& d, const dvector& x, const  int& ii);
02355 
02356 public:
02357 #if defined(USE_ADPVM)
02358   virtual void pvm_pack(void) { cerr << "Error" << endl; ad_exit(1);}
02359   virtual void pvm_unpack(void) { cerr << "Error" << endl; ad_exit(1);}
02360 #endif
02361 
02362   void allocate(const ad_integer& sl,const ad_integer& sh,
02363     const index_type& nrl,const index_type& nrh,
02364     const index_type& ncl,const index_type& nch,const char * s="UNNAMED");
02365 
02366   void allocate(const ad_integer& sl,const ad_integer& sh,
02367     const index_type& nrl,const index_type& nrh,
02368     const index_type& ncl,const index_type& nch,int phase_start=1,
02369     const char * s="UNNAMED");
02370 
02371   void allocate(int smin,int smax,int rmin,int rmax,int cmin,int cmax,
02372     int phase_start=1,const char * = "UNNAMED");
02373   void allocate(int smin,int smax,int rmin,int rmax,int cmin,int cmax,
02374     const char * = "UNNAMED");
02375   param_init_d3array(void);
02376 };
02377 
02382 class dll_param_init_d3array: public param_init_d3array
02383 {
02384   double * d;
02385 public:
02386   void allocate(double* _d,int hmin,int hmax,
02387     int rmin,int rmax,int cmin,int cmax,
02388     int phase_start=1,const char * = "UNNAMED");
02389   void allocate(double * _d,int hmin,int hmax,int rmin,int rmax,
02390     int cmin,int cmax,const char * = "UNNAMED");
02391   virtual ~dll_param_init_d3array();
02392   dll_param_init_d3array(){d=NULL;}
02393   dll_param_init_d3array& operator=(const d3_array&);
02394   dll_param_init_d3array& operator=(const dvar3_array&);
02395 };
02396 
02401 class dll_param_d3array: public named_dvar3_array
02402 {
02403   double * d;
02404 public:
02405   void allocate(double* _d,int hmin,int hmax,
02406     int rmin,int rmax,int cmin,int cmax,
02407     int phase_start=1,const char * = "UNNAMED");
02408   void allocate(double * _d,int hmin,int hmax,int rmin,int rmax,
02409     int cmin,int cmax,const char * = "UNNAMED");
02410   virtual ~dll_param_d3array();
02411   dll_param_d3array(){d=NULL;}
02412   dll_param_d3array& operator=(const d3_array&);
02413   dll_param_d3array& operator=(const dvar3_array&);
02414 };
02415 
02416 
02417 //double set_value_mc(const double& x, const double fmin, const double fmax);
02418 
02419 void set_value_mc(const dvar_vector& x, const dvar_vector& v, const int& ii,
02420   const double fmin, const double fmax);
02421 
02422 double set_value_inv_mc(double v,double fmin,double fmax);
02423 
02424 double set_value_inv_mc(const prevariable& v, double fmin, double fmax);
02425 
02426 void set_value_inv_mc(const dvar_vector& x, const dvector& v, const int& ii,
02427   const double fmin, const double fmax);
02428 
02429 //double set_value_inv_mc(const dvector&, const dvector& x, int ii, double minb,
02430 //  double maxb);
02431 
02432 double set_value_mc(double z,double min,double max);
02433 double ndfboundp( double x, double fmin, double fmax,const double& fpen);
02434 double ndfboundp_mc( double x, double fmin, double fmax,const double& fpen);
02435 
02436 void copy_value_from_vector(const double& _sd,const dvector& x,const int & _ii);
02437 void copy_value_from_vector(const dvector& _sd,const dvector& x,
02438   const int & _ii);
02439 void copy_value_from_vector(const dmatrix& _sd,const dvector& x,
02440   const int & _ii);
02441 
02442 dvector bounded_multivariate_normal(int nvar, const dvector& a1,
02443   const dvector& b1, dmatrix& ch, const double& _wght,
02444   const random_number_generator & rng);
02445 
02446 void bounded_multivariate_normal_mcmc(int nvar, const dvector& a1,
02447   const dvector& b1, dmatrix& ch, const double& wght, const dvector& y,
02448   const random_number_generator & rng);
02449 
02450 void bounded_multivariate_uniform_mcmc(int nvar, const dvector& a1,
02451   const dvector& b1, dmatrix& ch, const double& wght, const dvector& y,
02452   const random_number_generator & rng);
02453 
02454 dvector bounded_multivariate_normal(int nvar, const dvector& a1,
02455   const dvector& b1, dmatrix& ch, const double& lprob,
02456   random_number_generator &rng);
02457 
02458 dvector bounded_multivariate_normal_sobol(int nvar, const dvector& a1,
02459   const dvector& b1, dmatrix& ch, const double& lprob,
02460   const random_number_generator &rng);
02461 
02462 dvector probing_bounded_multivariate_normal(int nvar, const dvector& a1,
02463   const dvector& b1, dmatrix& ch, const double& lprob, double pprobe,
02464   const random_number_generator &rng);
02465 
02466 dvector bounded_multivariate_uniform(int nvar, const dvector& a1,
02467   const dvector& b1, dmatrix& ch, const double& lprob,
02468   random_number_generator &rng);
02469 
02470 dvector bounded_robust_multivariate_normal(int nvar, const dvector& a1,
02471   const dvector& b1, dmatrix& ch, const dmatrix& ch3, double contaminant,
02472   const double& _wght, random_number_generator &rng);
02473 
02474 void probing_bounded_multivariate_normal_mcmc(int nvar, const dvector& a1,
02475   const dvector& b1, dmatrix& ch, const double& wght, const dvector& v,
02476   double pprobe, const random_number_generator &rng);
02477 
02478 /*
02479 int option_match(int argc,char * argv[], const char * string);
02480 int option_match(int argc,char * argv[], const char * string, const int& nopt);
02481 int option_match(char * s, const char * string, const int& _nopt);
02482 int option_match(char * s, const char * string);
02483 */
02484 
02485 double inv_cumd_exp(double x);
02486 double cumd_exp(double x);
02487 
02488 double ffmax(double a,double b);
02489 double ffmin(double a,double b);
02490 
02491 void check_datafile_pointer(void * p);
02492 
02493 adstring get_reportfile_name(void);
02494 
02495 void ad_make_code_reentrant(void);
02496 
02497 char** parse_dll_options(char *pname, const int& argc,
02498   char * dll_options);
02499 
02500 void parse_dll_options(char *pname, const int& argc, char *dll_options,
02501   char *** pargv);
02502 
02503 char** no_dll_options(char *pname, const int& argc);
02504 
02505 void cleanup_argv(int nopt,char *** pa);
02506 
02507 void get_sp_printf(void);
02508 
02509 void do_dll_housekeeping(int argc,char ** argv);
02510 
02511 void adwait(double);
02512 
02513 int ad_get_commandline_option(const char *option_label, const int &option_value,
02514   const char * error_message);
02515 
02520 class param_init_vector_vector
02521 {
02522   param_init_vector * v;
02523   int index_min;
02524   int index_max;
02525   double_index_type * it;
02526 
02527 public:
02528   param_init_vector_vector();
02529   ~param_init_vector_vector();
02530 
02531   void set_scalefactor(double s);
02532   void set_scalefactor(const dvector& s);
02533   dvector get_scalefactor(void);
02534 
02535 #if defined(OPT_LIB)
02536   param_init_vector& operator [] (int i) { return v[i];}
02537   param_init_vector& operator () (int i) { return v[i];}
02538   prevariable operator () (int i,int j) { return v[i][j];}
02539 #else
02540   param_init_vector& operator [] (int i);
02541   param_init_vector& operator () (int i);
02542   prevariable operator () (int i,int j);
02543 #endif
02544 
02545   void allocate(int min1,int max1,const index_type& min,
02546      const index_type& max,const index_type& phase_start,
02547      const char * s);
02548 
02549   void allocate(int min1,int max1,const index_type& min,
02550      const index_type& max,const char * s);
02551 
02552   bool allocated() const { return v != NULL; }
02553   int indexmin(void) {return (index_min);}
02554   int indexmax(void) {return (index_max);}
02555   void set_initial_value(const double_index_type& it);
02556   void deallocate(void);
02557 };
02558 
02563 class param_init_bounded_vector_vector
02564 {
02565   param_init_bounded_vector* v;
02566   int index_min;
02567   int index_max;
02568   double_index_type* it;
02569 public:
02570   param_init_bounded_vector_vector();
02571   ~param_init_bounded_vector_vector();
02572 
02573   void set_scalefactor(double s);
02574   void set_scalefactor(const dvector& s);
02575   dvector get_scalefactor(void);
02576 #if defined(OPT_LIB)
02577   param_init_bounded_vector& operator [] (int i) { return v[i];}
02578   param_init_bounded_vector& operator () (int i) { return v[i];}
02579   prevariable operator () (int i,int j) { return v[i][j];}
02580 #else
02581   param_init_bounded_vector& operator [] (int i);
02582   param_init_bounded_vector& operator () (int i);
02583   prevariable operator () (int i,int j);
02584 #endif
02585 
02586   void allocate(int min1,int max1,
02587     const index_type& min,
02588     const index_type& max,
02589     const double_index_type& dmin,
02590     const double_index_type& dmax,
02591     const index_type& phase_start,
02592     const char* s);
02593 
02594   void allocate(int min1,int max1,
02595     const index_type& min,
02596     const index_type& max,
02597     const double_index_type& dmin,
02598     const double_index_type& dmax,
02599     const char* s);
02600 
02601   bool allocated() const { return v != NULL; }
02602   int indexmin() const {return index_min;}
02603   int indexmax() const {return index_max;}
02604   void deallocate(void);
02605   void set_initial_value(const double_index_type& it);
02606 };
02607 
02612 class param_init_matrix_vector
02613 {
02614   param_init_matrix* v;
02615   int index_min;
02616   int index_max;
02617   double_index_type* it;
02618 
02619 public:
02620   param_init_matrix_vector();
02621    ~param_init_matrix_vector();
02622 
02623   void set_scalefactor(double s);
02624   void set_scalefactor(const dvector& s);
02625   dvector get_scalefactor(void);
02626 
02627 #if defined(OPT_LIB)
02628   param_init_matrix& operator [] (int i) { return v[i];}
02629   param_init_matrix& operator () (int i) { return v[i];}
02630   dvar_vector& operator () (int i,int j) { return v[i][j];}
02631   prevariable operator () (int i,int j,int k) { return v[i](j,k);}
02632 #else
02633   param_init_matrix& operator [] (int i);
02634   param_init_matrix& operator () (int i);
02635   dvar_vector& operator () (int i,int j);
02636   prevariable operator () (int i,int j,int k);
02637 #endif
02638   void allocate(int min0,int max0,const index_type& min,
02639      const index_type& max,const index_type& min1,
02640      const index_type& max1,const index_type& phase_start,
02641      const char * s);
02642 
02643   void allocate(int min0,int max0,const index_type& min,
02644      const index_type& max,const index_type& min1,
02645      const index_type& max1,const char * s);
02646 
02647   bool allocated() const { return v != NULL; }
02648   int indexmin() const { return index_min;}
02649   int indexmax() const { return index_max;}
02650   void set_initial_value(const double_index_type& it);
02651   void deallocate(void);
02652 };
02653 
02658 class param_init_bounded_matrix_vector
02659 {
02660   param_init_bounded_matrix * v;
02661   int index_min;
02662   int index_max;
02663   double_index_type * it;
02664 public:
02665   param_init_bounded_matrix_vector();
02666   ~param_init_bounded_matrix_vector();
02667 
02668   void set_scalefactor(double s);
02669   void set_scalefactor(const dvector& s);
02670   dvector get_scalefactor(void);
02671 
02672   void allocate(int min1,int max1,
02673     const index_type& min, const index_type& max, const index_type& min2,
02674     const index_type& max2, const double_index_type& dmin2,
02675     const double_index_type& dmax2, const index_type& phase_start,
02676     const char * s);
02677 
02678   void allocate(int min1,int max1,
02679     const index_type& min, const index_type& max, const index_type& min2,
02680     const index_type& max2, const double_index_type& dmin2,
02681     const double_index_type& dmax2,const char * s);
02682 
02683 #if defined(OPT_LIB)
02684   param_init_bounded_matrix& operator [] (int i) { return v[i];}
02685   param_init_bounded_matrix& operator () (int i) { return v[i];}
02686   dvar_vector& operator () (int i,int j) { return v[i][j];}
02687   prevariable operator () (int i,int j,int k) { return v[i](j,k);}
02688 #else
02689   param_init_bounded_matrix& operator [] (int i);
02690   param_init_bounded_matrix& operator () (int i);
02691   dvar_vector& operator () (int i,int j);
02692   prevariable operator () (int i,int j,int k);
02693 #endif
02694 
02695   bool allocated() const { return v != NULL; }
02696   int indexmin() const { return index_min; }
02697   int indexmax() const { return index_max; }
02698   void set_initial_value(const double_index_type& it);
02699   void deallocate(void);
02700 };
02701 
02706 class param_init_number_vector
02707 {
02708   param_init_number* v;
02709   int index_min;
02710   int index_max;
02711   double_index_type* it;
02712 public:
02713   param_init_number_vector();
02714   ~param_init_number_vector();
02715 
02716   void set_scalefactor(double s);
02717   void set_scalefactor(const dvector& s);
02718   dvector get_scalefactor(void);
02719 
02720 #if defined(OPT_LIB)
02721   param_init_number& operator[](int i) { return v[i];}
02722   param_init_number& operator()(int i) { return v[i];}
02723 #else
02724   param_init_number& operator[](int i);
02725   param_init_number& operator()(int i);
02726 #endif
02727 
02728   void allocate(int min1,int max1,const index_type& phase_start,
02729      const char * s);
02730 
02731   void allocate(int min1,int max1,const char * s);
02732 
02733   bool allocated() const { return v != NULL; }
02734   int indexmin() const { return index_min; }
02735   int indexmax() const { return index_max; }
02736   void set_initial_value(const double_index_type& it);
02737   void deallocate(void);
02738 };
02739 
02747  class data_matrix;
02748 class param_init_bounded_number_vector
02749 {
02750   param_init_bounded_number* v;
02751   int index_min;
02752   int index_max;
02753   double_index_type* it;
02754 public:
02755   param_init_bounded_number_vector();
02756   ~param_init_bounded_number_vector();
02757 
02758   void set_scalefactor(double s);
02759   void set_scalefactor(const dvector& s);
02760   dvector get_scalefactor(void);
02761 
02762 #if defined(OPT_LIB)
02763   param_init_bounded_number& operator [] (int i) { return v[i];}
02764   param_init_bounded_number& operator () (int i) { return v[i];}
02765 #else
02766   param_init_bounded_number& operator [] (int i);
02767   param_init_bounded_number& operator () (int i);
02768 #endif
02769 
02770   void allocate(int min1,int max1,const double_index_type & bmin,
02771     const double_index_type & bmax,const index_type& phase_start,
02772     const char * s);
02773 
02774   void allocate(int min1,int max1,const double_index_type & bmin,
02775      const double_index_type & bmax,const char * s);
02776 
02777   // Added by Steve Martell, Jan 18, 2014.
02778   void allocate(const data_matrix &m, const char *s);
02779 
02780   bool allocated() const { return v != NULL; }
02781   int indexmin() const { return (index_min); }
02782   int indexmax() const { return (index_max); }
02783   void set_initial_value(const double_index_type& it);
02784   void deallocate(void);
02785 };
02786   extern int traceflag;
02787   void tracing_message(int traceflag,const char *s,int *pn);
02788   void tracing_message(int traceflag,const char *s,double *pn);
02789   void set_gauss_covariance_matrix(const dll_data_matrix& m);
02790   void set_gauss_covariance_matrix(const dmatrix& m);
02791   void set_covariance_matrix(const dll_data_matrix& m);
02792   void set_covariance_matrix(const dmatrix& m);
02793 
02794   //ostream& operator<<(const ostream&, const param_init_number_vector);
02795   //ostream& operator<<(const ostream&, const param_init_bounded_number_vector);
02796   //ostream& operator<<(const ostream&, const param_init_vector_vector);
02797   //ostream& operator<<(const ostream&, const param_init_bounded_vector_vector);
02798   //ostream& operator<<(const ostream&, const param_init_matrix_vector);
02799   //ostream& operator<<(const ostream&, const param_init_bounded_matrix_vector);
02800 
02805   class vector_kludge : public dvar_vector
02806   {
02807     public:
02808      vector_kludge(const param_init_number_vector &);
02809      vector_kludge(const param_init_bounded_number_vector &);
02810   };
02815   class matrix_kludge : public dvar_matrix
02816   {
02817     public:
02818      matrix_kludge(const param_init_vector_vector &);
02819      matrix_kludge(const param_init_bounded_vector_vector &);
02820   };
02821 
02822 void read_covariance_matrix(const dmatrix& S, int nvar, int& hbf,
02823   dvector& sscale);
02824 
02825 dvector value(const param_init_number_vector& t);
02826 dvector value(const param_init_bounded_number_vector& t);
02827 //dvector value(const param_init_bounded_number_matrix& t);
02828 //dvector value(const param_init_vector_vector& t);
02829 //dvector value(const param_init_bounded_vector_vector& t);
02830 
02831 dvector read_old_scale(int & old_nvar);
02832 
02833 int withinbound(int lb,int n,int ub);
02834 
02835 #if defined(_MSC_VER)
02836 #  if defined(min)
02837 #    undef min
02838 #  endif
02839 #  if defined(max)
02840 #    undef max
02841 #  endif
02842 #endif
02843 
02844 #include "param_init_bounded_number_matrix.h"
02845 
02846 #endif