Revision 692 branches/mergetrunkdavef/src/sparse/hs_sparse.cpp
hs_sparse.cpp (revision 692)  

245  245 
return tmp; 
246  246 
} 
247  247  
248 
dvar_matrix make_sdvar_matrix(dvar_compressed_triplet& M,int n,int m) 

249 
{ 

250 
dvar_matrix tmp(1,n,1,m); 

251 
int nz = M.indexmax() M.indexmin() + 1; 

252 
for (int i=1;i<=nz;i++) 

253 
{ 

254 
tmp(M(1,i),M(2,i))=M(i); 

255 
if (M(1,i) != M(2,i)) 

256 
{ 

257 
tmp(M(2,i),M(1,i))=M(i); 

258 
} 

259 
} 

260 
return tmp; 

261 
} 

248  262  
249  263 
hs_smatrix::hs_smatrix(int _n, XCONST dcompressed_triplet &_M) 
250  264 
{ 
...  ...  
281  295 
for (k = 0 ; k < nz ; k++) 
282  296 
lower_tri += Ti[k]>Tj[k]; 
283  297 
if(lower_tri) 
298 
{ 

284  299 
cout << "hs_smatrix::hs_smatrix: M must be upper triangular" << endl; 
300 
ad_exit(1); 

301 
} 

285  302  
286  303 
// Matrix in compressed format 
287  304 
p.allocate(0,n); 
...  ...  
2565  2582 
{ 
2566  2583 
allocate(mmin,mmax,_n,_m); 
2567  2584 
} 
2585 
dcompressed_triplet::dcompressed_triplet(void) 

2586 
{ 

2587 
n=0; 

2588 
m=1; 

2589 
} 

2568  2590  
2569  
2570  2591 
void dvar_compressed_triplet::allocate(int mmin,int mmax,int _n,int _m) 
2571  2592 
{ 
2572  2593 
n=_n; 
...  ...  
2809  2830 
return x; 
2810  2831 
} 
2811  2832 

2833 
dvector solve(const dcompressed_triplet & _st,const dvector& _grad) 

2834 
{ 

2835 
ADUNCONST(dcompressed_triplet,st) 

2836 
ADUNCONST(dvector,grad) 

2837 
int n=st.get_n(); 

2838 
//int n=Hess.indexmax(); 

2839 
// fill up compressed triplet with nonzero entries of the Hessian 

2840 
hs_symbolic S(st,1); // Fill reducing rowcol permutation 

2841 


2842 
hs_smatrix HS(n,st); // Convert triplet to working format 

2843 


2844 
hs_smatrix L(S); // Allocates cholesky factor 

2845 


2846 
ivector nxcount; 

2847 
chol(HS,S,L); // Does numeric factorization 

2848 


2849 
dvector x(0,n1); 

2850 
grad.shift(0); 

2851 
x = cs_ipvec(S.pinv, grad); 

2852 
grad.shift(1); 

2853 
x = cs_lsolve(L,x); 

2854 
x = cs_ltsolve(L,x); 

2855 
x = cs_pvec(S.pinv,x); 

2856 
x.shift(1); 

2857 
return x; 

2858 
} 

2859 


2812  2860 
dvector solve(const dcompressed_triplet & _st,const dvector& _grad,const hs_symbolic& S,int& ierr) 
2813  2861 
{ 
2814  2862 
ADUNCONST(dcompressed_triplet,st) 
...  ...  
2857  2905 
{ 
2858  2906 
int n=VM.get_n(); 
2859  2907 
dvar_hs_smatrix H(n,VM); 
2860 
hs_symbolic S(VM,1); // Fill reducing rowcol permutation 

2908 
cout << "Needs fixing ??? " << endl; 

2909 
//sleep(1); 

2910 
// what is difference between 0 and 1 

2911 
hs_symbolic S(VM,0); // Fill reducing rowcol permutation 

2912 
//hs_symbolic S(VM,1); // Fill reducing rowcol permutation 

2861  2913 
dvar_hs_smatrix L(S); // Allocates cholesky factor 
2862  2914 
int ierr=chol(H,S,L); // Does numeric factorization 
2863  2915 
if (ierr==0) 
...  ...  
3048  3100  
3049  3101 
static void dfcholeski_sparse(void); 
3050  3102 
Also available in: Unified diff