ConicBundle
CMlowrankss.hxx
Go to the documentation of this file.
1 
2 #ifndef CONICBUNDLE_CMLOWRANKSS_HXX
3 #define CONICBUNDLE_CMLOWRANKSS_HXX
4 
12 #include "Coeffmat.hxx"
13 #include "CMlowrankdd.hxx"
14 
15 namespace ConicBundle {
16 
20 
24  class CMlowrankss: public Coeffmat
25  {
26  private:
29  public:
32  {A=Ain;B=Bin;CM_type=CM_lowrankss;infop=cip;}
34  virtual ~CMlowrankss(){}
35 
37  virtual Coeffmat* clone() const
38  {return new CMlowrankss(A,B,ConicBundle::clone(infop));}
39 
41  virtual CH_Matrix_Classes::Integer dim() const { return A.rowdim(); }
42 
45  { return CH_Matrix_Classes::ip(A.row(i),B.row(j))+ CH_Matrix_Classes::ip(B.row(i),A.row(j));}
46 
50 
52  virtual CH_Matrix_Classes::Real norm(void) const
55  return sqrt(2.*CH_Matrix_Classes::ip(C,D)+d);}
56 
58  virtual Coeffmat* subspace(const CH_Matrix_Classes::Matrix& P) const
60  return new CMlowrankdd(C,D,ConicBundle::clone(infop)); }
61 
64  { A*=d; if (infop) infop->multiply(d);}
65 
69 
73  return 2.*CH_Matrix_Classes::ip(C,D); }
74 
77  {
80  const CH_Matrix_Classes::Integer brd=C.rowdim();
81  const CH_Matrix_Classes::Integer pcd=P.coldim();
82  const CH_Matrix_Classes::Integer prd=P.rowdim();
83  const CH_Matrix_Classes::Real *pp=P.get_store()+start_row;
85  const CH_Matrix_Classes::Real *vp=(A.get_colval()).get_store();
86  const CH_Matrix_Classes::Integer *ip=(A.get_colindex()).get_store();
87  for(CH_Matrix_Classes::Integer j=0;j<colinfo.rowdim();j++){
88  CH_Matrix_Classes::Integer indj=colinfo(j,0);
89  for(CH_Matrix_Classes::Integer i=colinfo(j,1);--i>=0;){
90  CH_Matrix_Classes::mat_xpeya(pcd,bp+indj,brd,pp+(*ip++),prd,(*vp++));
91  }
92  }
94  bp=D.get_store();
95  const CH_Matrix_Classes::Indexmatrix &colinfo2=B.get_colinfo();
96  vp=(B.get_colval()).get_store();
97  ip=(B.get_colindex()).get_store();
98  {for(CH_Matrix_Classes::Integer j=0;j<colinfo2.rowdim();j++){
99  CH_Matrix_Classes::Integer indj=colinfo2(j,0);
100  for(CH_Matrix_Classes::Integer i=colinfo2(j,1);--i>=0;){
101  CH_Matrix_Classes::mat_xpeya(pcd,bp+indj,brd,pp+(*ip++),prd,(*vp++));
102  }
103  }}
104  CH_Matrix_Classes::Real trval=0;
105  if (Lam==0) {
106  trval=CH_Matrix_Classes::ip(C,D);
107  }
108  else {
109  assert(Lam->dim()==P.coldim());
110  const CH_Matrix_Classes::Real *cp=C.get_store();
111  const CH_Matrix_Classes::Real *dp=D.get_store();
112  const CH_Matrix_Classes::Real *lp=Lam->get_store();
113  for (CH_Matrix_Classes::Integer i=0;i<C.coldim();i++,cp+=C.rowdim(),dp+=D.rowdim())
114  trval+=(*lp++)*CH_Matrix_Classes::mat_ip(C.rowdim(),cp,dp);
115  }
116  return 2.*trval;
117  }
118 
122 
126  CH_Matrix_Classes::genmult(B,CH_Matrix_Classes::genmult(A,C,E,1.,0.,1),D,d,1.);}
127 
131  CH_Matrix_Classes::genmult(B,CH_Matrix_Classes::genmult(A,C,E,1.,0.,1),D,d,1.);}
132 
135  {
136  CH_Matrix_Classes::Matrix tmp1; CH_Matrix_Classes::genmult(P,A,tmp1,1.,0.,1,0);
137  CH_Matrix_Classes::Matrix tmp2; CH_Matrix_Classes::genmult(B,Q,tmp2,1.,0.,1,0);
138  CH_Matrix_Classes::genmult(tmp1,tmp2,R,1.,0.,0,0);
139  CH_Matrix_Classes::genmult(P,B,tmp1,1.,0.,1,0);
140  CH_Matrix_Classes::genmult(A,Q,tmp2,1.,0.,1,0);
141  CH_Matrix_Classes::genmult(tmp1,tmp2,R,1.,1.,0,0);
142  }
143 
146  { return 4*A.nonzeros()+4*B.nonzeros(); }
147 
149  virtual int dense() const
150  {return 0;}
151 
153  virtual int sparse() const
154  { return 0;}
155 
159  CH_Matrix_Classes::Matrix& /* val */,
160  CH_Matrix_Classes::Real /* d=1. */)const
161  {return 0;}
162 
164  virtual int support_in(const CH_Matrix_Classes::Sparsesym& /* S */) const
165  {return 0;}
166 
169  {return 2.*CH_Matrix_Classes::ip(A,S*B);}
170 
174  CH_Matrix_Classes::rank2add(C,D,S,2.);}
175  // S=P^t*(*this)*P
176 
179  {
182  const CH_Matrix_Classes::Integer brd=C.rowdim();
183  const CH_Matrix_Classes::Integer pcd=P.coldim();
184  const CH_Matrix_Classes::Integer prd=P.rowdim();
185  const CH_Matrix_Classes::Real *pp=P.get_store()+start_row;
186  const CH_Matrix_Classes::Indexmatrix &colinfo=A.get_colinfo();
187  const CH_Matrix_Classes::Real *vp=(A.get_colval()).get_store();
188  const CH_Matrix_Classes::Integer *ip=(A.get_colindex()).get_store();
189  for(CH_Matrix_Classes::Integer j=0;j<colinfo.rowdim();j++){
190  CH_Matrix_Classes::Integer indj=colinfo(j,0);
191  for(CH_Matrix_Classes::Integer i=colinfo(j,1);--i>=0;){
192  CH_Matrix_Classes::mat_xpeya(pcd,bp+indj,brd,pp+(*ip++),prd,(*vp++));
193  }
194  }
196  bp=D.get_store();
197  const CH_Matrix_Classes::Indexmatrix &colinfo2=B.get_colinfo();
198  vp=(B.get_colval()).get_store();
199  ip=(B.get_colindex()).get_store();
200  {for(CH_Matrix_Classes::Integer j=0;j<colinfo2.rowdim();j++){
201  CH_Matrix_Classes::Integer indj=colinfo2(j,0);
202  for(CH_Matrix_Classes::Integer i=colinfo2(j,1);--i>=0;){
203  CH_Matrix_Classes::mat_xpeya(pcd,bp+indj,brd,pp+(*ip++),prd,(*vp++));
204  }
205  }}
206  CH_Matrix_Classes::rank2add(C,D,S,2.*alpha,1.,1);
207  }
208 
211  CH_Matrix_Classes::Real alpha=1.,CH_Matrix_Classes::Real beta=0.,int dtrans=0) const
212  {
214  CH_Matrix_Classes::genmult(A,CH_Matrix_Classes::genmult(B,D,E,1.,0.,1,dtrans),C,alpha,beta);
215  return CH_Matrix_Classes::genmult(B,CH_Matrix_Classes::genmult(A,D,E,1.,0.,1,dtrans),C,alpha,1.);
216  }
217 
220  CH_Matrix_Classes::Real alpha=1.,CH_Matrix_Classes::Real beta=0.,int dtrans=0) const
221  {
223  CH_Matrix_Classes::genmult(CH_Matrix_Classes::genmult(D,A,E,1.,0.,dtrans),B,C,alpha,beta,0,1);
224  return CH_Matrix_Classes::genmult(CH_Matrix_Classes::genmult(D,B,E,1.,0.,dtrans),A,C,alpha,1.,0,1);
225  }
226 
228  virtual int equal(const Coeffmat* p,double tol=1e-6) const
229  {
230  const CMlowrankss *pp=dynamic_cast<const CMlowrankss *>(p);
231  if (pp==0)
232  return 0;
233  if (CH_Matrix_Classes::equal(A,pp->A,tol) && CH_Matrix_Classes::equal(B,pp->B,tol))
234  return 1;
235  if (CH_Matrix_Classes::equal(A,pp->B,tol) && CH_Matrix_Classes::equal(B,pp->A,tol))
236  return 1;
237  return 0;
238  }
239 
241  virtual std::ostream& display(std::ostream& o) const
242  {o<<"CMlowrankss\n";A.display(o);B.display(o);return o;}
243 
245  virtual std::ostream& out(std::ostream& o) const
246  {return o<<"LOWRANK_SPARSE_SPARSE\n"<<A<<B;}
247 
249  virtual std::istream& in(std::istream& i)
250  {return i>>A>>B;
251  if((A.rowdim()!=B.rowdim())||(A.coldim()!=B.coldim())){
252  i.clear(std::ios::failbit);
253  if (CH_Matrix_Classes::materrout) (*CH_Matrix_Classes::materrout)<<"*** ERROR: CMlowrankss::in(): dimensions of A and B do not match"<<std::endl;
254  }
255  return i;
256  }
257 
259  CMlowrankss(std::istream& is,CoeffmatInfo* cip=0)
260  {CM_type=CM_lowrankss;infop=cip;in(is);}
261 
262  };
263 
265 
266 }
267 
268 #endif
269 
int Integer
all integer numbers in calculations and indexing are of this type
Definition: matop.hxx:40
virtual void multiply(CH_Matrix_Classes::Real d)
multiply constraint permanentely by d; this is to allow scaling or sign changes in the constraints ...
Definition: CMlowrankss.hxx:63
CH_Matrix_Classes::Sparsemat A
this is A in A*B^T+B*A^T
Definition: CMlowrankss.hxx:27
Header declaring the classes ConicBundle::Coeffmat, ConicBundle::CoeffmatPointer, ConicBundle::Coeffm...
virtual void make_symmatrix(CH_Matrix_Classes::Symmatrix &S) const
returns a dense symmetric constraint matrix
Definition: CMlowrankss.hxx:48
virtual int sparse(CH_Matrix_Classes::Indexmatrix &, CH_Matrix_Classes::Indexmatrix &, CH_Matrix_Classes::Matrix &, CH_Matrix_Classes::Real) const
returns 0 if not sparse. If it is sparse it returns 1 and the nonzero structure in I...
Definition: CMlowrankss.hxx:157
virtual CH_Matrix_Classes::Real norm(void) const
returns the Frobenius norm of the matrix
Definition: CMlowrankss.hxx:52
Integer rowdim() const
returns the row dimension
Definition: matrix.hxx:215
std::ostream * materrout
if not zero, this is the output stream for runtime error messages, by default it is set to &std::cout...
double Real
all real numbers in calculations are of this type
Definition: matop.hxx:50
allows to memorize the scalings applied to a Coeffmat and offers the basis for storing further user d...
Definition: Coeffmat.hxx:52
virtual CH_Matrix_Classes::Real ip(const CH_Matrix_Classes::Symmatrix &S) const
returns ip(*this,S)=trace(*this*S), the trace inner product
Definition: CMlowrankss.hxx:67
virtual Coeffmat * subspace(const CH_Matrix_Classes::Matrix &P) const
delivers a new object on the heap corresponding to the matrix P^T(*this)P, the caller is responsible ...
Definition: CMlowrankss.hxx:58
const Indexmatrix & get_colinfo() const
returns information on nonzero columns, k by 3, listing: index, %nonzeros, first index in colindex/co...
Definition: sparsmat.hxx:265
defines a base class for coefficient matrices in semidefinite programming, in particular for use with...
Definition: Coeffmat.hxx:125
void display(std::ostream &out, int precision=0, int width=0, int screenwidth=0) const
displays a matrix in a pretty way for bounded screen widths; for variables of value zero default valu...
Matrix class for integral values of type Integer
Definition: indexmat.hxx:195
CMlowrankss(std::istream &is, CoeffmatInfo *cip=0)
constructor with istream and possibly additional user information
Definition: CMlowrankss.hxx:259
virtual int sparse() const
returns 0 if not sparse, otherwise 1
Definition: CMlowrankss.hxx:153
virtual void add_projection(CH_Matrix_Classes::Symmatrix &S, const CH_Matrix_Classes::Matrix &P, CH_Matrix_Classes::Real alpha=1., CH_Matrix_Classes::Integer start_row=0) const
computes S+=Q^T(*this)Q for Q=P.rows(start_row,start_row+dim-1)
Definition: CMlowrankss.hxx:178
virtual int dense() const
returns 1 if its structure is as bad as its dense symmetric representation, otherwise 0 ...
Definition: CMlowrankss.hxx:149
virtual CH_Matrix_Classes::Integer prodvec_flops() const
returns an estimate of number of flops to compute addprodto for a vector
Definition: CMlowrankss.hxx:145
const Matrix & get_colval() const
returns the value vector of the column representation holding the value for each element ...
Definition: sparsmat.hxx:269
CMlowrankss(const CH_Matrix_Classes::Sparsemat &Ain, const CH_Matrix_Classes::Sparsemat &Bin, CoeffmatInfo *cip=0)
copy Ain, Bin and store the user information
Definition: CMlowrankss.hxx:31
Real * get_store()
returns the current address of the internal value array; use cautiously, do not use delete! ...
Definition: matrix.hxx:326
Matrix class of symmetric matrices with real values of type Real
Definition: symmat.hxx:43
virtual CH_Matrix_Classes::Real gramip(const CH_Matrix_Classes::Matrix &P) const
returns ip(*this,PP^T)=trace P^T(*this)P
Definition: CMlowrankss.hxx:71
virtual CH_Matrix_Classes::Real gramip(const CH_Matrix_Classes::Matrix &P, CH_Matrix_Classes::Integer start_row, const CH_Matrix_Classes::Matrix *Lam=0) const
returns ip(*this,QQ^T)=trace Q^T(*this)Q for Q=P.rows(start_row,start_row+dim-1)
Definition: CMlowrankss.hxx:76
virtual void left_right_prod(const CH_Matrix_Classes::Matrix &P, const CH_Matrix_Classes::Matrix &Q, CH_Matrix_Classes::Matrix &R) const
computes R=P^T*(*this)*Q
Definition: CMlowrankss.hxx:134
virtual void addprodto(CH_Matrix_Classes::Matrix &D, const CH_Matrix_Classes::Matrix &C, CH_Matrix_Classes::Real d=1.) const
computes D+=d*(*this)*C
Definition: CMlowrankss.hxx:124
conic bundle method solver for sum of convex functions. See the ConicBundle_Manual for a quick introd...
Definition: CBSolver.hxx:22
virtual CH_Matrix_Classes::Real ip(const CH_Matrix_Classes::Sparsesym &S) const
returns the inner product of the constraint matrix with S
Definition: CMlowrankss.hxx:168
Coeffmattype CM_type
in order to enable type identification
Definition: Coeffmat.hxx:134
virtual int equal(const Coeffmat *p, double tol=1e-6) const
returns 1, if p is the same derived class and entries differ by less than tol, otherwise zero ...
Definition: CMlowrankss.hxx:228
CoeffmatInfo * infop
allows the user to specify and output additional information
Definition: Coeffmat.hxx:135
virtual std::ostream & out(std::ostream &o) const
put entire contents onto ostream with the class type in the beginning so that the derived class can b...
Definition: CMlowrankss.hxx:245
virtual CH_Matrix_Classes::Integer dim() const
returns the order of the represented symmetric matrix
Definition: CMlowrankss.hxx:41
virtual void addprodto(CH_Matrix_Classes::Matrix &D, const CH_Matrix_Classes::Sparsemat &C, CH_Matrix_Classes::Real d=1.) const
computes D+=d*(*this)*C
Definition: CMlowrankss.hxx:129
Integer trace(const Indexmatrix &A)
returns the sum of the diagonal elements A(i,i) over all i
Matrix class of symmetric matrices with real values of type Real
Definition: sparssym.hxx:89
virtual CH_Matrix_Classes::Real operator()(CH_Matrix_Classes::Integer i, CH_Matrix_Classes::Integer j) const
returns the value of the matrix element (i,j)
Definition: CMlowrankss.hxx:44
void multiply(CH_Matrix_Classes::Real sf)
scales the scale factor
Definition: Coeffmat.hxx:70
CH_Matrix_Classes::Sparsemat B
this is B in A*B^T+B*A^T
Definition: CMlowrankss.hxx:28
implements a low rank matrix as Coeffmat with each a dense rectangular CH_Matrix_Classes::Matrix (f...
Definition: CMlowrankdd.hxx:24
Integer coldim() const
returns the column dimension
Definition: matrix.hxx:218
Matrix class for real values of type Real
Definition: matrix.hxx:74
virtual const CH_Matrix_Classes::Matrix & postgenmult(const CH_Matrix_Classes::Matrix &D, CH_Matrix_Classes::Matrix &C, CH_Matrix_Classes::Real alpha=1., CH_Matrix_Classes::Real beta=0., int dtrans=0) const
computes C= alpha*(*this)*D^(T if dtrans) + beta*C, C is also returned
Definition: CMlowrankss.hxx:210
virtual int support_in(const CH_Matrix_Classes::Sparsesym &) const
returns 0 if the support of the costraint matrix is not contained in the support of the sparse symmet...
Definition: CMlowrankss.hxx:164
virtual void project(CH_Matrix_Classes::Symmatrix &S, const CH_Matrix_Classes::Matrix &P) const
computes S=P^T*(*this)*P
Definition: CMlowrankss.hxx:172
Val mat_ip(Integer len, const Val *x, const Val *y, const Val *d=0)
return sum(x[i]*y[i]) summing over len elements of the arrays x and y.
Definition: matop.hxx:1096
Matrix class of sparse matrices with real values of type Real
Definition: sparsmat.hxx:74
for CMlowrankss
Definition: Coeffmat.hxx:35
void mat_xpeya(Integer len, Val *x, const Val *y, const Val a)
Set x[i]+=a*y[i] for len elements of the arrays x and y.
Definition: matop.hxx:543
const Indexmatrix & get_colindex() const
returns the index vector of the column representation holding the row index for each element ...
Definition: sparsmat.hxx:267
Integer rowdim() const
returns the row dimension
Definition: indexmat.hxx:321
virtual void addmeto(CH_Matrix_Classes::Symmatrix &S, CH_Matrix_Classes::Real d=1.) const
computes S+=d*(*this);
Definition: CMlowrankss.hxx:120
implements a low rank matrix as Coeffmat with each a sparse rectangular CH_Matrix_Classes::Sparsema...
Definition: CMlowrankss.hxx:24
Indexmatrix & genmult(const Indexmatrix &A, const Indexmatrix &B, Indexmatrix &C, Integer alpha=1, Integer beta=0, int atrans=0, int btrans=0)
returns C=beta*C+alpha*A*B, where A and B may be transposed; C must not be equal to A and B; if beta=...
Sparsemat row(Integer i) const
returns row i copied to a new sparse matrix
virtual std::istream & in(std::istream &i)
counterpart to out(), does not read the class type, though. This is assumed to have been read in orde...
Definition: CMlowrankss.hxx:249
Symmatrix & rank2add(const Matrix &A, const Matrix &B, Symmatrix &C, Real alpha=1., Real beta=0., int trans=0)
returns C=beta*C+alpha*(A*B^T+B*A^T)/2 [or for transposed (A^T*B+B^T*A)/2]. If beta==0. then C is initiliazed to the correct size.
bool equal(const Matrix &A, const Matrix &B)
returns true if both matrices have the same size and the same elements
Definition: matrix.hxx:1433
Integer coldim() const
returns the column dimension
Definition: sparsmat.hxx:202
virtual Coeffmat * clone() const
makes an explicit copy of itself and returns a pointer to it
Definition: CMlowrankss.hxx:37
Integer rowdim() const
returns the row dimension
Definition: sparsmat.hxx:199
double sqrt(int a)
return sqrt for int a
Definition: mymath.hxx:121
CoeffmatInfo * clone(const CoeffmatInfo *cip)
if cip is not zero, it calls and returns cip->clone() and 0 otherwise
Definition: Coeffmat.hxx:106
virtual const CH_Matrix_Classes::Matrix & pregenmult(const CH_Matrix_Classes::Matrix &D, CH_Matrix_Classes::Matrix &C, CH_Matrix_Classes::Real alpha=1., CH_Matrix_Classes::Real beta=0., int dtrans=0) const
computes C= alpha*D^(T if dtrans)*(*this) + beta*C, C is also returned
Definition: CMlowrankss.hxx:219
virtual std::ostream & display(std::ostream &o) const
display constraint information
Definition: CMlowrankss.hxx:241
Real ip(const Matrix &A, const Matrix &B)
returns the usual inner product of A and B, i.e., the sum of A(i,j)*B(i,j) over all i...
Definition: matrix.hxx:1165
Header declaring the class ConicBundle::CMlowrankdd (needed for ConicBundle::AffineMatrixFunction) ...
Integer nonzeros() const
returns the number of nonzeros
Definition: sparsmat.hxx:205