ROOT_528-00b_version: math/quadp/src/TQpProbSparse.cxx Source File

00001 // @(#)root/quadp:$Id: TQpProbSparse.cxx 20882 2007-11-19 11:31:26Z rdm $
00002 // Author: Eddy Offermann   May 2004
00003 
00004 /*************************************************************************
00005  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers.               *
00006  * All rights reserved.                                                  *
00007  *                                                                       *
00008  * For the licensing terms see $ROOTSYS/LICENSE.                         *
00009  * For the list of contributors see $ROOTSYS/README/CREDITS.             *
00010  *************************************************************************/
00011 
00012 /*************************************************************************
00013  * Parts of this file are copied from the OOQP distribution and          *
00014  * are subject to the following license:                                 *
00015  *                                                                       *
00016  * COPYRIGHT 2001 UNIVERSITY OF CHICAGO                                  *
00017  *                                                                       *
00018  * The copyright holder hereby grants you royalty-free rights to use,    *
00019  * reproduce, prepare derivative works, and to redistribute this software*
00020  * to others, provided that any changes are clearly documented. This     *
00021  * software was authored by:                                             *
00022  *                                                                       *
00023  *   E. MICHAEL GERTZ      gertz@mcs.anl.gov                             *
00024  *   Mathematics and Computer Science Division                           *
00025  *   Argonne National Laboratory                                         *
00026  *   9700 S. Cass Avenue                                                 *
00027  *   Argonne, IL 60439-4844                                              *
00028  *                                                                       *
00029  *   STEPHEN J. WRIGHT     swright@cs.wisc.edu                           *
00030  *   Computer Sciences Department                                        *
00031  *   University of Wisconsin                                             *
00032  *   1210 West Dayton Street                                             *
00033  *   Madison, WI 53706   FAX: (608)262-9777                              *
00034  *                                                                       *
00035  * Any questions or comments may be directed to one of the authors.      *
00036  *                                                                       *
00037  * ARGONNE NATIONAL LABORATORY (ANL), WITH FACILITIES IN THE STATES OF   *
00038  * ILLINOIS AND IDAHO, IS OWNED BY THE UNITED STATES GOVERNMENT, AND     *
00039  * OPERATED BY THE UNIVERSITY OF CHICAGO UNDER PROVISION OF A CONTRACT   *
00040  * WITH THE DEPARTMENT OF ENERGY.                                        *
00041  *************************************************************************/
00042 
00043 #include "TQpProbSparse.h"
00044 #include "TMatrixD.h"
00045 
00046 //////////////////////////////////////////////////////////////////////////
00047 //                                                                      //
00048 // TQpProbSparse                                                        //
00049 //                                                                      //
00050 // dense matrix problem formulation                                     //
00051 //                                                                      //
00052 //////////////////////////////////////////////////////////////////////////
00053 
00054 ClassImp(TQpProbSparse)
00055 
00056 //______________________________________________________________________________
00057 TQpProbSparse::TQpProbSparse(Int_t nx,Int_t my,Int_t mz) :
00058                TQpProbBase(nx,my,mz)
00059 {
00060 // Constructor
00061 
00062    // We do not want more constrains than variables
00063    R__ASSERT(nx-my-mz > 0);
00064 }
00065 
00066 
00067 //______________________________________________________________________________
00068 TQpProbSparse::TQpProbSparse(const TQpProbSparse &another) : TQpProbBase(another)
00069 {
00070 // Copy constructor
00071 
00072    *this = another;
00073 }
00074 
00075 
00076 //______________________________________________________________________________
00077 TQpDataBase *TQpProbSparse::MakeData(Double_t *c,
00078                                      Int_t nnzQ,Int_t *irowQ,Int_t *icolQ,Double_t *Q,
00079                                      Double_t *xlo,Bool_t *ixlo,
00080                                      Double_t *xup,Bool_t *ixup,
00081                                      Int_t nnzA,Int_t *irowA,Int_t *icolA,Double_t *A,
00082                                      Double_t *bA,
00083                                      Int_t nnzC,Int_t *irowC,Int_t *icolC,Double_t *C,
00084                                      Double_t *clo,Bool_t *iclo,
00085                                      Double_t *cup,Bool_t *icup)
00086 {
00087 // Setup the data
00088 
00089    TVectorD       vc  ; vc  .Use(fNx,c);
00090    TMatrixDSparse mQ  ; mQ  .Use(fNx,fNx,nnzQ,irowQ,icolQ,Q);
00091    TVectorD       vxlo; vxlo.Use(fNx,xlo);
00092    TVectorD       vxup; vxup.Use(fNx,xup);
00093    TMatrixDSparse mA  ;
00094    TVectorD       vbA ;
00095    if (fMy > 0) {
00096       mA  .Use(fMy,fNx,nnzA,irowA,icolA,A);
00097       vbA .Use(fMy,bA);
00098    }
00099    TMatrixDSparse mC  ;
00100    TVectorD       vclo;
00101    TVectorD       vcup;
00102    if (fMz > 0) {
00103       mC  .Use(fMz,fNx,nnzC,irowC,icolC,C);
00104       vclo.Use(fMz,clo);
00105       vcup.Use(fMz,cup);
00106    }
00107 
00108    TVectorD vixlo(fNx);
00109    TVectorD vixup(fNx);
00110    for (Int_t ix = 0; ix < fNx; ix++) {
00111       vixlo[ix] = (ixlo[ix]) ? 1.0 : 0.0;
00112       vixup[ix] = (ixup[ix]) ? 1.0 : 0.0;
00113    }
00114 
00115    TVectorD viclo(fMz);
00116    TVectorD vicup(fMz);
00117    for (Int_t ic = 0; ic < fMz; ic++) {
00118       viclo[ic] = (iclo[ic]) ? 1.0 : 0.0;
00119       vicup[ic] = (icup[ic]) ? 1.0 : 0.0;
00120    }
00121 
00122    TQpDataSparse *data = new TQpDataSparse(vc,mQ,vxlo,vixlo,vxup,vixup,mA,vbA,mC,vclo,
00123       viclo,vcup,vicup);
00124 
00125    return data;
00126 }
00127 
00128 
00129 //______________________________________________________________________________
00130 TQpDataBase *TQpProbSparse::MakeData(TVectorD     &c,
00131                                      TMatrixDBase &Q_in,
00132                                      TVectorD     &xlo, TVectorD &ixlo,
00133                                      TVectorD     &xup, TVectorD &ixup,
00134                                      TMatrixDBase &A_in,TVectorD &bA,
00135                                      TMatrixDBase &C_in,
00136                                      TVectorD     &clo, TVectorD &iclo,
00137                                      TVectorD     &cup, TVectorD &icup)
00138 {
00139 // Setup the data
00140 
00141    TMatrixDSparse &mQ = (TMatrixDSparse &) Q_in;
00142    TMatrixDSparse &mA = (TMatrixDSparse &) A_in;
00143    TMatrixDSparse &mC = (TMatrixDSparse &) C_in;
00144 
00145    R__ASSERT(mQ.GetNrows() == fNx && mQ.GetNcols() == fNx);
00146    if (fMy > 0) R__ASSERT(mA.GetNrows() == fMy && mA.GetNcols() == fNx);
00147    else         R__ASSERT(mA.GetNrows() == fMy);
00148    if (fMz > 0) R__ASSERT(mC.GetNrows() == fMz && mC.GetNcols() == fNx);
00149    else         R__ASSERT(mC.GetNrows() == fMz);
00150 
00151    R__ASSERT(c.GetNrows()    == fNx);
00152    R__ASSERT(xlo.GetNrows()  == fNx);
00153    R__ASSERT(ixlo.GetNrows() == fNx);
00154    R__ASSERT(xup.GetNrows()  == fNx);
00155    R__ASSERT(ixup.GetNrows() == fNx);
00156 
00157    R__ASSERT(bA.GetNrows()   == fMy);
00158    R__ASSERT(clo.GetNrows()  == fMz);
00159    R__ASSERT(iclo.GetNrows() == fMz);
00160    R__ASSERT(cup.GetNrows()  == fMz);
00161    R__ASSERT(icup.GetNrows() == fMz);
00162 
00163    TQpDataSparse *data = new TQpDataSparse(c,mQ,xlo,ixlo,xup,ixup,mA,bA,mC,clo,iclo,cup,icup);
00164 
00165    return data;
00166 }
00167 
00168 
00169 //______________________________________________________________________________
00170 TQpResidual* TQpProbSparse::MakeResiduals(const TQpDataBase *data_in)
00171 {
00172 // Setup the residuals
00173 
00174    TQpDataSparse *data = (TQpDataSparse *) data_in;
00175    return new TQpResidual(fNx,fMy,fMz,data->fXloIndex,data->fXupIndex,data->fCloIndex,data->fCupIndex);
00176 }
00177 
00178 
00179 //______________________________________________________________________________
00180 TQpVar* TQpProbSparse::MakeVariables(const TQpDataBase *data_in)
00181 {
00182 // Setup the variables
00183 
00184    TQpDataSparse *data = (TQpDataSparse *) data_in;
00185 
00186    return new TQpVar(fNx,fMy,fMz,data->fXloIndex,data->fXupIndex,data->fCloIndex,data->fCupIndex);
00187 }
00188 
00189 
00190 //______________________________________________________________________________
00191 TQpLinSolverBase* TQpProbSparse::MakeLinSys(const TQpDataBase *data_in)
00192 {
00193 // Setup the linear solver
00194 
00195    TQpDataSparse *data = (TQpDataSparse *) data_in;
00196    return new TQpLinSolverSparse(this,data);
00197 }
00198 
00199 
00200 //______________________________________________________________________________
00201 void TQpProbSparse::JoinRHS(TVectorD &rhs,TVectorD &rhs1_in,TVectorD &rhs2_in,TVectorD &rhs3_in)
00202 {
00203 // Assembles a single vector object from three given vectors .
00204 //     rhs_out (output) final joined vector
00205 //     rhs1_in (input) first part of rhs
00206 //     rhs2_in (input) middle part of rhs
00207 //     rhs3_in (input) last part of rhs .
00208 
00209    rhs.SetSub(0,rhs1_in);
00210    if (fMy > 0) rhs.SetSub(fNx,    rhs2_in);
00211    if (fMz > 0) rhs.SetSub(fNx+fMy,rhs3_in);
00212 }
00213 
00214 
00215 //______________________________________________________________________________
00216 void TQpProbSparse::SeparateVars(TVectorD &x_in,TVectorD &y_in,TVectorD &z_in,TVectorD &vars_in)
00217 {
00218 // Extracts three component vectors from a given aggregated vector.
00219 //     vars_in  (input) aggregated vector
00220 //     x_in (output) first part of vars
00221 //     y_in (output) middle part of vars
00222 //     z_in (output) last part of vars
00223 
00224    x_in = vars_in.GetSub(0,fNx-1);
00225    if (fMy > 0) y_in = vars_in.GetSub(fNx,    fNx+fMy-1);
00226    if (fMz > 0) z_in = vars_in.GetSub(fNx+fMy,fNx+fMy+fMz-1);
00227 }
00228 
00229 
00230 //______________________________________________________________________________
00231 void TQpProbSparse::MakeRandomData(TQpDataSparse *&data,TQpVar *&soln,Int_t nnzQ,Int_t nnzA,Int_t nnzC)
00232 {
00233 // Create a random QP problem
00234 
00235    data = new TQpDataSparse(fNx,fMy,fMz);
00236    soln = this->MakeVariables(data);
00237    data->SetNonZeros(nnzQ,nnzA,nnzC);
00238    data->DataRandom(soln->fX,soln->fY,soln->fZ,soln->fS);
00239 }
00240 
00241 
00242 //______________________________________________________________________________
00243 TQpProbSparse &TQpProbSparse::operator=(const TQpProbSparse &source)
00244 {
00245 // Assignment operator
00246 
00247    if (this != &source) {
00248       TQpProbBase::operator=(source);
00249    }
00250    return *this;
00251 }