COIN-OR::LEMON - Graph Library

source: glpk-cmake/src/glpnpp01.c @ 1:c445c931472f

Last change on this file since 1:c445c931472f was 1:c445c931472f, checked in by Alpar Juttner <alpar@…>, 10 years ago

Import glpk-4.45

  • Generated files and doc/notes are removed
File size: 28.4 KB
Line 
1/* glpnpp01.c */
2
3/***********************************************************************
4*  This code is part of GLPK (GNU Linear Programming Kit).
5*
6*  Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008,
7*  2009, 2010 Andrew Makhorin, Department for Applied Informatics,
8*  Moscow Aviation Institute, Moscow, Russia. All rights reserved.
9*  E-mail: <mao@gnu.org>.
10*
11*  GLPK is free software: you can redistribute it and/or modify it
12*  under the terms of the GNU General Public License as published by
13*  the Free Software Foundation, either version 3 of the License, or
14*  (at your option) any later version.
15*
16*  GLPK is distributed in the hope that it will be useful, but WITHOUT
17*  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
18*  or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
19*  License for more details.
20*
21*  You should have received a copy of the GNU General Public License
22*  along with GLPK. If not, see <http://www.gnu.org/licenses/>.
23***********************************************************************/
24
25#include "glpnpp.h"
26
27NPP *npp_create_wksp(void)
28{     /* create LP/MIP preprocessor workspace */
29      NPP *npp;
30      npp = xmalloc(sizeof(NPP));
31      npp->orig_dir = 0;
32      npp->orig_m = npp->orig_n = npp->orig_nnz = 0;
33      npp->pool = dmp_create_pool();
34      npp->name = npp->obj = NULL;
35      npp->c0 = 0.0;
36      npp->nrows = npp->ncols = 0;
37      npp->r_head = npp->r_tail = NULL;
38      npp->c_head = npp->c_tail = NULL;
39      npp->stack = dmp_create_pool();
40      npp->top = NULL;
41#if 0 /* 16/XII-2009 */
42      memset(&npp->count, 0, sizeof(npp->count));
43#endif
44      npp->m = npp->n = npp->nnz = 0;
45      npp->row_ref = npp->col_ref = NULL;
46      npp->sol = npp->scaling = 0;
47      npp->p_stat = npp->d_stat = npp->t_stat = npp->i_stat = 0;
48      npp->r_stat = NULL;
49      /*npp->r_prim =*/ npp->r_pi = NULL;
50      npp->c_stat = NULL;
51      npp->c_value = /*npp->c_dual =*/ NULL;
52      return npp;
53}
54
55void npp_insert_row(NPP *npp, NPPROW *row, int where)
56{     /* insert row to the row list */
57      if (where == 0)
58      {  /* insert row to the beginning of the row list */
59         row->prev = NULL;
60         row->next = npp->r_head;
61         if (row->next == NULL)
62            npp->r_tail = row;
63         else
64            row->next->prev = row;
65         npp->r_head = row;
66      }
67      else
68      {  /* insert row to the end of the row list */
69         row->prev = npp->r_tail;
70         row->next = NULL;
71         if (row->prev == NULL)
72            npp->r_head = row;
73         else
74            row->prev->next = row;
75         npp->r_tail = row;
76      }
77      return;
78}
79
80void npp_remove_row(NPP *npp, NPPROW *row)
81{     /* remove row from the row list */
82      if (row->prev == NULL)
83         npp->r_head = row->next;
84      else
85         row->prev->next = row->next;
86      if (row->next == NULL)
87         npp->r_tail = row->prev;
88      else
89         row->next->prev = row->prev;
90      return;
91}
92
93void npp_activate_row(NPP *npp, NPPROW *row)
94{     /* make row active */
95      if (!row->temp)
96      {  row->temp = 1;
97         /* move the row to the beginning of the row list */
98         npp_remove_row(npp, row);
99         npp_insert_row(npp, row, 0);
100      }
101      return;
102}
103
104void npp_deactivate_row(NPP *npp, NPPROW *row)
105{     /* make row inactive */
106      if (row->temp)
107      {  row->temp = 0;
108         /* move the row to the end of the row list */
109         npp_remove_row(npp, row);
110         npp_insert_row(npp, row, 1);
111      }
112      return;
113}
114
115void npp_insert_col(NPP *npp, NPPCOL *col, int where)
116{     /* insert column to the column list */
117      if (where == 0)
118      {  /* insert column to the beginning of the column list */
119         col->prev = NULL;
120         col->next = npp->c_head;
121         if (col->next == NULL)
122            npp->c_tail = col;
123         else
124            col->next->prev = col;
125         npp->c_head = col;
126      }
127      else
128      {  /* insert column to the end of the column list */
129         col->prev = npp->c_tail;
130         col->next = NULL;
131         if (col->prev == NULL)
132            npp->c_head = col;
133         else
134            col->prev->next = col;
135         npp->c_tail = col;
136      }
137      return;
138}
139
140void npp_remove_col(NPP *npp, NPPCOL *col)
141{     /* remove column from the column list */
142      if (col->prev == NULL)
143         npp->c_head = col->next;
144      else
145         col->prev->next = col->next;
146      if (col->next == NULL)
147         npp->c_tail = col->prev;
148      else
149         col->next->prev = col->prev;
150      return;
151}
152
153void npp_activate_col(NPP *npp, NPPCOL *col)
154{     /* make column active */
155      if (!col->temp)
156      {  col->temp = 1;
157         /* move the column to the beginning of the column list */
158         npp_remove_col(npp, col);
159         npp_insert_col(npp, col, 0);
160      }
161      return;
162}
163
164void npp_deactivate_col(NPP *npp, NPPCOL *col)
165{     /* make column inactive */
166      if (col->temp)
167      {  col->temp = 0;
168         /* move the column to the end of the column list */
169         npp_remove_col(npp, col);
170         npp_insert_col(npp, col, 1);
171      }
172      return;
173}
174
175NPPROW *npp_add_row(NPP *npp)
176{     /* add new row to the current problem */
177      NPPROW *row;
178      row = dmp_get_atom(npp->pool, sizeof(NPPROW));
179      row->i = ++(npp->nrows);
180      row->name = NULL;
181      row->lb = -DBL_MAX, row->ub = +DBL_MAX;
182      row->ptr = NULL;
183      row->temp = 0;
184      npp_insert_row(npp, row, 1);
185      return row;
186}
187
188NPPCOL *npp_add_col(NPP *npp)
189{     /* add new column to the current problem */
190      NPPCOL *col;
191      col = dmp_get_atom(npp->pool, sizeof(NPPCOL));
192      col->j = ++(npp->ncols);
193      col->name = NULL;
194#if 0
195      col->kind = GLP_CV;
196#else
197      col->is_int = 0;
198#endif
199      col->lb = col->ub = col->coef = 0.0;
200      col->ptr = NULL;
201      col->temp = 0;
202      npp_insert_col(npp, col, 1);
203      return col;
204}
205
206NPPAIJ *npp_add_aij(NPP *npp, NPPROW *row, NPPCOL *col, double val)
207{     /* add new element to the constraint matrix */
208      NPPAIJ *aij;
209      aij = dmp_get_atom(npp->pool, sizeof(NPPAIJ));
210      aij->row = row;
211      aij->col = col;
212      aij->val = val;
213      aij->r_prev = NULL;
214      aij->r_next = row->ptr;
215      aij->c_prev = NULL;
216      aij->c_next = col->ptr;
217      if (aij->r_next != NULL)
218         aij->r_next->r_prev = aij;
219      if (aij->c_next != NULL)
220         aij->c_next->c_prev = aij;
221      row->ptr = col->ptr = aij;
222      return aij;
223}
224
225int npp_row_nnz(NPP *npp, NPPROW *row)
226{     /* count number of non-zero coefficients in row */
227      NPPAIJ *aij;
228      int nnz;
229      xassert(npp == npp);
230      nnz = 0;
231      for (aij = row->ptr; aij != NULL; aij = aij->r_next)
232         nnz++;
233      return nnz;
234}
235
236int npp_col_nnz(NPP *npp, NPPCOL *col)
237{     /* count number of non-zero coefficients in column */
238      NPPAIJ *aij;
239      int nnz;
240      xassert(npp == npp);
241      nnz = 0;
242      for (aij = col->ptr; aij != NULL; aij = aij->c_next)
243         nnz++;
244      return nnz;
245}
246
247void *npp_push_tse(NPP *npp, int (*func)(NPP *npp, void *info),
248      int size)
249{     /* push new entry to the transformation stack */
250      NPPTSE *tse;
251      tse = dmp_get_atom(npp->stack, sizeof(NPPTSE));
252      tse->func = func;
253      tse->info = dmp_get_atom(npp->stack, size);
254      tse->link = npp->top;
255      npp->top = tse;
256      return tse->info;
257}
258
259#if 1 /* 23/XII-2009 */
260void npp_erase_row(NPP *npp, NPPROW *row)
261{     /* erase row content to make it empty */
262      NPPAIJ *aij;
263      while (row->ptr != NULL)
264      {  aij = row->ptr;
265         row->ptr = aij->r_next;
266         if (aij->c_prev == NULL)
267            aij->col->ptr = aij->c_next;
268         else
269            aij->c_prev->c_next = aij->c_next;
270         if (aij->c_next == NULL)
271            ;
272         else
273            aij->c_next->c_prev = aij->c_prev;
274         dmp_free_atom(npp->pool, aij, sizeof(NPPAIJ));
275      }
276      return;
277}
278#endif
279
280void npp_del_row(NPP *npp, NPPROW *row)
281{     /* remove row from the current problem */
282#if 0 /* 23/XII-2009 */
283      NPPAIJ *aij;
284#endif
285      if (row->name != NULL)
286         dmp_free_atom(npp->pool, row->name, strlen(row->name)+1);
287#if 0 /* 23/XII-2009 */
288      while (row->ptr != NULL)
289      {  aij = row->ptr;
290         row->ptr = aij->r_next;
291         if (aij->c_prev == NULL)
292            aij->col->ptr = aij->c_next;
293         else
294            aij->c_prev->c_next = aij->c_next;
295         if (aij->c_next == NULL)
296            ;
297         else
298            aij->c_next->c_prev = aij->c_prev;
299         dmp_free_atom(npp->pool, aij, sizeof(NPPAIJ));
300      }
301#else
302      npp_erase_row(npp, row);
303#endif
304      npp_remove_row(npp, row);
305      dmp_free_atom(npp->pool, row, sizeof(NPPROW));
306      return;
307}
308
309void npp_del_col(NPP *npp, NPPCOL *col)
310{     /* remove column from the current problem */
311      NPPAIJ *aij;
312      if (col->name != NULL)
313         dmp_free_atom(npp->pool, col->name, strlen(col->name)+1);
314      while (col->ptr != NULL)
315      {  aij = col->ptr;
316         col->ptr = aij->c_next;
317         if (aij->r_prev == NULL)
318            aij->row->ptr = aij->r_next;
319         else
320            aij->r_prev->r_next = aij->r_next;
321         if (aij->r_next == NULL)
322            ;
323         else
324            aij->r_next->r_prev = aij->r_prev;
325         dmp_free_atom(npp->pool, aij, sizeof(NPPAIJ));
326      }
327      npp_remove_col(npp, col);
328      dmp_free_atom(npp->pool, col, sizeof(NPPCOL));
329      return;
330}
331
332void npp_del_aij(NPP *npp, NPPAIJ *aij)
333{     /* remove element from the constraint matrix */
334      if (aij->r_prev == NULL)
335         aij->row->ptr = aij->r_next;
336      else
337         aij->r_prev->r_next = aij->r_next;
338      if (aij->r_next == NULL)
339         ;
340      else
341         aij->r_next->r_prev = aij->r_prev;
342      if (aij->c_prev == NULL)
343         aij->col->ptr = aij->c_next;
344      else
345         aij->c_prev->c_next = aij->c_next;
346      if (aij->c_next == NULL)
347         ;
348      else
349         aij->c_next->c_prev = aij->c_prev;
350      dmp_free_atom(npp->pool, aij, sizeof(NPPAIJ));
351      return;
352}
353
354void npp_load_prob(NPP *npp, glp_prob *orig, int names, int sol,
355      int scaling)
356{     /* load original problem into the preprocessor workspace */
357      int m = orig->m;
358      int n = orig->n;
359      NPPROW **link;
360      int i, j;
361      double dir;
362      xassert(names == GLP_OFF || names == GLP_ON);
363      xassert(sol == GLP_SOL || sol == GLP_IPT || sol == GLP_MIP);
364      xassert(scaling == GLP_OFF || scaling == GLP_ON);
365      if (sol == GLP_MIP) xassert(!scaling);
366      npp->orig_dir = orig->dir;
367      if (npp->orig_dir == GLP_MIN)
368         dir = +1.0;
369      else if (npp->orig_dir == GLP_MAX)
370         dir = -1.0;
371      else
372         xassert(npp != npp);
373      npp->orig_m = m;
374      npp->orig_n = n;
375      npp->orig_nnz = orig->nnz;
376      if (names && orig->name != NULL)
377      {  npp->name = dmp_get_atom(npp->pool, strlen(orig->name)+1);
378         strcpy(npp->name, orig->name);
379      }
380      if (names && orig->obj != NULL)
381      {  npp->obj = dmp_get_atom(npp->pool, strlen(orig->obj)+1);
382         strcpy(npp->obj, orig->obj);
383      }
384      npp->c0 = dir * orig->c0;
385      /* load rows */
386      link = xcalloc(1+m, sizeof(NPPROW *));
387      for (i = 1; i <= m; i++)
388      {  GLPROW *rrr = orig->row[i];
389         NPPROW *row;
390         link[i] = row = npp_add_row(npp);
391         xassert(row->i == i);
392         if (names && rrr->name != NULL)
393         {  row->name = dmp_get_atom(npp->pool, strlen(rrr->name)+1);
394            strcpy(row->name, rrr->name);
395         }
396         if (!scaling)
397         {  if (rrr->type == GLP_FR)
398               row->lb = -DBL_MAX, row->ub = +DBL_MAX;
399            else if (rrr->type == GLP_LO)
400               row->lb = rrr->lb, row->ub = +DBL_MAX;
401            else if (rrr->type == GLP_UP)
402               row->lb = -DBL_MAX, row->ub = rrr->ub;
403            else if (rrr->type == GLP_DB)
404               row->lb = rrr->lb, row->ub = rrr->ub;
405            else if (rrr->type == GLP_FX)
406               row->lb = row->ub = rrr->lb;
407            else
408               xassert(rrr != rrr);
409         }
410         else
411         {  double rii = rrr->rii;
412            if (rrr->type == GLP_FR)
413               row->lb = -DBL_MAX, row->ub = +DBL_MAX;
414            else if (rrr->type == GLP_LO)
415               row->lb = rrr->lb * rii, row->ub = +DBL_MAX;
416            else if (rrr->type == GLP_UP)
417               row->lb = -DBL_MAX, row->ub = rrr->ub * rii;
418            else if (rrr->type == GLP_DB)
419               row->lb = rrr->lb * rii, row->ub = rrr->ub * rii;
420            else if (rrr->type == GLP_FX)
421               row->lb = row->ub = rrr->lb * rii;
422            else
423               xassert(rrr != rrr);
424         }
425      }
426      /* load columns and constraint coefficients */
427      for (j = 1; j <= n; j++)
428      {  GLPCOL *ccc = orig->col[j];
429         GLPAIJ *aaa;
430         NPPCOL *col;
431         col = npp_add_col(npp);
432         xassert(col->j == j);
433         if (names && ccc->name != NULL)
434         {  col->name = dmp_get_atom(npp->pool, strlen(ccc->name)+1);
435            strcpy(col->name, ccc->name);
436         }
437         if (sol == GLP_MIP)
438#if 0
439            col->kind = ccc->kind;
440#else
441            col->is_int = (char)(ccc->kind == GLP_IV);
442#endif
443         if (!scaling)
444         {  if (ccc->type == GLP_FR)
445               col->lb = -DBL_MAX, col->ub = +DBL_MAX;
446            else if (ccc->type == GLP_LO)
447               col->lb = ccc->lb, col->ub = +DBL_MAX;
448            else if (ccc->type == GLP_UP)
449               col->lb = -DBL_MAX, col->ub = ccc->ub;
450            else if (ccc->type == GLP_DB)
451               col->lb = ccc->lb, col->ub = ccc->ub;
452            else if (ccc->type == GLP_FX)
453               col->lb = col->ub = ccc->lb;
454            else
455               xassert(ccc != ccc);
456            col->coef = dir * ccc->coef;
457            for (aaa = ccc->ptr; aaa != NULL; aaa = aaa->c_next)
458               npp_add_aij(npp, link[aaa->row->i], col, aaa->val);
459         }
460         else
461         {  double sjj = ccc->sjj;
462            if (ccc->type == GLP_FR)
463               col->lb = -DBL_MAX, col->ub = +DBL_MAX;
464            else if (ccc->type == GLP_LO)
465               col->lb = ccc->lb / sjj, col->ub = +DBL_MAX;
466            else if (ccc->type == GLP_UP)
467               col->lb = -DBL_MAX, col->ub = ccc->ub / sjj;
468            else if (ccc->type == GLP_DB)
469               col->lb = ccc->lb / sjj, col->ub = ccc->ub / sjj;
470            else if (ccc->type == GLP_FX)
471               col->lb = col->ub = ccc->lb / sjj;
472            else
473               xassert(ccc != ccc);
474            col->coef = dir * ccc->coef * sjj;
475            for (aaa = ccc->ptr; aaa != NULL; aaa = aaa->c_next)
476               npp_add_aij(npp, link[aaa->row->i], col,
477                  aaa->row->rii * aaa->val * sjj);
478         }
479      }
480      xfree(link);
481      /* keep solution indicator and scaling option */
482      npp->sol = sol;
483      npp->scaling = scaling;
484      return;
485}
486
487void npp_build_prob(NPP *npp, glp_prob *prob)
488{     /* build resultant (preprocessed) problem */
489      NPPROW *row;
490      NPPCOL *col;
491      NPPAIJ *aij;
492      int i, j, type, len, *ind;
493      double dir, *val;
494      glp_erase_prob(prob);
495      glp_set_prob_name(prob, npp->name);
496      glp_set_obj_name(prob, npp->obj);
497      glp_set_obj_dir(prob, npp->orig_dir);
498      if (npp->orig_dir == GLP_MIN)
499         dir = +1.0;
500      else if (npp->orig_dir == GLP_MAX)
501         dir = -1.0;
502      else
503         xassert(npp != npp);
504      glp_set_obj_coef(prob, 0, dir * npp->c0);
505      /* build rows */
506      for (row = npp->r_head; row != NULL; row = row->next)
507      {  row->temp = i = glp_add_rows(prob, 1);
508         glp_set_row_name(prob, i, row->name);
509         if (row->lb == -DBL_MAX && row->ub == +DBL_MAX)
510            type = GLP_FR;
511         else if (row->ub == +DBL_MAX)
512            type = GLP_LO;
513         else if (row->lb == -DBL_MAX)
514            type = GLP_UP;
515         else if (row->lb != row->ub)
516            type = GLP_DB;
517         else
518            type = GLP_FX;
519         glp_set_row_bnds(prob, i, type, row->lb, row->ub);
520      }
521      /* build columns and the constraint matrix */
522      ind = xcalloc(1+prob->m, sizeof(int));
523      val = xcalloc(1+prob->m, sizeof(double));
524      for (col = npp->c_head; col != NULL; col = col->next)
525      {  j = glp_add_cols(prob, 1);
526         glp_set_col_name(prob, j, col->name);
527#if 0
528         glp_set_col_kind(prob, j, col->kind);
529#else
530         glp_set_col_kind(prob, j, col->is_int ? GLP_IV : GLP_CV);
531#endif
532         if (col->lb == -DBL_MAX && col->ub == +DBL_MAX)
533            type = GLP_FR;
534         else if (col->ub == +DBL_MAX)
535            type = GLP_LO;
536         else if (col->lb == -DBL_MAX)
537            type = GLP_UP;
538         else if (col->lb != col->ub)
539            type = GLP_DB;
540         else
541            type = GLP_FX;
542         glp_set_col_bnds(prob, j, type, col->lb, col->ub);
543         glp_set_obj_coef(prob, j, dir * col->coef);
544         len = 0;
545         for (aij = col->ptr; aij != NULL; aij = aij->c_next)
546         {  len++;
547            ind[len] = aij->row->temp;
548            val[len] = aij->val;
549         }
550         glp_set_mat_col(prob, j, len, ind, val);
551      }
552      xfree(ind);
553      xfree(val);
554      /* resultant problem has been built */
555      npp->m = prob->m;
556      npp->n = prob->n;
557      npp->nnz = prob->nnz;
558      npp->row_ref = xcalloc(1+npp->m, sizeof(int));
559      npp->col_ref = xcalloc(1+npp->n, sizeof(int));
560      for (row = npp->r_head, i = 0; row != NULL; row = row->next)
561         npp->row_ref[++i] = row->i;
562      for (col = npp->c_head, j = 0; col != NULL; col = col->next)
563         npp->col_ref[++j] = col->j;
564      /* transformed problem segment is no longer needed */
565      dmp_delete_pool(npp->pool), npp->pool = NULL;
566      npp->name = npp->obj = NULL;
567      npp->c0 = 0.0;
568      npp->r_head = npp->r_tail = NULL;
569      npp->c_head = npp->c_tail = NULL;
570      return;
571}
572
573void npp_postprocess(NPP *npp, glp_prob *prob)
574{     /* postprocess solution from the resultant problem */
575      GLPROW *row;
576      GLPCOL *col;
577      NPPTSE *tse;
578      int i, j, k;
579      double dir;
580      xassert(npp->orig_dir == prob->dir);
581      if (npp->orig_dir == GLP_MIN)
582         dir = +1.0;
583      else if (npp->orig_dir == GLP_MAX)
584         dir = -1.0;
585      else
586         xassert(npp != npp);
587      xassert(npp->m == prob->m);
588      xassert(npp->n == prob->n);
589      xassert(npp->nnz == prob->nnz);
590      /* copy solution status */
591      if (npp->sol == GLP_SOL)
592      {  npp->p_stat = prob->pbs_stat;
593         npp->d_stat = prob->dbs_stat;
594      }
595      else if (npp->sol == GLP_IPT)
596         npp->t_stat = prob->ipt_stat;
597      else if (npp->sol == GLP_MIP)
598         npp->i_stat = prob->mip_stat;
599      else
600         xassert(npp != npp);
601      /* allocate solution arrays */
602      if (npp->sol == GLP_SOL)
603      {  if (npp->r_stat == NULL)
604            npp->r_stat = xcalloc(1+npp->nrows, sizeof(char));
605         for (i = 1; i <= npp->nrows; i++)
606            npp->r_stat[i] = 0;
607         if (npp->c_stat == NULL)
608            npp->c_stat = xcalloc(1+npp->ncols, sizeof(char));
609         for (j = 1; j <= npp->ncols; j++)
610            npp->c_stat[j] = 0;
611      }
612#if 0
613      if (npp->r_prim == NULL)
614         npp->r_prim = xcalloc(1+npp->nrows, sizeof(double));
615      for (i = 1; i <= npp->nrows; i++)
616         npp->r_prim[i] = DBL_MAX;
617#endif
618      if (npp->c_value == NULL)
619         npp->c_value = xcalloc(1+npp->ncols, sizeof(double));
620      for (j = 1; j <= npp->ncols; j++)
621         npp->c_value[j] = DBL_MAX;
622      if (npp->sol != GLP_MIP)
623      {  if (npp->r_pi == NULL)
624            npp->r_pi = xcalloc(1+npp->nrows, sizeof(double));
625         for (i = 1; i <= npp->nrows; i++)
626            npp->r_pi[i] = DBL_MAX;
627#if 0
628         if (npp->c_dual == NULL)
629            npp->c_dual = xcalloc(1+npp->ncols, sizeof(double));
630         for (j = 1; j <= npp->ncols; j++)
631            npp->c_dual[j] = DBL_MAX;
632#endif
633      }
634      /* copy solution components from the resultant problem */
635      if (npp->sol == GLP_SOL)
636      {  for (i = 1; i <= npp->m; i++)
637         {  row = prob->row[i];
638            k = npp->row_ref[i];
639            npp->r_stat[k] = (char)row->stat;
640            /*npp->r_prim[k] = row->prim;*/
641            npp->r_pi[k] = dir * row->dual;
642         }
643         for (j = 1; j <= npp->n; j++)
644         {  col = prob->col[j];
645            k = npp->col_ref[j];
646            npp->c_stat[k] = (char)col->stat;
647            npp->c_value[k] = col->prim;
648            /*npp->c_dual[k] = dir * col->dual;*/
649         }
650      }
651      else if (npp->sol == GLP_IPT)
652      {  for (i = 1; i <= npp->m; i++)
653         {  row = prob->row[i];
654            k = npp->row_ref[i];
655            /*npp->r_prim[k] = row->pval;*/
656            npp->r_pi[k] = dir * row->dval;
657         }
658         for (j = 1; j <= npp->n; j++)
659         {  col = prob->col[j];
660            k = npp->col_ref[j];
661            npp->c_value[k] = col->pval;
662            /*npp->c_dual[k] = dir * col->dval;*/
663         }
664      }
665      else if (npp->sol == GLP_MIP)
666      {
667#if 0
668         for (i = 1; i <= npp->m; i++)
669         {  row = prob->row[i];
670            k = npp->row_ref[i];
671            /*npp->r_prim[k] = row->mipx;*/
672         }
673#endif
674         for (j = 1; j <= npp->n; j++)
675         {  col = prob->col[j];
676            k = npp->col_ref[j];
677            npp->c_value[k] = col->mipx;
678         }
679      }
680      else
681         xassert(npp != npp);
682      /* perform postprocessing to construct solution to the original
683         problem */
684      for (tse = npp->top; tse != NULL; tse = tse->link)
685      {  xassert(tse->func != NULL);
686         xassert(tse->func(npp, tse->info) == 0);
687      }
688      return;
689}
690
691void npp_unload_sol(NPP *npp, glp_prob *orig)
692{     /* store solution to the original problem */
693      GLPROW *row;
694      GLPCOL *col;
695      int i, j;
696      double dir;
697      xassert(npp->orig_dir == orig->dir);
698      if (npp->orig_dir == GLP_MIN)
699         dir = +1.0;
700      else if (npp->orig_dir == GLP_MAX)
701         dir = -1.0;
702      else
703         xassert(npp != npp);
704      xassert(npp->orig_m == orig->m);
705      xassert(npp->orig_n == orig->n);
706      xassert(npp->orig_nnz == orig->nnz);
707      if (npp->sol == GLP_SOL)
708      {  /* store basic solution */
709         orig->valid = 0;
710         orig->pbs_stat = npp->p_stat;
711         orig->dbs_stat = npp->d_stat;
712         orig->obj_val = orig->c0;
713         orig->some = 0;
714         for (i = 1; i <= orig->m; i++)
715         {  row = orig->row[i];
716            row->stat = npp->r_stat[i];
717            if (!npp->scaling)
718            {  /*row->prim = npp->r_prim[i];*/
719               row->dual = dir * npp->r_pi[i];
720            }
721            else
722            {  /*row->prim = npp->r_prim[i] / row->rii;*/
723               row->dual = dir * npp->r_pi[i] * row->rii;
724            }
725            if (row->stat == GLP_BS)
726               row->dual = 0.0;
727            else if (row->stat == GLP_NL)
728            {  xassert(row->type == GLP_LO || row->type == GLP_DB);
729               row->prim = row->lb;
730            }
731            else if (row->stat == GLP_NU)
732            {  xassert(row->type == GLP_UP || row->type == GLP_DB);
733               row->prim = row->ub;
734            }
735            else if (row->stat == GLP_NF)
736            {  xassert(row->type == GLP_FR);
737               row->prim = 0.0;
738            }
739            else if (row->stat == GLP_NS)
740            {  xassert(row->type == GLP_FX);
741               row->prim = row->lb;
742            }
743            else
744               xassert(row != row);
745         }
746         for (j = 1; j <= orig->n; j++)
747         {  col = orig->col[j];
748            col->stat = npp->c_stat[j];
749            if (!npp->scaling)
750            {  col->prim = npp->c_value[j];
751               /*col->dual = dir * npp->c_dual[j];*/
752            }
753            else
754            {  col->prim = npp->c_value[j] * col->sjj;
755               /*col->dual = dir * npp->c_dual[j] / col->sjj;*/
756            }
757            if (col->stat == GLP_BS)
758               col->dual = 0.0;
759#if 1
760            else if (col->stat == GLP_NL)
761            {  xassert(col->type == GLP_LO || col->type == GLP_DB);
762               col->prim = col->lb;
763            }
764            else if (col->stat == GLP_NU)
765            {  xassert(col->type == GLP_UP || col->type == GLP_DB);
766               col->prim = col->ub;
767            }
768            else if (col->stat == GLP_NF)
769            {  xassert(col->type == GLP_FR);
770               col->prim = 0.0;
771            }
772            else if (col->stat == GLP_NS)
773            {  xassert(col->type == GLP_FX);
774               col->prim = col->lb;
775            }
776            else
777               xassert(col != col);
778#endif
779            orig->obj_val += col->coef * col->prim;
780         }
781#if 1
782         /* compute primal values of inactive rows */
783         for (i = 1; i <= orig->m; i++)
784         {  row = orig->row[i];
785            if (row->stat == GLP_BS)
786            {  GLPAIJ *aij;
787               double temp;
788               temp = 0.0;
789               for (aij = row->ptr; aij != NULL; aij = aij->r_next)
790                  temp += aij->val * aij->col->prim;
791               row->prim = temp;
792            }
793         }
794         /* compute reduced costs of active columns */
795         for (j = 1; j <= orig->n; j++)
796         {  col = orig->col[j];
797            if (col->stat != GLP_BS)
798            {  GLPAIJ *aij;
799               double temp;
800               temp = col->coef;
801               for (aij = col->ptr; aij != NULL; aij = aij->c_next)
802                  temp -= aij->val * aij->row->dual;
803               col->dual = temp;
804            }
805         }
806#endif
807      }
808      else if (npp->sol == GLP_IPT)
809      {  /* store interior-point solution */
810         orig->ipt_stat = npp->t_stat;
811         orig->ipt_obj = orig->c0;
812         for (i = 1; i <= orig->m; i++)
813         {  row = orig->row[i];
814            if (!npp->scaling)
815            {  /*row->pval = npp->r_prim[i];*/
816               row->dval = dir * npp->r_pi[i];
817            }
818            else
819            {  /*row->pval = npp->r_prim[i] / row->rii;*/
820               row->dval = dir * npp->r_pi[i] * row->rii;
821            }
822         }
823         for (j = 1; j <= orig->n; j++)
824         {  col = orig->col[j];
825            if (!npp->scaling)
826            {  col->pval = npp->c_value[j];
827               /*col->dval = dir * npp->c_dual[j];*/
828            }
829            else
830            {  col->pval = npp->c_value[j] * col->sjj;
831               /*col->dval = dir * npp->c_dual[j] / col->sjj;*/
832            }
833            orig->ipt_obj += col->coef * col->pval;
834         }
835#if 1
836         /* compute row primal values */
837         for (i = 1; i <= orig->m; i++)
838         {  row = orig->row[i];
839            {  GLPAIJ *aij;
840               double temp;
841               temp = 0.0;
842               for (aij = row->ptr; aij != NULL; aij = aij->r_next)
843                  temp += aij->val * aij->col->pval;
844               row->pval = temp;
845            }
846         }
847         /* compute column dual values */
848         for (j = 1; j <= orig->n; j++)
849         {  col = orig->col[j];
850            {  GLPAIJ *aij;
851               double temp;
852               temp = col->coef;
853               for (aij = col->ptr; aij != NULL; aij = aij->c_next)
854                  temp -= aij->val * aij->row->dval;
855               col->dval = temp;
856            }
857         }
858#endif
859      }
860      else if (npp->sol == GLP_MIP)
861      {  /* store MIP solution */
862         xassert(!npp->scaling);
863         orig->mip_stat = npp->i_stat;
864         orig->mip_obj = orig->c0;
865#if 0
866         for (i = 1; i <= orig->m; i++)
867         {  row = orig->row[i];
868            /*row->mipx = npp->r_prim[i];*/
869         }
870#endif
871         for (j = 1; j <= orig->n; j++)
872         {  col = orig->col[j];
873            col->mipx = npp->c_value[j];
874            if (col->kind == GLP_IV)
875               xassert(col->mipx == floor(col->mipx));
876            orig->mip_obj += col->coef * col->mipx;
877         }
878#if 1
879         /* compute row primal values */
880         for (i = 1; i <= orig->m; i++)
881         {  row = orig->row[i];
882            {  GLPAIJ *aij;
883               double temp;
884               temp = 0.0;
885               for (aij = row->ptr; aij != NULL; aij = aij->r_next)
886                  temp += aij->val * aij->col->mipx;
887               row->mipx = temp;
888            }
889         }
890#endif
891      }
892      else
893         xassert(npp != npp);
894      return;
895}
896
897void npp_delete_wksp(NPP *npp)
898{     /* delete LP/MIP preprocessor workspace */
899      if (npp->pool != NULL)
900         dmp_delete_pool(npp->pool);
901      if (npp->stack != NULL)
902         dmp_delete_pool(npp->stack);
903      if (npp->row_ref != NULL)
904         xfree(npp->row_ref);
905      if (npp->col_ref != NULL)
906         xfree(npp->col_ref);
907      if (npp->r_stat != NULL)
908         xfree(npp->r_stat);
909#if 0
910      if (npp->r_prim != NULL)
911         xfree(npp->r_prim);
912#endif
913      if (npp->r_pi != NULL)
914         xfree(npp->r_pi);
915      if (npp->c_stat != NULL)
916         xfree(npp->c_stat);
917      if (npp->c_value != NULL)
918         xfree(npp->c_value);
919#if 0
920      if (npp->c_dual != NULL)
921         xfree(npp->c_dual);
922#endif
923      xfree(npp);
924      return;
925}
926
927/* eof */
Note: See TracBrowser for help on using the repository browser.