Revision 545 trunk/src/nh99/mod_hess.cpp
mod_hess.cpp (revision 545)  

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 
#if defined(USE_LAPLACE) 
8  8 
# include <df1b2fun.h> 
...  ...  
88  88 
if (ad_comm::wd_flag) 
89  89 
tmpstring = ad_comm::adprogram_name + ".hes"; 
90  90 
uostream ofs((char*)tmpstring); 
91 


91  
92  92 
ofs << nvar; 
93  93 
{ 
94  94 
{ 
...  ...  
138  138 
gradcalc(nvar,g2); 
139  139 
x(i)=xsave; 
140  140 
hess1=(g1g2)/(sdelta1sdelta2); 
141 


141  
142  142 
sdelta1=x(i)+eps*delta; 
143  143 
useless(sdelta1); 
144  144 
sdelta1=x(i); 
...  ...  
170  170 
double eps2=eps*eps; 
171  171 
hess2=(g1g2)/(sdelta1sdelta2); 
172  172 
hess=(eps2*hess1hess2) /(eps21.); 
173 


173  
174  174 
ofs << hess; 
175  175 
//if (adjm_ptr) ad_update_hess_stats_report(nvar,i); 
176  176 
} 
...  ...  
475  475 
dvector ggg(1,1); 
476  476 
gradcalc(0,ggg); 
477  477 
gradient_structure::set_YES_DERIVATIVES(); 
478 
initial_params::restore_start_phase();


478 
initial_params::restore_start_phase(); 

479  479 
int nvar=initial_params::nvarcalc(); // get the number of active parameters 
480  480 
int ndvar=stddev_params::num_stddev_calc(); 
481  481 
independent_variables x(1,nvar); 
...  ...  
492  492 
lapprox>no_function_component_flag=1; 
493  493 
} 
494  494 
#endif 
495 


495  
496  496 
dvariable vf; 
497  497 
vf=initial_params::reset(dvar_vector(x)); 
498  498 
*objective_function_value::pobjfun=0.0; 
...  ...  
523  523 
// symmetrize and invert the hessian 
524  524 
void function_minimizer::hess_inv(void) 
525  525 
{ 
526 
initial_params::set_inactive_only_random_effects();


526 
initial_params::set_inactive_only_random_effects(); 

527  527 
int nvar=initial_params::nvarcalc(); // get the number of active parameters 
528  528 
independent_variables x(1,nvar); 
529  529  
...  ...  
555  555 
ifs >> sscale; 
556  556 
if (!ifs) 
557  557 
{ 
558 
cerr << "Error reading sscale"


558 
cerr << "Error reading sscale" 

559  559 
<< " in routine hess_inv()" << endl; 
560  560 
} 
561  561  
...  ...  
584  584 
int zero_switch=0; 
585  585 
for (int j=1;j<=nvar;j++) 
586  586 
{ 
587 
if (hess(i,j)!=0.0)


587 
if (hess(i,j)!=0.0) 

588  588 
{ 
589  589 
zero_switch=1; 
590  590 
} 
...  ...  
601  601 
double llss=ln_det(hess,ssggnn); 
602  602 
int on1=0; 
603  603 
useless(llss); 
604 
{


604 
{ 

605  605 
ofstream ofs3((char*)(ad_comm::adprogram_name + adstring(".eva"))); 
606  606 
{ 
607 
dvector se=sort(eigenvalues(hess));


607 
dvector se=eigenvalues(hess);


608  608 
ofs3 << setshowpoint() << setw(14) << setprecision(10) 
609 
<< se << endl; 

610 
if (se(se.indexmin())<=0.0) 

609 
<< "unsorted:\t" << se << endl; 

610 
se=sort(se); 

611 
ofs3 << setshowpoint() << setw(14) << setprecision(10) 
Also available in: Unified diff