00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #include <cstdio>
00011 #include <math.h>
00012 #include "CoinTime.hpp"
00013
00014 #include "misc_util.hpp"
00015
00016 Stat::Stat() { reset(); }
00017 Stat::Stat(const double *vector,int n) {
00018 reset();
00019 for(int i=0;i<n;i++)
00020 addEntry(vector[i]);
00021 }
00022 Stat::Stat(const int *vector,int n) {
00023 reset();
00024 for(int i=0;i<n;i++)
00025 addEntry((double)vector[i]);
00026 }
00027 void Stat::reset() { _n = 0; _nz_entries = 0; _sum = 0; _sum_squares = 0;}
00028 void Stat::addEntry(double value) {
00029 if (fabs(value) > 1e-8)
00030 _nz_entries++;
00031 _sum += value;
00032 _sum_squares += value*value;
00033 _n++;
00034 if (_n == 1) {
00035 _min = value;
00036 _max = value;
00037 _minIndex = 1;
00038 _maxIndex = 1;
00039 } else {
00040 if (value < _min) {
00041 _min = value;
00042 _minIndex = _n;
00043 }
00044 if (value > _max) {
00045 _max = value;
00046 _maxIndex = _n;
00047 }
00048 }
00049 }
00050 int Stat::numEntries() const {return _n;}
00051 int Stat::numNZEntries() const {return _nz_entries;}
00052 double Stat::mean() const {
00053 if (_n>0)
00054 return _sum/_n;
00055 else
00056 return 0.0;
00057 }
00058 double Stat::stdDev() const {
00059 if (_n>0)
00060 return sqrt( fabs((_sum_squares / _n) - (mean()*mean())) );
00061 else
00062 return 0.0;
00063 }
00064 double Stat::sum() const {return _sum;}
00065 double Stat::min() const {
00066 if (_n>0)
00067 return _min;
00068 else
00069 return 0.0;
00070 }
00071 double Stat::max() const {
00072 if (_n>0)
00073 return _max;
00074 else
00075 return 0.0;
00076 }
00077 int Stat::minIndex() const {
00078 if (_n>0)
00079 return _minIndex;
00080 else
00081 return -1;
00082 }
00083 int Stat::maxIndex() const {
00084 if (_n>0)
00085 return _maxIndex;
00086 else
00087 return -1;
00088 }
00089
00090 Timer::Timer() { _pause = false; _starttime = -99999999;}
00091 Timer::~Timer() { if (_pause) delete _pausetimer; }
00092 void Timer::start() {
00093 _starttime = CoinCpuTime ();
00094 }
00095
00096 double Timer::time() {
00097 if (_starttime == -99999999)
00098 return 0.0;
00099 if (_pause) {
00100 return fabs( _pausetimer->starttime() - _starttime );
00101 }
00102
00103 return (fabs (CoinCpuTime () - _starttime));
00104 }
00105
00106 void Timer::pause() {
00107 if (_pause) {
00108 fprintf(stderr,"timer already in pause state\n");
00109 fprintf(stdout,"timer already in pause state\n");
00110 } else {
00111 _pause = true;
00112 _pausetimer = new Timer();
00113 _pausetimer->start();
00114 }
00115 }
00116 void Timer::restore() {
00117 if (_pause) {
00118 _pause = false;
00119 _starttime += _pausetimer->time();
00120 delete _pausetimer;
00121 } else {
00122 fprintf(stderr,"timer not in pause state\n");
00123 fprintf(stdout,"timer not in pause state\n");
00124 }
00125 }
00126
00127 double Timer::starttime() {
00128 return _starttime;
00129 }
00130
00131 void cpp_fprintvecINT(FILE *file,
00132 char const *vecstr, const int *x, const int n,
00133 const int numberAcross)
00134 {
00135 int num, fromto, upto, j, i;
00136
00137 num = (n/numberAcross) + 1;
00138 if(vecstr != NULL) {
00139 fprintf(file, "%s : \n", vecstr);
00140 }
00141 for(j=0; j<num; j++){
00142 fromto = numberAcross * j;
00143 upto = numberAcross * (j+1);
00144 if(n <= upto) upto = n;
00145 for(i=fromto; i<upto; i++)
00146 fprintf(file, " %4d", x[i]);
00147 fprintf(file, "\n");
00148 }
00149 fprintf(file, "\n");
00150 }
00151
00152
00153 void cpp_fprintvecINT(FILE *file,
00154 char const *vecstr, const int *x, const int n)
00155 {
00156 cpp_fprintvecINT(file, vecstr, x, n, 10);
00157 }
00158
00159
00160 void cpp_printvecINT(char const *vecstr, const int *x, const int n)
00161 {
00162 cpp_fprintvecINT(stdout, vecstr, x, n, 10);
00163 }
00164
00165
00166 void cpp_fprintmatINT(FILE *file, char const *vecstr, const int **x,
00167 const int m,
00168 const int n)
00169 {
00170 int i, j;
00171
00172 fprintf(file, "%s :\n", vecstr);
00173
00174 for(i=0; i<m; i++){
00175 for(j=0; j<n; j++){
00176 fprintf(file, " %4d", x[i][j]);
00177 }
00178 fprintf(file, "\n");
00179 }
00180 fprintf(file, "\n");
00181 }
00182
00183
00184 void cpp_fprintmatINT(FILE *file, char const *vecstr, int **x, const int m,
00185 const int n)
00186 {
00187 int i, j;
00188
00189 fprintf(file, "%s :\n", vecstr);
00190
00191 for(i=0; i<m; i++){
00192 for(j=0; j<n; j++){
00193 fprintf(file, " %4d", x[i][j]);
00194 }
00195 fprintf(file, "\n");
00196 }
00197 fprintf(file, "\n");
00198 }
00199
00200
00201 void cpp_printmatINT(char const *vecstr, const int **x, const int m, const int n)
00202 {
00203 int i, j;
00204
00205 printf("%s :\n", vecstr);
00206
00207 for(i=0; i<m; i++){
00208 for(j=0; j<n; j++){
00209 printf(" %4d", x[i][j]);
00210 }
00211 printf("\n");
00212 }
00213 printf("\n");
00214 }
00215
00216
00217 void cpp_printmatINT(char const *vecstr, int **x, const int m, const int n)
00218 {
00219 int i, j;
00220
00221 printf("%s :\n", vecstr);
00222
00223 for(i=0; i<m; i++){
00224 for(j=0; j<n; j++){
00225 printf(" %4d", x[i][j]);
00226 }
00227 printf("\n");
00228 }
00229 printf("\n");
00230 }
00231
00232
00233 void cpp_printvecCHAR(char const *vecstr, const char *x, const int n)
00234 {
00235 int num, fromto, upto, j, i;
00236
00237 num = (n/10) + 1;
00238 printf("%s :\n", vecstr);
00239 for(j=0; j<num; j++){
00240 fromto = 10 * j;
00241 upto = 10 * (j+1);
00242 if(n <= upto) upto = n;
00243 for(i = fromto; i<upto; i++)
00244 printf(" %c", x[i]);
00245 printf("\n");
00246 }
00247 printf("\n");
00248 }
00249
00250
00251 void cpp_printvecCHAR(char const *vecstr, char *x, const int n)
00252 {
00253 int num, fromto, upto, j, i;
00254
00255 num = (n/10) + 1;
00256 printf("%s :\n", vecstr);
00257 for(j=0; j<num; j++){
00258 fromto = 10 * j;
00259 upto = 10 * (j+1);
00260 if(n <= upto) upto = n;
00261 for(i = fromto; i<upto; i++)
00262 printf(" %c", x[i]);
00263 printf("\n");
00264 }
00265 printf("\n");
00266 }
00267
00268
00269 void cpp_fprintvecCHAR(FILE *file, char const *vecstr, const char *x, const int n)
00270 {
00271 int num, fromto, upto, j, i;
00272
00273 num = n/10 + 1;
00274 fprintf(file, "%s :\n", vecstr);
00275 for(j=0; j<num; j++){
00276 fromto = 10 * j;
00277 upto = 10 * (j+1);
00278 if(n <= upto) upto = n;
00279 for(i=fromto; i<upto; i++)
00280 fprintf(file, " %c", x[i]);
00281 fprintf(file, "\n");
00282 }
00283 fprintf(file, "\n");
00284 }
00285
00286
00287 void cpp_fprintvecCHAR(FILE *file, char const *vecstr, char *x, const int n)
00288 {
00289 int num, fromto, upto, j, i;
00290
00291 num = n/10 + 1;
00292 fprintf(file, "%s :\n", vecstr);
00293 for(j=0; j<num; j++){
00294 fromto = 10 * j;
00295 upto = 10 * (j+1);
00296 if(n <= upto) upto = n;
00297 for(i=fromto; i<upto; i++)
00298 fprintf(file, " %c", x[i]);
00299 fprintf(file, "\n");
00300 }
00301 fprintf(file, "\n");
00302 }
00303
00304
00305 void cpp_fprintvecDBL(FILE *file, char const *vecstr, const double *x,
00306 const int n, const int numberAcross,
00307 char *form)
00308 {
00309 int num, fromto, upto, j, i;
00310
00311 num = (n/numberAcross) + 1;
00312 if(vecstr != NULL) {
00313 fprintf(file,"%s:\n",vecstr);
00314 }
00315 for(j=0;j<num; j++){
00316 fromto = numberAcross * j;
00317 upto = numberAcross * (j + 1);
00318 if(n <= upto) upto = n;
00319 for(i = fromto; i < upto; i++)
00320 fprintf(file, form, x[i]);
00321 fprintf(file,"\n");
00322 }
00323 fprintf(file,"\n");
00324 }
00325
00326
00327 void cpp_fprintvecDBL(FILE *file, char const *vecstr, const double *x,
00328 const int n, const int numberAcross,
00329 const int print_pos, const int decimals)
00330 {
00331 char form[15];
00332 sprintf(form, " %%%d.%df", print_pos, decimals);
00333 cpp_fprintvecDBL(file, vecstr, x, n, numberAcross, form);
00334 }
00335
00336
00337 void cpp_fprintvecDBLg(FILE *file, char const *vecstr, const double *x,
00338 const int n, const int numberAcross,
00339 const int print_pos, const int print_digits)
00340 {
00341 char form[15];
00342 sprintf(form, " %%%d.%dg", print_pos, print_digits);
00343 cpp_fprintvecDBL(file, vecstr, x, n, numberAcross, form);
00344 }
00345
00346
00347 void cpp_fprintvecDBL(FILE *file, char const *vecstr, const double *x, const int n)
00348 {
00349 cpp_fprintvecDBL(file, vecstr, x, n, 10, 7, 3);
00350 }
00351
00352
00353 void cpp_printvecDBL(char const *vecstr, const double *x, const int n)
00354 {
00355 cpp_fprintvecDBL(stdout, vecstr, x, n, 10, 7, 3);
00356 }
00357
00358
00359 void cpp_printmatDBL(char const *vecstr, const double **x, const int m, const int n)
00360 {
00361 int i, j;
00362
00363 printf("%s :\n", vecstr);
00364
00365 for(i=0; i<m; i++){
00366 for(j=0; j<n; j++){
00367 printf(" %4.3f", x[i][j]);
00368 }
00369 printf("\n");
00370 }
00371 printf("\n");
00372 }
00373
00374
00375 void cpp_printmatDBL(char const *vecstr, double **x, const int m, const int n)
00376 {
00377 int i, j;
00378
00379 printf("%s :\n", vecstr);
00380
00381 for(i=0; i<m; i++){
00382 for(j=0; j<n; j++){
00383 printf(" %4.3f", x[i][j]);
00384 }
00385 printf("\n");
00386 }
00387 printf("\n");
00388 }
00389
00390
00391 void cpp_fprintmatDBL(FILE *file, char const *vecstr, const double **x, const int m,
00392 const int n)
00393 {
00394 int i, j;
00395
00396 fprintf(file, "%s :\n", vecstr);
00397
00398 for(i=0; i<m; i++){
00399 for(j=0; j<n; j++){
00400 fprintf(file, " %4.3f", x[i][j]);
00401 }
00402 fprintf(file, "\n");
00403 }
00404 fprintf(file, "\n");
00405 }
00406
00407
00408 void cpp_fprintmatDBL(FILE *file, char const *vecstr, double **x, const int m,
00409 const int n)
00410 {
00411 int i, j;
00412
00413 fprintf(file, "%s :\n", vecstr);
00414
00415 for(i=0; i<m; i++){
00416 for(j=0; j<n; j++){
00417 fprintf(file, " %4.3f", x[i][j]);
00418 }
00419 fprintf(file, "\n");
00420 }
00421 fprintf(file, "\n");
00422 }
00423
00424
00425
00426
00427
00428
00429
00430
00431
00432
00433
00434
00435
00436 double cpp_genalea (int *x0)
00437 {
00438 int m = 2147483647;
00439 int a = 16807 ;
00440 int b = 127773 ;
00441 int c = 2836 ;
00442 int x1, k;
00443
00444 k = (int) ((*x0)/b) ;
00445 x1 = a*(*x0 - k*b) - k*c ;
00446 if(x1 < 0) x1 = x1 + m;
00447 *x0 = x1;
00448
00449 return((double)x1/(double)m);
00450
00451
00452
00453
00454
00455 }
00456
00457
00458
00459 int cpp_partition_dec(const int q, const int r, int *key, double *value)
00460 {
00461 double v, t;
00462 int i, j, index;
00463
00464 v = value[r];
00465 i = q-1;
00466 j = r;
00467 do{
00468 i++;j--;
00469 while ( value[i]> v) i++;
00470 while ( (value[j]< v)&&(j>q) ) j--;
00471 t = value[i];
00472 value[i] = value[j];
00473 value[j] = t;
00474 index = key[i];
00475 key[i] = key[j];
00476 key[j] = index;
00477 }while (j>i);
00478
00479 value[j] = value[i];
00480 value[i] = value[r];
00481 value[r] = t;
00482 key[j] = key[i];
00483 key[i] = key[r];
00484 key[r] = index;
00485 return(i);
00486 }
00487
00488
00489 void cpp_quicksort_dec(const int k, const int l, int *key, double *value)
00490
00491 {
00492 int p, i, q, r;
00493 int stack[1000];
00494
00495 if(l == 1) {
00496 return;
00497 }
00498
00499 stack[0] = stack[1] = 0;
00500 q = k;
00501 r = l-1;
00502 p = 2;
00503 do{
00504 if (r>q){
00505 i = cpp_partition_dec(q,r,key,value);
00506 if (i-q>r-i){
00507 stack[p] = q;
00508 stack[p+1] = i - 1;
00509 q = i + 1;
00510 }
00511 else{
00512 stack[p] = i + 1;
00513 stack[p+1] = r;
00514 r = i - 1;
00515 }
00516 p += 2;
00517 }
00518 else{
00519 p -= 2;
00520 q = stack[p];
00521 r = stack[p+1];
00522 }
00523 }while (p);
00524 }
00525
00526
00527 int cpp_partitionINT_dec(const int q, const int r, int *key, int *value)
00528 {
00529 int v, t;
00530 int i, j, index;
00531
00532 v = value[r];
00533 i = q-1;
00534 j = r;
00535 do{
00536 i++;j--;
00537 while ( value[i]> v) i++;
00538 while ( (value[j]< v)&&(j>q) ) j--;
00539 t = value[i];
00540 value[i] = value[j];
00541 value[j] = t;
00542 index = key[i];
00543 key[i] = key[j];
00544 key[j] = index;
00545 }while (j>i);
00546
00547 value[j] = value[i];
00548 value[i] = value[r];
00549 value[r] = t;
00550 key[j] = key[i];
00551 key[i] = key[r];
00552 key[r] = index;
00553 return(i);
00554 }
00555
00556
00557 void cpp_quicksortINT_dec(const int k, const int l, int *key, int *value)
00558
00559 {
00560 int p, i, q, r;
00561 int stack[1000];
00562
00563 if(l == 1) {
00564 return;
00565 }
00566
00567 stack[0] = stack[1] = 0;
00568 q = k;
00569 r = l-1;
00570 p = 2;
00571 do{
00572 if (r>q){
00573 i = cpp_partitionINT_dec(q,r,key,value);
00574 if (i-q>r-i){
00575 stack[p] = q;
00576 stack[p+1] = i - 1;
00577 q = i + 1;
00578 }
00579 else{
00580 stack[p] = i + 1;
00581 stack[p+1] = r;
00582 r = i - 1;
00583 }
00584 p += 2;
00585 }
00586 else{
00587 p -= 2;
00588 q = stack[p];
00589 r = stack[p+1];
00590 }
00591 }while (p);
00592 }
00593
00594
00595 int cpp_partition_inc(const int q, const int r, int *key, double *value)
00596 {
00597 double v, t;
00598 int i, j, index;
00599
00600 v = value[r];
00601 i = q-1;
00602 j = r;
00603 do{
00604 i++;j--;
00605 while ( value[i]< v) i++;
00606 while ( (value[j]> v)&&(j>q) ) j--;
00607 t = value[i];
00608 value[i] = value[j];
00609 value[j] = t;
00610 index = key[i];
00611 key[i] = key[j];
00612 key[j] = index;
00613 }while (j>i);
00614
00615 value[j] = value[i];
00616 value[i] = value[r];
00617 value[r] = t;
00618 key[j] = key[i];
00619 key[i] = key[r];
00620 key[r] = index;
00621 return(i);
00622 }
00623
00624
00625 void cpp_quicksort_inc(const int k, const int l, int *key, double *value)
00626
00627 {
00628 int p, i, q, r;
00629 int stack[1000];
00630
00631
00632 if(l == 1) {
00633 return;
00634 }
00635
00636 stack[0] = stack[1] = 0;
00637 q = k;
00638 r = l-1;
00639 p = 2;
00640 do{
00641 if (r>q){
00642 i = cpp_partition_inc(q,r,key,value);
00643 if (i-q>r-i){
00644 stack[p] = q;
00645 stack[p+1] = i - 1;
00646 q = i + 1;
00647 }
00648 else{
00649 stack[p] = i + 1;
00650 stack[p+1] = r;
00651 r = i - 1;
00652 }
00653 p += 2;
00654 }
00655 else{
00656 p -= 2;
00657 q = stack[p];
00658 r = stack[p+1];
00659 }
00660 }while (p);
00661 }
00662
00663
00664 int cpp_partitionINT_inc(const int q, const int r, int *key, int *value)
00665 {
00666 int v, t;
00667 int i, j, index;
00668
00669 v = value[r];
00670 i = q-1;
00671 j = r;
00672 do{
00673 i++;j--;
00674 while ( value[i]< v) i++;
00675 while ( (value[j]> v)&&(j>q) ) j--;
00676 t = value[i];
00677 value[i] = value[j];
00678 value[j] = t;
00679 index = key[i];
00680 key[i] = key[j];
00681 key[j] = index;
00682 }while (j>i);
00683
00684 value[j] = value[i];
00685 value[i] = value[r];
00686 value[r] = t;
00687 key[j] = key[i];
00688 key[i] = key[r];
00689 key[r] = index;
00690 return(i);
00691 }
00692
00693
00694 void cpp_quicksortINT_inc(const int k, const int l, int *key, int *value)
00695
00696 {
00697 int p, i, q, r;
00698 int stack[1000];
00699
00700
00701 if(l == 1) {
00702 return;
00703 }
00704
00705 stack[0] = stack[1] = 0;
00706 q = k;
00707 r = l-1;
00708 p = 2;
00709 do{
00710 if (r>q){
00711 i = cpp_partitionINT_inc(q,r,key,value);
00712 if (i-q>r-i){
00713 stack[p] = q;
00714 stack[p+1] = i - 1;
00715 q = i + 1;
00716 }
00717 else{
00718 stack[p] = i + 1;
00719 stack[p+1] = r;
00720 r = i - 1;
00721 }
00722 p += 2;
00723 }
00724 else{
00725 p -= 2;
00726 q = stack[p];
00727 r = stack[p+1];
00728 }
00729 }while (p);
00730 }
00731
00732
00733
00734
00735
00736 void get_barQ(double **Q, const double *b, const int n, double **mat) {
00737
00738 int np = n+1;
00739 mat[0][0] = 0;
00740
00741 for(int i=1; i<np; i++) {
00742 mat[0][i] = b[i-1];
00743 mat[i][0] = b[i-1];
00744 }
00745
00746 for(int i=1; i<np; i++) {
00747 for(int j=i; j<np; j++) {
00748 mat[i][j] = Q[i-1][j-1];
00749 mat[j][i] = Q[i-1][j-1];
00750 }
00751 }
00752 }
00753
00754
00755 void print_barQ(FILE *f, char *header,
00756 double **Q, const double *b, const int n) {
00757
00758 int i, np = n+1;
00759 double **mat = (double **) malloc (np * sizeof (double *));
00760
00761 for (i=0; i<np; i++) {
00762 mat[i] = (double *) malloc (np * sizeof (double));
00763 }
00764
00765 get_barQ(Q, b, n, mat);
00766 cpp_fprintmatDBL(f, header, mat, np, np);
00767
00768 for (i=0; i<np; i++) {
00769 free (mat[i]);
00770 }
00771 free (mat);
00772
00773
00774 }
00775
00776
00777 void get_vec_from_matbar(double **mat, const int dim, double *v,
00778 const int include_entry_00) {
00779
00780 int ind = 0;
00781
00782 if(include_entry_00) {
00783 v[ind] = mat[0][0];
00784 ind++;
00785 }
00786
00787 for(int j=1; j<dim; j++) {
00788 v[ind] = mat[0][j];
00789 ind++;
00790 }
00791 for(int i=1; i<dim; i++) {
00792 for(int j=i; j<dim; j++) {
00793 v[ind] = mat[i][j];
00794 ind++;
00795 }
00796 }
00797 }
00798
00799
00800 void get_mat_from_vec(const double *v, const int n, const double entry_00,
00801 double **mat) {
00802
00803 mat[0][0] = entry_00;
00804
00805 for(int i=1; i<n+1; i++) {
00806 mat[0][i] = v[i-1];
00807 mat[i][0] = v[i-1];
00808 }
00809
00810 for(int i=1; i<n+1; i++) {
00811 for(int j=i; j<n+1; j++) {
00812 int ind = indexQ(i-1, j-1, n);
00813 mat[i][j] = v[ind];
00814 mat[j][i] = v[ind];
00815 }
00816 }
00817 }
00818
00819
00820 void print_mat_from_vec(FILE *f, char *header,
00821 const double *v, const int n, const double entry_00) {
00822
00823 int i, np = n+1;
00824 double **mat = (double **) malloc (np * sizeof (double *));
00825
00826 for (i=0; i<np; i++) {
00827 mat[i] = (double *) malloc (np * sizeof (double));
00828 }
00829
00830 get_mat_from_vec(v, n, entry_00, mat);
00831 cpp_fprintmatDBL(f, header, mat, np, np);
00832
00833 for (i=0; i<np; i++) {
00834 free (mat[i]);
00835 }
00836 free (mat);
00837 }
00838
00839
00840 void get_LPsol_vec_from_vvT(const double *v, const int n, double *vec,
00841 const int include_entry_00) {
00842
00843 int i, j, np = n+1;
00844 double **mat = (double **) malloc (np * sizeof (double *));
00845
00846 for (i=0; i<np; i++) {
00847 mat[i] = (double *) malloc (np * sizeof (double));
00848 }
00849
00850 mat[0][0] = 1;
00851
00852 for(i=1; i<np; i++) {
00853 mat[0][i] = v[i-1];
00854 mat[i][0] = v[i-1];
00855 }
00856
00857 for(i=1; i<np; i++) {
00858 for(j=1; j<np; j++) {
00859 mat[i][j] = v[i-1] * v[j-1];
00860 }
00861 }
00862 get_vec_from_matbar(mat, np, vec, include_entry_00);
00863
00864 for (i=0; i<np; i++) {
00865 free(mat [i]);
00866 }
00867 free(mat);
00868 }
00869
00870
00871 void get_mat_from_vvT(const double *v, const int n, double **mat) {
00872
00873 for(int i=0; i<n+1; i++) {
00874 for(int j=0; j<n+1; j++) {
00875 mat[i][j] = v[i] * v[j];
00876 }
00877 }
00878 }
00879
00880
00881 void print_mat_from_vvT(FILE *f, char *header,
00882 const double *v, const int n) {
00883
00884 int i, np = n+1;
00885 double **mat = (double **) malloc (np * sizeof (double *));
00886
00887 for (i=0; i<np; i++) {
00888 mat[i] = (double *) malloc (np * sizeof (double));
00889 }
00890
00891 get_mat_from_vvT(v, n, mat);
00892 cpp_fprintmatDBL(f, header, mat, np, np);
00893
00894 for (i=0; i<np; i++) {
00895 free (mat[i]);
00896 }
00897 free (mat);
00898 }
00899
00900
00901 void fprintvecmat(FILE *f, char *header, const double *v, const int n,
00902 const double entry_00) {
00903
00904 int i, np = n+1;
00905 double **mat = (double **) malloc (np * sizeof (double *));
00906
00907 for (i=0; i<np; i++) {
00908 mat[i] = (double *) malloc (np * sizeof (double));
00909 }
00910
00911 get_mat_from_vec(v, n, entry_00, mat);
00912 cpp_fprintmatDBL(f, header, mat, np, np);
00913
00914 for (i=0; i<np; i++) {
00915 free (mat[i]);
00916 }
00917 free (mat);
00918
00919 }
00920
00921
00922
00923 void check_prod_row_barQ(double **Q, const double *b, const int n,
00924 const OsiSolverInterface *si,
00925 const int from) {
00926 int i, j, k, np = n+1;
00927 double **mat = (double **) malloc (np * sizeof (double *));
00928 double **mat_row = (double **) malloc (np * sizeof (double *));
00929
00930 for (i=0; i<np; i++) {
00931 mat[i] = (double *) malloc (np * sizeof (double));
00932 mat_row[i] = (double *) malloc (np * sizeof (double));
00933 }
00934
00935 get_barQ(Q, b, n, mat);
00936
00937 int curr_nrows = si->getNumRows();
00938 const double *y = si->getRowPrice();
00939 const CoinPackedMatrix *byRow = si->getMatrixByRow();
00940 const double *rhs = si->getRightHandSide();
00941
00942 int row_size = np * np;
00943
00944 double *row_v = new double[row_size];
00945
00946 for(int i = from; i<curr_nrows; i++) {
00947 if(fabs(y[i]) > 1e-5) {
00948
00949 const int *rowind = byRow->getVector(i).getIndices();
00950 int card_row = byRow->getVector(i).getNumElements();
00951 const double *rowelem = byRow->getVector(i).getElements();
00952
00953 for(j=0; j<row_size; j++) {
00954 row_v[j] = 0;
00955 }
00956
00957 for(j=0; j<card_row; j++) {
00958 row_v[rowind[j]] = rowelem[j];
00959 }
00960
00961 get_mat_from_vec(row_v, n, -rhs[i], mat_row);
00962
00963 double ckmuQ = 0;
00964
00965 for (j=0; j<np; j++) {
00966 for (k=0; k<np; k++) {
00967 ckmuQ += mat[j][j] * mat_row[j][k];
00968 }
00969 }
00970
00971 if(ckmuQ > 0) {
00972 printf("### WARNING: SdpCutGen::check_prod_row_barQ(): row %d with dual: %8.4f muQ: %8.4f\n", i, y[i], ckmuQ);
00973
00974 #ifdef TRACE_ALL
00975 cpp_printmatDBL("mat_row", mat_row, np, np);
00976 cpp_printmatDBL("mat_barQ", mat, np, np);
00977 #endif
00978
00979 exit(1);
00980 }
00981 if(ckmuQ > 1e30) {
00982 exit(1);
00983 }
00984 }
00985 }
00986
00987 for (i=0; i<np; i++) {
00988 free (mat[i]);
00989 free (mat_row[i]);
00990 }
00991 free (mat);
00992 free (mat_row);
00993 delete[] row_v;
00994
00995 }
00996
00997