1//
2// This file is part of Gambit
3// Copyright (c) 1994-2016, The Gambit Project (http://www.gambit-project.org)
4//
5// FILE: src/liblinear/btableau.imp
6// Implementation of base tableau classes
7//
8// This program is free software; you can redistribute it and/or modify
9// it under the terms of the GNU General Public License as published by
10// the Free Software Foundation; either version 2 of the License, or
11// (at your option) any later version.
12//
13// This program is distributed in the hope that it will be useful,
14// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
16// GNU General Public License for more details.
17//
18// You should have received a copy of the GNU General Public License
19// along with this program; if not, write to the Free Software
20// Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
21//
22
23// ---------------------------------------------------------------------------
24//                BaseTableau method definitions
25// ---------------------------------------------------------------------------
26
27#include "btableau.h"
28
29namespace Gambit {
30
31namespace linalg {
32
33template <class T>
34bool BaseTableau<T>::ColIndex(int x) const
35{ return MinCol()<=x && x<=MaxCol(); }
36
37template <class T>
38bool BaseTableau<T>::RowIndex(int x) const
39{ return MinRow()<=x && x<=MaxRow(); }
40
41template <class T>
42bool BaseTableau<T>::ValidIndex(int x) const
43{ return (ColIndex(x) || RowIndex(-x)); }
44
45template <class T>
46void BaseTableau<T>::CompPivot(int outlabel, int col)
47{
48  Pivot(Find(outlabel),col);
49  Pivot(Find(-col),-outlabel);
50}
51
52
53// ---------------------------------------------------------------------------
54//            TableauInterface method definitions
55// ---------------------------------------------------------------------------
56
57// Constructors and Destructor
58
59template <class T>
60TableauInterface<T>::TableauInterface(const Gambit::Matrix<T> &A, const Gambit::Vector<T> &b)
61  : A(&A), b(&b), basis(A.MinRow(),A.MaxRow(),A.MinCol(),A.MaxCol()),
62    solution(A.MinRow(),A.MaxRow()), npivots(0),
63    artificial(A.MaxCol()+1,A.MaxCol())
64
65{
66  // These are the values recommended by Murtagh (1981) for 15 digit
67  // accuracy in LP problems
68  // Note: for Rational, eps1 and eps2 resolve to 0
69  Gambit::linalg::epsilon(eps1, 5);
70  Gambit::linalg::epsilon(eps2);
71}
72
73template <class T>
74TableauInterface<T>::TableauInterface(const Gambit::Matrix<T> &A,
75				      const Gambit::Array<int> &art, const Gambit::Vector<T> &b)
76  : A(&A), b(&b),
77    basis(A.MinRow(),A.MaxRow(),A.MinCol(),A.MaxCol()+art.Length()),
78    solution(A.MinRow(),A.MaxRow()), npivots(0),
79    artificial(A.MaxCol()+1,A.MaxCol()+art.Length())
80{
81  Gambit::linalg::epsilon(eps1, 5);
82  Gambit::linalg::epsilon(eps2);
83  for(int i = 0;i<art.Length();i++)
84    artificial[A.MaxCol()+1+i] = art[art.First()+i];
85}
86
87template <class T>
88TableauInterface<T>::TableauInterface(const TableauInterface<T> &orig)
89  : A(orig.A), b(orig.b), basis(orig.basis), solution(orig.solution),
90    npivots(orig.npivots), eps1(orig.eps1), eps2(orig.eps2),
91    artificial(orig.artificial)
92{ }
93
94template <class T>
95TableauInterface<T>::~TableauInterface()
96{ }
97
98template <class T>
99TableauInterface<T>& TableauInterface<T>::operator=(const TableauInterface<T> &orig)
100{
101  if(this!= &orig) {
102    A = orig.A;
103    b = orig.b;
104    basis= orig.basis;
105    solution= orig.solution;
106    npivots = orig.npivots;
107    artificial = orig.artificial;
108  }
109  return *this;
110}
111
112// getting information
113
114template <class T>
115int TableauInterface<T>::MinRow() const { return A->MinRow(); }
116
117template <class T>
118int TableauInterface<T>::MaxRow() const { return A->MaxRow(); }
119
120template <class T>
121int TableauInterface<T>::MinCol() const { return basis.MinCol(); }
122
123template <class T>
124int TableauInterface<T>::MaxCol() const { return basis.MaxCol(); }
125
126template <class T>
127Basis & TableauInterface<T>::GetBasis(void) {return basis; }
128
129template <class T>
130const Gambit::Matrix<T> & TableauInterface<T>::Get_A(void) const {return *A; }
131
132template <class T>
133const Gambit::Vector<T> & TableauInterface<T>::Get_b(void) const {return *b;}
134
135template <class T>
136bool TableauInterface<T>::Member(int i) const
137{ return basis.Member(i);}
138
139template <class T>
140int TableauInterface<T>::Label(int i) const
141{ return basis.Label(i);}
142
143template <class T>
144int TableauInterface<T>::Find(int i) const
145{ return basis.Find(i);}
146
147template <class T>
148long TableauInterface<T>::NumPivots() const
149{ return npivots; }
150
151template <class T>
152long &TableauInterface<T>::NumPivots()
153{ return npivots; }
154
155template <class T>
156void TableauInterface<T>::Mark(int label)
157{basis.Mark(label);}
158
159template <class T>
160void TableauInterface<T>::UnMark(int label)
161{basis.UnMark(label);}
162
163template <class T>
164bool TableauInterface<T>::IsBlocked(int label) const
165{
166  return basis.IsBlocked(label);
167}
168
169template <class T>
170void TableauInterface<T>::GetColumn(int col, Gambit::Vector<T> &ret) const
171{
172  if(IsArtifColumn(col)) {
173    ret = (T) 0;
174    ret[artificial[col]] = (T)1;
175  }
176  else if(basis.IsRegColumn(col))
177    A->GetColumn(col, ret);
178  else if (basis.IsSlackColumn(col)) {
179    ret = (T) 0;
180    ret[-col] = (T) 1;
181  }
182}
183
184template <class T>
185void TableauInterface<T>::GetBasis(Basis &out) const
186{
187  out= basis;
188}
189
190template <class T>
191BFS<T> TableauInterface<T>::GetBFS()
192{
193  Gambit::Vector<T> sol(basis.First(),basis.Last());
194  BasisVector(sol);
195
196  BFS<T> cbfs;
197  for(int i=MinCol();i<=MaxCol();i++)  {
198    if(Member(i)) {
199      cbfs.insert(i,sol[basis.Find(i)]);
200    }
201  }
202  return cbfs;
203}
204
205template <class T>
206BFS<T> TableauInterface<T>::GetBFS1() const
207{
208  Gambit::Vector<T> sol(basis.First(),basis.Last());
209  BasisVector(sol);
210
211  BFS<T> cbfs;
212  int i;
213  for(i=-MaxRow();i<=-MinRow();i++)  {
214    if(Member(i)) {
215      cbfs.insert(i,sol[basis.Find(i)]);
216    }
217  }
218  for(i=MinCol();i<=MaxCol();i++)  {
219    if(Member(i)) {
220      cbfs.insert(i,sol[basis.Find(i)]);
221    }
222  }
223  return cbfs;
224}
225
226// miscellaneous functions
227
228template <class T>
229bool TableauInterface<T>::EqZero(T x) const
230{
231  return (LeZero(x) && GeZero(x));
232}
233
234template <class T>
235bool TableauInterface<T>::LtZero(T x) const
236{
237  return !GeZero(x);
238}
239
240template <class T>
241bool TableauInterface<T>::GtZero(T x) const
242{
243  return !LeZero(x);
244}
245
246template <class T>
247bool TableauInterface<T>::LeZero(T x) const
248{
249  if(x <=eps2) return 1;
250  return 0;
251}
252
253template <class T>
254bool TableauInterface<T>::GeZero(T x) const
255{
256  if(x >= -eps2) return 1;
257  return 0;
258}
259
260template <class T>
261T TableauInterface<T>::Epsilon(int i) const
262{
263  if(i!=1 && i!=2) {
264    throw Gambit::DimensionException();
265  }
266  if(i==1) return eps1;
267  return eps2;
268}
269
270template <class T>
271bool TableauInterface<T>::IsArtifColumn(int col) const
272{
273  return col>=artificial.First() && col<=artificial.Last();
274}
275
276}  // end namespace Gambit::linalg
277
278}  // end namespace Gambit
279
280