Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TQpLinSolverDens.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////////////////////////////////////////////////////////////////////////////////
44///
45/// \class TQpLinSolverDens
46///
47/// Implements the aspects of the solvers for dense general QP
48/// formulation that are specific to the dense case.
49///
50////////////////////////////////////////////////////////////////////////////////
51
52#include "TQpLinSolverDens.h"
53
55
56////////////////////////////////////////////////////////////////////////////////
57/// Constructor
58
60 TQpLinSolverBase(factory,data)
61{
62 const Int_t n = factory->fNx+factory->fMy+factory->fMz;
64
65 data->PutQIntoAt(fKkt,0, 0);
66 if (fMy > 0) data->PutAIntoAt(fKkt,fNx, 0);
67 if (fMz > 0) data->PutCIntoAt(fKkt,fNx+fMy,0);
68 for (Int_t ix = fNx; ix < fNx+fMy+fMz; ix++) {
69 for (Int_t iy = fNx; iy < fNx+fMy+fMz; iy++)
70 fKkt(ix,iy) = 0.0;
71 }
72
74}
75
76
77////////////////////////////////////////////////////////////////////////////////
78/// Copy constructor
79
81{
82 *this = another;
83}
84
85
86////////////////////////////////////////////////////////////////////////////////
87/// Sets up the matrix for the main linear system in "augmented system" form.
88
90{
91 TQpLinSolverBase::Factor(prob,vars);
93}
94
95
96////////////////////////////////////////////////////////////////////////////////
97/// Places the diagonal resulting from the bounds on x into the augmented system matrix
98
100{
101 TMatrixDDiag diag(fKkt);
102 for (Int_t i = 0; i < xdiag.GetNrows(); i++)
103 diag[i] = xdiag[i];
104}
105
106
107////////////////////////////////////////////////////////////////////////////////
108/// Places the diagonal resulting from the bounds on Cx into the augmented system matrix
109
111{
112 TMatrixDDiag diag(fKkt);
113 for (Int_t i = 0; i < zdiag.GetNrows(); i++)
114 diag[i+fNx+fMy] = zdiag[i];
115}
116
117
118////////////////////////////////////////////////////////////////////////////////
119/// Perform the actual solve using the factors produced in factor.
120/// rhs on input contains the aggregated right-hand side of the augmented system;
121/// on output contains the solution in aggregated form .
122
124{
125 fSolveLU.Solve(compressedRhs);
126}
127
128
129////////////////////////////////////////////////////////////////////////////////
130/// Assignment operator
131
133{
134 if (this != &source) {
136 fKkt.ResizeTo(source.fKkt); fKkt = source.fKkt;
137 fSolveLU = source.fSolveLU;
138 }
139 return *this;
140}
#define ClassImp(name)
Definition Rtypes.h:377
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
LU Decomposition class.
Definition TDecompLU.h:24
virtual void SetMatrix(const TMatrixD &a)
Set matrix to be decomposed.
Bool_t Solve(TVectorD &b) override
Solve Ax=b assuming the LU form of A is stored in fLU, but assume b has not been transformed.
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
Data for the dense QP formulation.
Definition TQpDataDens.h:63
Implementation of main solver for linear systems.
TQpLinSolverBase & operator=(const TQpLinSolverBase &source)
Assignment operator.
virtual void Factor(TQpDataBase *prob, TQpVar *vars)
Sets up the matrix for the main linear system in "augmented system" form.
Implements the aspects of the solvers for dense general QP formulation that are specific to the dense...
void PutZDiagonal(TVectorD &zdiag) override
Places the diagonal resulting from the bounds on Cx into the augmented system matrix.
void Factor(TQpDataBase *prob, TQpVar *vars) override
Sets up the matrix for the main linear system in "augmented system" form.
void SolveCompressed(TVectorD &rhs) override
Perform the actual solve using the factors produced in factor.
void PutXDiagonal(TVectorD &xdiag) override
Places the diagonal resulting from the bounds on x into the augmented system matrix.
TQpLinSolverDens & operator=(const TQpLinSolverDens &source)
Assignment operator.
dense matrix problem formulation
Definition TQpProbDens.h:61
Class containing the variables for the general QP formulation.
Definition TQpVar.h:60
Int_t GetNrows() const
Definition TVectorT.h:75
const Int_t n
Definition legend1.C:16