Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TQpDataDens.cxx
Go to the documentation of this file.
1// @(#)root/quadp:$Id$
2// Author: Eddy Offermann May 2004
3
4/*************************************************************************
5 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12/*************************************************************************
13 * Parts of this file are copied from the OOQP distribution and *
14 * are subject to the following license: *
15 * *
16 * COPYRIGHT 2001 UNIVERSITY OF CHICAGO *
17 * *
18 * The copyright holder hereby grants you royalty-free rights to use, *
19 * reproduce, prepare derivative works, and to redistribute this software*
20 * to others, provided that any changes are clearly documented. This *
21 * software was authored by: *
22 * *
23 * E. MICHAEL GERTZ gertz@mcs.anl.gov *
24 * Mathematics and Computer Science Division *
25 * Argonne National Laboratory *
26 * 9700 S. Cass Avenue *
27 * Argonne, IL 60439-4844 *
28 * *
29 * STEPHEN J. WRIGHT swright@cs.wisc.edu *
30 * Computer Sciences Department *
31 * University of Wisconsin *
32 * 1210 West Dayton Street *
33 * Madison, WI 53706 FAX: (608)262-9777 *
34 * *
35 * Any questions or comments may be directed to one of the authors. *
36 * *
37 * ARGONNE NATIONAL LABORATORY (ANL), WITH FACILITIES IN THE STATES OF *
38 * ILLINOIS AND IDAHO, IS OWNED BY THE UNITED STATES GOVERNMENT, AND *
39 * OPERATED BY THE UNIVERSITY OF CHICAGO UNDER PROVISION OF A CONTRACT *
40 * WITH THE DEPARTMENT OF ENERGY. *
41 *************************************************************************/
42
43#include "TQpDataDens.h"
44
45////////////////////////////////////////////////////////////////////////////////
46///
47/// \class TQpDataDens
48///
49/// Data for the dense QP formulation
50///
51////////////////////////////////////////////////////////////////////////////////
52
54
55////////////////////////////////////////////////////////////////////////////////
56/// Constructor
57
59: TQpDataBase(nx,my,mz)
60{
64}
65
66
67////////////////////////////////////////////////////////////////////////////////
68/// Constructor
69
71 TVectorD &xlow_in,TVectorD &ixlow_in,
72 TVectorD &xupp_in,TVectorD &ixupp_in,
73 TMatrixD &A_in, TVectorD &bA_in,
74 TMatrixD &C_in,
75 TVectorD &clow_in,TVectorD &iclow_in,
76 TVectorD &cupp_in,TVectorD &icupp_in)
77{
78 fG .ResizeTo(c_in) ; fG = c_in;
79 fBa .ResizeTo(bA_in) ; fBa = bA_in;
80 fXloBound.ResizeTo(xlow_in) ; fXloBound = xlow_in;
81 fXloIndex.ResizeTo(ixlow_in); fXloIndex = ixlow_in;
82 fXupBound.ResizeTo(xupp_in) ; fXupBound = xupp_in;
83 fXupIndex.ResizeTo(ixupp_in); fXupIndex = ixupp_in;
84 fCloBound.ResizeTo(clow_in) ; fCloBound = clow_in;
85 fCloIndex.ResizeTo(iclow_in); fCloIndex = iclow_in;
86 fCupBound.ResizeTo(cupp_in) ; fCupBound = cupp_in;
87 fCupIndex.ResizeTo(icupp_in); fCupIndex = icupp_in;
88
89 fNx = fG.GetNrows();
90 fQ.Use(Q_in);
91
92 if (A_in.GetNrows() > 0) {
93 fA.Use(A_in);
94 fMy = fA.GetNrows();
95 } else
96 fMy = 0;
97
98 if (C_in.GetNrows() > 0) {
99 fC.Use(C_in);
100 fMz = fC.GetNrows();
101 } else
102 fMz = 0;
103}
104
105
106////////////////////////////////////////////////////////////////////////////////
107/// Copy constructor
108
110{
111 *this = another;
112}
113
114
115////////////////////////////////////////////////////////////////////////////////
116/// calculate y = beta*y + alpha*(fQ*x)
117
119{
120 y *= beta;
121 if (fQ.GetNoElements() > 0)
122 y += alpha*(fQ*x);
123}
124
125
126////////////////////////////////////////////////////////////////////////////////
127/// calculate y = beta*y + alpha*(fA*x)
128
130{
131 y *= beta;
132 if (fA.GetNoElements() > 0)
133 y += alpha*(fA*x);
134}
135
136
137////////////////////////////////////////////////////////////////////////////////
138/// calculate y = beta*y + alpha*(fC*x)
139
141{
142 y *= beta;
143 if (fC.GetNoElements() > 0)
144 y += alpha*(fC*x);
145}
146
147
148////////////////////////////////////////////////////////////////////////////////
149/// calculate y = beta*y + alpha*(fA^T*x)
150
152{
153 y *= beta;
154 if (fA.GetNoElements() > 0)
155 y += alpha*(TMatrixD(TMatrixD::kTransposed,fA)*x);
156}
157
158
159////////////////////////////////////////////////////////////////////////////////
160/// calculate y = beta*y + alpha*(fC^T*x)
161
163{
164 y *= beta;
165 if (fC.GetNoElements() > 0)
166 y += alpha*(TMatrixD(TMatrixD::kTransposed,fC)*x);
167}
168
169
170////////////////////////////////////////////////////////////////////////////////
171/// Return the largest component of several vectors in the data class
172
174{
175 Double_t norm = 0.0;
176
177 Double_t componentNorm = fG.NormInf();
178 if (componentNorm > norm) norm = componentNorm;
179
180 TMatrixDSym fQ_abs(fQ);
181 componentNorm = (fQ_abs.Abs()).Max();
182 if (componentNorm > norm) norm = componentNorm;
183
184 componentNorm = fBa.NormInf();
185 if (componentNorm > norm) norm = componentNorm;
186
187 TMatrixD fA_abs(fQ);
188 componentNorm = (fA_abs.Abs()).Max();
189 if (componentNorm > norm) norm = componentNorm;
190
191 TMatrixD fC_abs(fQ);
192 componentNorm = (fC_abs.Abs()).Max();
193 if (componentNorm > norm) norm = componentNorm;
194
196 componentNorm = fXloBound.NormInf();
197 if (componentNorm > norm) norm = componentNorm;
198
200 componentNorm = fXupBound.NormInf();
201 if (componentNorm > norm) norm = componentNorm;
202
204 componentNorm = fCloBound.NormInf();
205 if (componentNorm > norm) norm = componentNorm;
206
208 componentNorm = fCupBound.NormInf();
209 if (componentNorm > norm) norm = componentNorm;
210
211 return norm;
212}
213
214
215////////////////////////////////////////////////////////////////////////////////
216/// Print all class members
217
218void TQpDataDens::Print(Option_t * /*opt*/) const
219{
220 fQ.Print("Q");
221 fG.Print("c");
222
223 fXloBound.Print("xlow");
224 fXloIndex.Print("ixlow");
225
226 fXupBound.Print("xupp");
227 fXupIndex.Print("ixupp");
228
229 fA.Print("A");
230 fBa.Print("b");
231 fC.Print("C");
232
233 fCloBound.Print("clow");
234 fCloIndex.Print("iclow");
235
236 fCupBound.Print("cupp");
237 fCupIndex.Print("icupp");
238}
239
240
241////////////////////////////////////////////////////////////////////////////////
242/// Insert the Hessian Q into the matrix M at index (row,col) for the fundamental
243/// linear system
244
246{
247 m.SetSub(row,col,fQ);
248}
249
250
251////////////////////////////////////////////////////////////////////////////////
252/// Insert the constraint matrix A into the matrix M at index (row,col) for the fundamental
253/// linear system
254
256{
257 m.SetSub(row,col,fA);
258}
259
260
261////////////////////////////////////////////////////////////////////////////////
262/// Insert the constraint matrix C into the matrix M at index (row,col) for the fundamental
263/// linear system
264
266{
267 m.SetSub(row,col,fC);
268}
269
270
271////////////////////////////////////////////////////////////////////////////////
272/// Return in vector dq the diagonal of matrix fQ (Quadratic part of Objective function)
273
275{
276 const Int_t n = TMath::Min(fQ.GetNrows(),fQ.GetNcols());
277 dq.ResizeTo(n);
278 dq = TMatrixDDiag(fQ);
279}
280
281
282////////////////////////////////////////////////////////////////////////////////
283/// Return value of the objective function
284
286{
287 TVectorD tmp(fG);
288 this->Qmult(1.0,tmp,0.5,vars->fX);
289
290 return tmp*vars->fX;
291}
292
293
294////////////////////////////////////////////////////////////////////////////////
295/// Choose randomly a QP problem
296
298{
299 Double_t ix = 3074.20374;
300
301 TVectorD xdual(fNx);
303 TVectorD sprime(fMz);
305
306 fQ.RandomizePD(0.0,1.0,ix);
307 fA.Randomize(-10.0,10.0,ix);
308 fC.Randomize(-10.0,10.0,ix);
309 y .Randomize(-10.0,10.0,ix);
310
311 fG = xdual;
312 fG -= fQ*x;
313
316
317 fBa = fA*x;
318 s = fC*x;
319
320 // Now compute the real q = s-sprime
321 const TVectorD q = s-sprime;
322
323 // Adjust fCloBound and fCupBound appropriately
324 Add(fCloBound,1.0,q);
325 Add(fCupBound,1.0,q);
326
329}
330
331
332////////////////////////////////////////////////////////////////////////////////
333/// Assignment operator
334
336{
337 if (this != &source) {
339 fQ.ResizeTo(source.fQ); fQ = source.fQ;
340 fA.ResizeTo(source.fA); fA = source.fA;
341 fC.ResizeTo(source.fC); fC = source.fC;
342 }
343 return *this;
344}
const char Option_t
Definition RtypesCore.h:66
#define ClassImp(name)
Definition Rtypes.h:377
#define R__ASSERT(e)
Definition TError.h:118
float * q
TMatrixTDiag< Double_t > TMatrixDDiag
TMatrixT< Double_t > TMatrixD
Definition TMatrixDfwd.h:23
virtual TMatrixTBase< Element > & Randomize(Element alpha, Element beta, Double_t &seed)
Randomize matrix element values.
Int_t GetNrows() const
void Print(Option_t *name="") const override
Print the matrix as a table of elements.
Int_t GetNoElements() const
Int_t GetNcols() const
virtual TMatrixTBase< Element > & Abs()
Take an absolute value of a matrix, i.e. apply Abs() to each element.
TMatrixTBase< Element > & ResizeTo(Int_t nrows, Int_t ncols, Int_t=-1) override
Set size of the matrix to nrows x ncols New dynamic elements are created, the overlapping part of the...
virtual TMatrixTSym< Element > & RandomizePD(Element alpha, Element beta, Double_t &seed)
randomize matrix element values but keep matrix symmetric positive definite
TMatrixTSym< Element > & Use(Int_t row_lwb, Int_t row_upb, Element *data)
TMatrixT< Element > & Use(Int_t row_lwb, Int_t row_upb, Int_t col_lwb, Int_t col_upb, Element *data)
Use the array data to fill the matrix ([row_lwb..row_upb] x [col_lwb..col_upb])
TMatrixTBase< Element > & ResizeTo(Int_t nrows, Int_t ncols, Int_t=-1) override
Set size of the matrix to nrows x ncols New dynamic elements are created, the overlapping part of the...
Data for the general QP formulation.
Definition TQpDataBase.h:61
TQpDataBase & operator=(const TQpDataBase &source)
Assignment operator.
TVectorD fXupIndex
Definition TQpDataBase.h:80
TVectorD fXloBound
Definition TQpDataBase.h:81
TVectorD fCloIndex
Definition TQpDataBase.h:86
TVectorD fCupIndex
Definition TQpDataBase.h:84
TVectorD fG
Definition TQpDataBase.h:77
TVectorD fXupBound
Definition TQpDataBase.h:79
TVectorD fXloIndex
Definition TQpDataBase.h:82
TVectorD fCupBound
Definition TQpDataBase.h:83
TVectorD fCloBound
Definition TQpDataBase.h:85
TVectorD fBa
Definition TQpDataBase.h:78
static void RandomlyChooseBoundedVariables(TVectorD &x, TVectorD &dualx, TVectorD &blx, TVectorD &ixlow, TVectorD &bux, TVectorD &ixupp, Double_t &ix, Double_t percentLowerOnly, Double_t percentUpperOnly, Double_t percentBound)
Randomly choose x and its boundaries.
Data for the dense QP formulation.
Definition TQpDataDens.h:63
Double_t DataNorm() override
Return the largest component of several vectors in the data class.
void Qmult(Double_t beta, TVectorD &y, Double_t alpha, const TVectorD &x) override
calculate y = beta*y + alpha*(fQ*x)
void Cmult(Double_t beta, TVectorD &y, Double_t alpha, const TVectorD &x) override
calculate y = beta*y + alpha*(fC*x)
void GetDiagonalOfQ(TVectorD &dQ) override
Return in vector dq the diagonal of matrix fQ (Quadratic part of Objective function)
void DataRandom(TVectorD &x, TVectorD &y, TVectorD &z, TVectorD &s) override
Choose randomly a QP problem.
TMatrixD fC
Definition TQpDataDens.h:70
void Amult(Double_t beta, TVectorD &y, Double_t alpha, const TVectorD &x) override
calculate y = beta*y + alpha*(fA*x)
TMatrixDSym fQ
Definition TQpDataDens.h:68
Double_t ObjectiveValue(TQpVar *vars) override
Return value of the objective function.
TMatrixD fA
Definition TQpDataDens.h:69
void PutCIntoAt(TMatrixDBase &M, Int_t row, Int_t col) override
Insert the constraint matrix C into the matrix M at index (row,col) for the fundamental linear system...
void Print(Option_t *opt="") const override
Print all class members.
void PutQIntoAt(TMatrixDBase &M, Int_t row, Int_t col) override
Insert the Hessian Q into the matrix M at index (row,col) for the fundamental linear system.
void PutAIntoAt(TMatrixDBase &M, Int_t row, Int_t col) override
Insert the constraint matrix A into the matrix M at index (row,col) for the fundamental linear system...
void ATransmult(Double_t beta, TVectorD &y, Double_t alpha, const TVectorD &x) override
calculate y = beta*y + alpha*(fA^T*x)
TQpDataDens & operator=(const TQpDataDens &source)
Assignment operator.
void CTransmult(Double_t beta, TVectorD &y, Double_t alpha, const TVectorD &x) override
calculate y = beta*y + alpha*(fC^T*x)
Class containing the variables for the general QP formulation.
Definition TQpVar.h:60
TVectorD fX
Definition TQpVar.h:91
Element NormInf() const
Compute the infinity-norm of the vector MAX{ |v[i]| }.
Definition TVectorT.cxx:603
Bool_t MatchesNonZeroPattern(const TVectorT< Element > &select)
Check if vector elements as selected through array select are non-zero.
TVectorT< Element > & ResizeTo(Int_t lwb, Int_t upb)
Resize the vector to [lwb:upb] .
Definition TVectorT.cxx:294
Int_t GetNrows() const
Definition TVectorT.h:73
TVectorT< Element > & SelectNonZeros(const TVectorT< Element > &select)
Keep only element as selected through array select non-zero.
Definition TVectorT.cxx:544
void Print(Option_t *option="") const override
Print the vector as a list of elements.
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
const Int_t n
Definition legend1.C:16
Short_t Min(Short_t a, Short_t b)
Returns the smallest of a and b.
Definition TMathBase.h:198
TMatrixT< Element > & Add(TMatrixT< Element > &target, Element scalar, const TMatrixT< Element > &source)
Modify addition: target += scalar * source.
TMarker m
Definition textangle.C:8