Annotation of OpenXM/src/ox_cdd/cddlib.c, Revision 1.2
1.2 ! noro 1: /* $OpenXM: OpenXM/src/ox_cdd/cddlib.c,v 1.1 2005/05/25 04:42:20 noro Exp $ */
1.1 noro 2:
3: /* This program is free software; you can redistribute it and/or modify
4: it under the terms of the GNU General Public License as published by
5: the Free Software Foundation; either version 2 of the License, or
6: (at your option) any later version.
7:
8: This program is distributed in the hope that it will be useful,
9: but WITHOUT ANY WARRANTY; without even the implied warranty of
10: MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11: GNU General Public License for more details.
12:
13: You should have received a copy of the GNU General Public License
14: along with this program; if not, write to the Free Software
15: Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
16: */
17:
18:
19: #include "setoper.h"
20: #include "cdd.h"
21: #include <stdio.h>
22: #include <stdlib.h>
23: #include <time.h>
24: #include <math.h>
25: #include <string.h>
26:
27: /* for gc */
28: #include <gc/gc.h>
29: #define MALLOC(x) GC_MALLOC((x))
30: #define MALLOC_ATOMIC(x) GC_MALLOC_ATOMIC((x))
31: #define ALLOCA(x) alloca((x))
32: /* #define FREE(x) free((x)) */
33: #define FREE(x)
34:
35: int debug_print = 1;
36: int mytype2int( mytype , int );
37:
38:
39: void init_cdd(void)
40: {
41: dd_set_global_constants(); /* First, this must be called. */
42: debug_print = 0;
43: }
44:
45: int **redcheck(int row,int col,int **matrix,int *presult_row)
46: {
47: dd_MatrixPtr M=NULL,M2=NULL,M3=NULL;
48: dd_ErrorType err=dd_NoError;
49: dd_rowset redrows,linrows,ignoredrows, basisrows;
50: dd_colset ignoredcols, basiscols;
51: long rank;
52:
53: time_t starttime, endtime;
54:
55: int i,j;
56: int **result=NULL;
57:
58: M = dd_CreateMatrix( row, col );
59: for(i=0;i<row;i++){
60: for(j=0;j<col;j++){
61: dd_set_si(M->matrix[i][j],matrix[i][j]);
62: }
63: }
64:
65: M->representation=dd_Generator;
66:
67: if(debug_print>=2){
68: printf("Inputed Matrix : \n");
69: dd_WriteMatrix(stdout, M);
70:
71: fprintf(stdout, "\nredundant rows:\n");
72: time(&starttime);
73: }
74:
75: redrows=dd_RedundantRows(M, &err);
76:
77: if(debug_print>=2){
78: time(&endtime);
79: set_fwrite(stdout, redrows);
80: dd_WriteTimes(stdout,starttime,endtime);
81: }
82:
83: M2=dd_MatrixSubmatrix(M, redrows);
84:
85: if(debug_print>=2){
86: printf("Implicit linearity (after removal of redundant rows): ");
87: }
88:
89: linrows=dd_ImplicitLinearityRows(M2, &err);
90:
91: if(debug_print>=2){
92: fprintf(stdout," %ld ", set_card(linrows));
93: set_fwrite(stdout,linrows);
94: }
95: set_uni(M2->linset, M2->linset, linrows);
96: /* add the implicit linrows to the explicit linearity rows */
97:
98: if(debug_print>=2){
99: printf("\nNonredundant representation (except possibly for the linearity part):\n");
100: dd_WriteMatrix(stdout, M2);
101: }
102:
103: /* To remove redundancy of the linearity part,
104: we need to compute the rank and a basis of the linearity part. */
105: set_initialize(&ignoredrows, M2->rowsize);
106: set_initialize(&ignoredcols, M2->colsize);
107: set_compl(ignoredrows, M2->linset);
108: rank=dd_MatrixRank(M2,ignoredrows,ignoredcols, &basisrows, &basiscols);
109: set_diff(ignoredrows, M2->linset, basisrows);
110: M3=dd_MatrixSubmatrix(M2, ignoredrows);
111: if (rank>0){
112: if(debug_print>=2){
113: fprintf(stdout,"\nThe following %ld linearity rows are dependent and unnecessary:", set_card(ignoredrows));
114: set_fwrite(stdout,ignoredrows);
115: }
116: } else{
117: if(debug_print>=2){
118: fprintf(stdout,"\nThe linearity rows are independent and thus minimal\n");
119: }
120: }
121:
122: if(debug_print>=2){
123: printf("Nonredundant representation (= minimal representation):\n");
124: dd_WriteMatrix(stdout, M3);
125: }
126:
127: *presult_row = M3->rowsize;
128:
129: result = MALLOC( *presult_row * sizeof(int*) );
130: for(i=0;i<*presult_row;i++){
131: result[i] = MALLOC( col * sizeof(int) );
132: }
133:
134: for(i=0;i<*presult_row;i++){
135: for(j=0;j<col;j++){
136: result[i][j] = mytype2int( M3->matrix[i][j], 0 );
137: }
138: }
139:
140: set_free(linrows);
141: set_free(basisrows);
142: set_free(basiscols);
143: set_free(ignoredrows);
144: set_free(ignoredcols);
145: set_free(redrows);
146:
147: dd_FreeMatrix(M);
148: dd_FreeMatrix(M2);
149: dd_FreeMatrix(M3);
150:
151: if (err!=dd_NoError) dd_WriteErrorMessages(stderr,err);
152: return result;
153: }
154:
155: mytype *lpsolve(dd_LPObjectiveType type,int row,int col,int **matrix,int *obj)
156: {
157: dd_ErrorType error=dd_NoError;
158: dd_LPSolverType solver; /* either DualSimplex or CrissCross */
159: dd_LPPtr lp;
160:
161: dd_MatrixPtr A;
162: dd_ErrorType err;
163: int i,j;
164: static mytype result;
165:
166: A = dd_CreateMatrix( row, col );
167: for(i=0;i<row;i++){
168: for(j=0;j<col;j++){
169: dd_set_si(A->matrix[i][j],matrix[i][j]);
170: }
171: }
172:
173: for(j=0;j<col;j++){
174: dd_set_si(A->rowvec[j],obj[j]);
175: }
176:
177: A->objective = type;
178: lp=dd_Matrix2LP(A, &err); /* load an LP */
179: if (lp==NULL) goto _L99;
180:
181: if(debug_print>=2){
182: /* Print the LP. */
183: printf("\n--- LP to be solved ---\n");
184: dd_WriteLP(stdout, lp);
185:
186: /* Solve the LP by cdd LP solver. */
187: printf("\n--- Running dd_LPSolve ---\n");
188: }
189:
190: solver=dd_DualSimplex;
191: dd_LPSolve(lp, solver, &error); /* Solve the LP */
192: if (error!=dd_NoError) goto _L99;
193:
194: if(debug_print>=2){
195: /* Write the LP solutions by cdd LP reporter. */
196: dd_WriteLPResult(stdout, lp, error);
197: }
198:
199: dd_init( result );
200: dd_set( result, lp->optvalue );
201:
202: /* Free allocated spaces. */
203: dd_FreeLPData(lp);
204: dd_FreeMatrix(A);
205:
206: _L99:
207: if (error!=dd_NoError) dd_WriteErrorMessages(stdout, error);
208: return &result;
209: }
1.2 ! noro 210:
! 211: mytype *lpintpt(int row,int col,int **matrix,mytype *ptp)
! 212: {
! 213: dd_ErrorType error=dd_NoError;
! 214: dd_rowset ImL, Lbasis;
! 215: dd_LPSolutionPtr lps1;
! 216:
! 217: dd_MatrixPtr A;
! 218: dd_ErrorType err;
! 219: int i,j;
! 220: mytype *ret;
! 221:
! 222: A = dd_CreateMatrix( row, col );
! 223: for(i=0;i<row;i++){
! 224: for(j=0;j<col;j++){
! 225: dd_set_si(A->matrix[i][j],matrix[i][j]);
! 226: }
! 227: }
! 228:
! 229: dd_FindRelativeInterior(A,&ImL,&Lbasis,&lps1,&error);
! 230: if (error!=dd_NoError) goto _L99;
! 231: if ( dd_Positive(lps1->optvalue) ) {
! 232: ret = (mytype *)MALLOC(lps1->d*sizeof(mytype));
! 233: for ( i = 0; i < lps1->d; i++ ) {
! 234: dd_init(ret[i]);
! 235: dd_set(ret[i],lps1->sol[i]);
! 236: }
! 237: } else {
! 238: ret = 0;
! 239: }
! 240: dd_FreeLPSolution(lps1);
! 241: set_free(ImL);
! 242: set_free(Lbasis);
! 243: _L99:
! 244: if (error!=dd_NoError) dd_WriteErrorMessages(stdout, error);
! 245: return ret;
! 246: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>