ADMB Documentation  11.1.2495
 All Classes Files Functions Variables Typedefs Friends Defines
adrndeff.h
Go to the documentation of this file.
00001 /*
00002  * $Id: adrndeff.h 1919 2014-04-22 22:02:01Z 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  */
00046 #if !defined(__RANDEFFECTS__)
00047 #  define __RANDEFFECTS__
00048 
00049 #  include <admodel.h>
00050 
00051 class dcompressed_triplet;
00052 class hs_symbolic;
00053 class dvar_compressed_triplet;
00054 
00058 class random_effects_number : public param_init_number
00059 {
00060   virtual void set_random_effects_active();
00061   virtual void set_random_effects_inactive();
00062   virtual void set_only_random_effects_active();
00063   virtual void set_only_random_effects_inactive();
00064 };
00065 
00069 class random_effects_bounded_number : public param_init_bounded_number
00070 {
00071   virtual void set_random_effects_active();
00072   virtual void set_random_effects_inactive();
00073   virtual void set_only_random_effects_active();
00074   virtual void set_only_random_effects_inactive();
00075 };
00076 
00080 class random_effects_vector : public param_init_vector
00081 {
00082   virtual void set_random_effects_active();
00083   virtual void set_random_effects_inactive();
00084   virtual void set_only_random_effects_active();
00085   virtual void set_only_random_effects_inactive();
00086 };
00087 
00088 
00089 class random_effects_bounded_vector;
00090 
00094 class random_effects_bounded_vector : public param_init_bounded_vector
00095 {
00096   virtual void set_random_effects_active();
00097   virtual void set_random_effects_inactive();
00098   virtual void set_only_random_effects_active();
00099   virtual void set_only_random_effects_inactive();
00100 };
00101 
00105 class random_effects_matrix : public param_init_matrix
00106 {
00107   virtual void set_random_effects_active();
00108   virtual void set_random_effects_inactive();
00109   virtual void set_only_random_effects_active();
00110   virtual void set_only_random_effects_inactive();
00111 };
00112 
00116 class random_effects_bounded_matrix : public param_init_bounded_matrix
00117 {
00118   virtual void set_random_effects_active();
00119   virtual void set_random_effects_inactive();
00120   virtual void set_only_random_effects_active();
00121   virtual void set_only_random_effects_inactive();
00122 };
00123 
00124 class gauss_hermite_stuff;
00125 
00126 class nested_calls_shape;
00127 //class sparse_symbolic;
00128 
00132 class nested_calls_indices
00133 {
00134   imatrix * ptr1;
00135   i3_array * ptr2;
00136   i4_array * ptr3;
00137   i5_array * ptr4;
00138 public:
00139   nested_calls_indices(void) : ptr1(0),ptr2(0),ptr3(0),ptr4(0){}
00140   ivector & level1(int i) {return (*ptr1)(i);}
00141   ivector & level2(int i,int j) {return (*ptr2)(i,j);}
00142   ivector & level3(int i,int j,int k) {return (*ptr3)(i,j,k);}
00143   ivector & level4(int i,int j,int k,int l) {return (*ptr4)(i,j,k,l);}
00144   void allocate(const nested_calls_shape& nsc);
00145 };
00146 
00150 class nested_calls_shape
00151 {
00152   ivector * ptr1;
00153   imatrix * ptr2;
00154   i3_array * ptr3;
00155   i4_array * ptr4;
00156 public:
00157   nested_calls_shape(void) : ptr1(0),ptr2(0),ptr3(0),ptr4(0){}
00158   ivector * get_ptr1(void);
00159   imatrix * get_ptr2(void);
00160   i3_array * get_ptr3(void);
00161   i4_array * get_ptr4(void);
00162   void trim(void);
00163   void allocate(int);
00164   void allocate(int,int);
00165   void allocate(int,int,int);
00166   void allocate(int,int,int,int);
00167   void initialize(void);
00168   ~nested_calls_shape();
00169   int & level1(int i) {return (*ptr1)(i);}
00170   int & operator () (int i) {return (*ptr1)(i);}
00171   int & level2(int i,int j) {return (*ptr2)(i,j);}
00172   int & operator () (int i,int j) {return (*ptr2)(i,j);}
00173   int & level3(int i,int j,int k) {return (*ptr3)(i,j,k);}
00174   int & operator () (int i,int j,int k) {return (*ptr3)(i,j,k);}
00175   int & level4(int i,int j,int k,int l) {return (*ptr4)(i,j,k,l);}
00176   int & operator () (int i,int j,int k,int l) {return (*ptr4)(i,j,k,l);}
00177 };
00178 
00182 class laplace_approximation_calculator
00183 {
00184 public:
00185   int init_switch;
00186   nested_calls_indices nested_indices;
00187   nested_calls_shape nested_shape;
00188   int separable_call_level;
00189   int dd_nr_flag;
00190   dmatrix * antiepsilon;
00191   i3_array * triplet_information;
00192   imatrix * compressed_triplet_information;
00193   imatrix * calling_set;
00194   dvector *  importance_sampling_values;
00195   dvector *  importance_sampling_weights;
00196   int is_diagnostics_flag;
00197   static int saddlepointflag;
00198   static int sparse_hessian_flag;
00199   static int antiflag;
00200   int rseed;
00201   static int print_importance_sampling_weights_flag;
00202   static dvar_vector * variance_components_vector;
00203   static int where_are_we_flag; // 1 in inner opt 2 in newton-raphson 3 in
00204                                 // laplace approximation
00205   int num_separable_calls; // 1 in inner opt 2 in newton-raphson 3 in
00206   int separable_calls_counter; // 1 in inner opt 2 in newton-raphson 3 in
00207   ivector nested_separable_calls_counter;
00208   ivector nested_tree_position;
00209   ivector * num_local_re_array; // 1 in inner opt 2 in newton-raphson 3 in
00210   ivector * num_local_fixed_array; // 1 in inner opt 2 in newton-raphson 3 in
00211   int use_outliers;
00212   int use_gauss_hermite;
00213   int in_gauss_hermite_phase;
00214   int multi_random_effects;
00215   int no_function_component_flag;
00216   dvector * separable_function_difference;
00217   gauss_hermite_stuff * gh;
00218   dvar_matrix * importance_sampling_components;
00219   int importance_sampling_counter;
00220   d3_array * block_diagonal_hessian;
00221   d3_array * block_diagonal_ch;
00222   dvar3_array * block_diagonal_vch;
00223   d3_array * block_diagonal_Dux;
00224   imatrix * block_diagonal_re_list;
00225   imatrix * block_diagonal_fe_list;
00226   d3_array * block_diagonal_vhessianadjoint;
00227   dvar3_array * block_diagonal_vhessian;
00228   dvector local_dtemp;
00229   double inner_crit;
00230   double max_separable_g;
00231   double nr_crit;
00232   ivector used_flags;
00233   int have_users_hesstype;
00234   int inner_maxfn;
00235   int nr_debug;
00236   int inner_lmnflag;
00237   int inner_lmnsteps;
00238   int inner_iprint;
00239   int inner_noprintx;
00240   function_minimizer * pmin;
00241   int block_diagonal_flag;
00242   int bw;
00243   int xsize;
00244   int usize;
00245   int nvariables;
00246   int nvar;
00247   ivector minder;
00248   ivector maxder;
00249   int num_der_blocks;
00250   int num_nr_iters;
00251   int num_importance_samples;
00252   dmatrix epsilon;
00253   int hesstype;
00254   int var_flag;
00255   int have_bounded_random_effects;
00256   int isfunnel_flag;
00257   int nfunnelblocks;
00258   int bHess_pd_flag;
00259   banded_symmetric_dmatrix * bHess;
00260   banded_symmetric_dmatrix * bHessadjoint;
00261   dcompressed_triplet * sparse_triplet;
00262   ivector * sparse_iterator;
00263   int sparse_count;
00264   int sparse_count_adjoint;
00265   dcompressed_triplet * sparse_triplet2;
00266   dvar_compressed_triplet * vsparse_triplet;
00267   dcompressed_triplet * vsparse_triplet_adjoint;
00268   hs_symbolic * sparse_symbolic;
00269   hs_symbolic * sparse_symbolic2;
00270 
00271   void make_sparse_triplet(void);
00272   void check_for_need_to_reallocate(int ip);
00273   void get_newton_raphson_info(function_minimizer * pmin);
00274   void get_newton_raphson_info_slave(function_minimizer * pmin);
00275   void get_newton_raphson_info_master(function_minimizer * pmin);
00276   dvector get_newton_raphson_info_block_diagonal(function_minimizer * pmin);
00277   dvector get_newton_raphson_info_banded(function_minimizer * pmin);
00278   double do_one_feval(const dvector& x,function_minimizer * pfmin);
00279   void test_trust_region_method(function_minimizer * pmin);
00280   void get_complete_hessian(dmatrix& H,function_minimizer * pfmin);
00281   void get_complete_hessian(dmatrix& H,dvector& g,function_minimizer * pfmin);
00282   dvector lincg(dvector& x,dvector& c, dmatrix& H,double tol,double Delta,
00283     function_minimizer * pfmin,double& truef,double& e,double& f,
00284     double& fbest,int& iflag,int& iter,int maxfn);
00285   dvector test_trust_region_method(const dvector& _x,const double& _f,
00286     function_minimizer * pfmin);
00287   dmatrix get_gradient_for_hessian_calcs(const dmatrix& local_Hess,
00288     double & f);
00289   fmm fmc1;
00290   //fmmt1 fmc1;
00291   fmm fmc;
00292   dvector scale;
00293   dvector curv;
00294   dvector xadjoint;
00295   dvector check_local_uadjoint;
00296   dvector check_local_uadjoint2;
00297   dvector check_local_xadjoint;
00298   dvector check_local_xadjoint2;
00299   dvector uadjoint;
00300   dvector uhat;
00301   dvector ubest;
00302   dvector grad;
00303   dvector * grad_x_u;
00304   dvector * grad_x;
00305   dvector step;
00306   dmatrix Hess;
00307   d3_array * Hess_components;
00308   dmatrix * pHess_non_quadprior_part;
00309   imatrix * derindex;
00310   dmatrix Hessadjoint;
00311   dmatrix Dux;
00312   init_df1b2vector y;
00313   dvector get_uhat_quasi_newton_block_diagonal(const dvector& x,
00314     function_minimizer * pfmin);
00315   dvector get_uhat_quasi_newton(const dvector& x,function_minimizer * pfmin);
00316   dvector get_uhat_quasi_newton_qd(const dvector& x,function_minimizer * pfmin);
00317   void set_u_dot(int i);
00318 
00319   void do_separable_stuff_laplace_approximation_banded_adjoint
00320     (const df1b2variable& ff);
00321 
00322   dvector get_uhat_lm_newton2(const dvector& x,function_minimizer * pfmin);
00323 
00324   dvector get_uhat_lm_newton(const dvector& x,function_minimizer * pfmin);
00325   laplace_approximation_calculator(int _xsize,int _usize,int _minder,
00326     int _maxder,function_minimizer * pfmin);
00327 
00328   laplace_approximation_calculator(int _xsize,int _usize,ivector _minder,
00329     ivector _maxder, function_minimizer * pfmin);
00330 
00331   dvector operator () (const dvector& _x,const double& _f,
00332     function_minimizer * pfmin);
00333 
00334   dvector get_uhat(const dvector& x,function_minimizer * pfmin);
00335 
00336   ~laplace_approximation_calculator();
00337 
00338   void  do_separable_stuff(void);
00339   void  do_separable_stuff_newton_raphson_block_diagonal(df1b2variable&);
00340   void  do_separable_stuff_newton_raphson_banded(df1b2variable&);
00341   void  do_separable_stuff_laplace_approximation_block_diagonal(df1b2variable&);
00342   void  do_separable_stuff_laplace_approximation_banded(df1b2variable&);
00343   dvector default_calculations_check_derivatives(const dvector& _x,
00344     function_minimizer * pfmin,const double& f);
00345   dvector default_calculations(const dvector& _x,const double& _f,
00346     function_minimizer * pfmin);
00347   dvector banded_calculations(const dvector& _x,const double& _f,
00348     function_minimizer * pfmin);
00349   dvector block_diagonal_calculations(const dvector& _x,const double& _f,
00350     function_minimizer * pfmin);
00351   dvector default_calculations_parallel_master(const dvector& _x,
00352     const double& _f,function_minimizer * pfmin);
00353   void default_calculations_parallel_slave(const dvector& _x,
00354     const double& _f,function_minimizer * pfmin);
00355   void pvm_slave_function_evaluation_random_effects(void);
00356   dvector banded_calculations_trust_region_approach(const dvector& _uhat,
00357     function_minimizer * pmin);
00358   void do_newton_raphson_banded(function_minimizer * pmin,double,int&);
00359   double inner_optimization_banded(/*dvector& uhat,*/ dvector& x,
00360     function_minimizer * pfmin,int& no_converge_flag);
00361   void set_default_hessian_type(void);
00362   void check_hessian_type(const dvector& _x,function_minimizer * );
00363   void check_hessian_type(function_minimizer * pfmin);
00364   void check_hessian_type2(function_minimizer * pfmin);
00365   void do_separable_stuff_hessian_type_information(void);
00366   void get_block_diagonal_hessian(df1b2variable&);
00367   void allocate_block_diagonal_stuff(void);
00368 
00369   void do_separable_stuff_laplace_approximation_importance_sampling_adjoint
00370     (df1b2variable&);
00371   void get_hessian_components_banded_lme(function_minimizer * pfmin);
00372   dvar_matrix get_hessian_from_components_lme(function_minimizer * pfmin);
00373   dvector banded_calculations_lme(const dvector& _x,const double& _f,
00374     function_minimizer * pfmin);
00375   dvector get_gradient_lme(const dvector& x,function_minimizer * pfmin);
00376   dvector get_gradient_lme(function_minimizer * pfmin);
00377   dvector get_gradient_lme_hp(const double& x,function_minimizer * pfmin);
00378   void check_pool_size(void);
00379   void check_derivatives(const dvector&,function_minimizer * pfmin,
00380     double fval1);
00381   dvector local_minimization_routine(dvector& s,dmatrix& Hess,
00382     dvector& grad,double lambda);
00383   dvector local_minimization(dvector& s,dmatrix& Hess,
00384     dvector& grad,double lambda);
00385   imatrix check_sparse_matrix_structure(void);
00386   double get_fx_fu(function_minimizer * pfmin);
00387   void do_separable_stuff_x_u_block_diagonal(df1b2variable& ff);
00388   void generate_antithetical_rvs();
00389   double standard_type3_loop(int no_converge_flag);
00390   void do_newton_raphson_state_space
00391     (function_minimizer * pfmin,double f_from_1,int& no_converge_flag);
00392   void begin_separable_call_stuff(void);
00393   void end_separable_call_stuff(void);
00394   void build_up_nested_shape(void);
00395 };
00396 
00400 class gauss_hermite_stuff
00401 {
00402 public:
00403   dvar_matrix gauss_hermite_values;
00404   dvector x;
00405   dvector w;
00406   int is;
00407   multi_index * mi;
00408 
00409   gauss_hermite_stuff(laplace_approximation_calculator * lapprox,
00410     int use_gauss_hermite,int num_separable_calls ,const ivector& itmp);
00411 
00412   friend class laplace_approximation_calculator;
00413 };
00414 
00415  ostream& operator <<  (const ostream &,const nested_calls_shape &);
00416 
00417 #endif  //  #if !defined(__RANDEFFECTS__)