Revision 1107 trunk/src/df1b2separable/test_trust.cpp
test_trust.cpp (revision 1107)  

2  2 
* $Id$ 
3  3 
* 
4  4 
* Author: David Fournier 
5 
* Copyright (c) 20082012 Regents of the University of California


5 
* Copyright (c) 20082012 Regents of the University of California 

6  6 
*/ 
7  7 
/** 
8  8 
* \file 
...  ...  
20  20 
double calculate_laplace_approximation(const dvector& x,const dvector& u0, 
21  21 
const dmatrix& Hess,const dvector& _xadjoint,const dvector& _uadjoint, 
22  22 
const dmatrix& _Hessadjoint,function_minimizer * pmin); 
23 


23  
24  24 
double calculate_importance_sample(const dvector& x,const dvector& u0, 
25  25 
const dmatrix& Hess,const dvector& _xadjoint,const dvector& _uadjoint, 
26  26 
const dmatrix& _Hessadjoint,function_minimizer * pmin); 
27 


27  
28  28 
dmatrix choleski_decomp_positive(const dmatrix& M,double b); 
29  29  
30  30 
//dvector laplace_approximation_calculator::default_calculations 
...  ...  
43  43 
//int i,j; 
44  44 
int i; 
45  45  
46 
initial_params::set_inactive_only_random_effects();


46 
initial_params::set_inactive_only_random_effects(); 

47  47 
gradient_structure::set_NO_DERIVATIVES(); 
48  48 
initial_params::reset(x); // get current x values into the model 
49  49 
gradient_structure::set_YES_DERIVATIVES(); 
50  50  
51  51 
/* 
52 
initial_params::set_active_only_random_effects();


52 
initial_params::set_active_only_random_effects(); 

53  53 
int lmn_flag=0; 
54  54 
if (ad_comm::time_flag) 
55  55 
{ 
...  ...  
61  61 
{ 
62  62 
ad_comm::ptm>get_elapsed_time_and_reset(); 
63  63 
} 
64 
}


64 
} 

65  65 
if (ad_comm::time_flag) 
66  66 
{ 
67  67 
if (ad_comm::ptm) 
...  ...  
74  74 
} 
75  75 
} 
76  76 
} 
77 


77  
78  78 
double maxg; 
79  79 
double maxg_save; 
80  80 
dvector uhat_old(1,usize); 
...  ...  
138  138 
//int niters=0; 
139  139 
dvector g(1,n); 
140  140 
dmatrix H(1,n,1,n); 
141 


141  
142  142 
// test newton raphson 
143  143 
get_complete_hessian(H,g,pfmin); 
144  144  
145 


146  145 
double tol=1.e6; 
147  146 
//int itmax=1000; 
148  147 
//int itol=1; 
149  148 
//int iter=0; 
150  149 
//double err=0; 
151  150 
//lambda=1; 
152 


151  
153  152 
//cout << "input Delta" << endl; 
154  153 
//cin >> Delta; 
155  154 
//cout << "input lambda" << endl; 
156  155 
//cin >> lambda; 
157 


156  
158  157 
int i; 
159 


158  
160  159 
for (i=1;i<=n;i++) 
161  160 
{ 
162  161 
H(i,i)+=lambda; 
163  162 
} 
164 


163  
165  164 
cout << "initial fun value  " << double(0.5)*xx*(H*xx)+xx*g; 
166  165 
double truef=do_one_feval(xx,pfmin); 
167  166 
cout << "real function value = " << truef << endl; 
...  ...  
173  172 
int maxfn=15; 
174  173 
dvector xret=lincg(xx,g,H,tol,Delta,pfmin,truef,estdiff, 
175  174 
truediff,bestf,iflag,inner_iter,maxfn); 
Also available in: Unified diff