Kernel Quantum Probability Library
The KQP library aims at providing tools for working with quantums probabilities
feature_matrix.hpp
1 /*
2  This file is part of the Kernel Quantum Probability library (KQP).
3 
4  KQP is free software: you can redistribute it and/or modify
5  it under the terms of the GNU General Public License as published by
6  the Free Software Foundation, either version 3 of the License, or
7  (at your option) any later version.
8 
9  KQP is distributed in the hope that it will be useful,
10  but WITHOUT ANY WARRANTY; without even the implied warranty of
11  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12  GNU General Public License for more details.
13 
14  You should have received a copy of the GNU General Public License
15  along with KQP. If not, see <http://www.gnu.org/licenses/>.
16  */
17 
18 #ifndef __KQP_FEATURE_MATRIX_H__
19 #define __KQP_FEATURE_MATRIX_H__
20 
21 #include <kqp/kqp.hpp>
22 #include <kqp/logging.hpp>
23 #include <kqp/picojson.hpp>
24 #include <boost/shared_ptr.hpp>
25 #include <kqp/alt_matrix.hpp>
26 
27 
28 # include <kqp/define_header_logger.hpp>
29  DEFINE_KQP_HLOGGER("kqp.feature_matrix");
30 
31 namespace kqp
32 {
33 
34 #ifndef SWIG
35 template <typename Scalar>
37 {
38  typedef typename Eigen::NumTraits<Scalar>::Real Real;
40  typedef kqp::AltMatrix< typename kqp::AltVector<Real>::VectorType, typename kqp::AltVector<Real>::ConstantVectorType > RealAltVector;
41 
42  static inline bool isDense(const ScalarAltMatrix &m) {
43  return m.isT1();
44  }
45  static inline bool isIdentity(const ScalarAltMatrix &m) {
46  return m.isT2();
47  }
48 
49  static inline bool isDense(const RealAltVector &m) {
50  return m.isT1();
51  }
52  static inline bool isConstant(const RealAltVector &m) {
53  return m.isT2();
54  }
55 };
56 #endif
57 
58 #define KQP_SCALAR_TYPEDEFS(Scalar) \
59  typedef typename Eigen::NumTraits<Scalar>::Real Real; \
60  /* Scalar Matrices */ \
61  typedef Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> ScalarMatrix; \
62  typedef Eigen::Matrix<Scalar, Eigen::Dynamic, 1> ScalarVector; \
63  typedef Eigen::Matrix<Real, Eigen::Dynamic, Dynamic> RealMatrix; \
64  typedef Eigen::Matrix<Real, Eigen::Dynamic, 1> RealVector; \
65  typedef typename ScalarDefinitions<Scalar>::ScalarAltMatrix ScalarAltMatrix; \
66  typedef typename ScalarDefinitions<Scalar>::RealAltVector RealAltVector; \
67  typedef typename AltVector<Real>::ConstantVectorType ConstantRealVector; \
68  typedef typename AltDense<Scalar>::IdentityType IdentityScalarMatrix; \
69  /* Feature matrices */ \
70  typedef FeatureMatrixBase<Scalar> FMatrixBase; \
71  typedef boost::shared_ptr< FeatureMatrixBase<Scalar> > FMatrixPtr; \
72  typedef boost::shared_ptr< const FeatureMatrixBase<Scalar> > FMatrixCPtr; \
73  /* Feature spaces */ \
74  typedef SpaceBase<Scalar> FSpaceBase; \
75  typedef boost::shared_ptr< SpaceBase<Scalar> > FSpacePtr; \
76  typedef boost::shared_ptr< const SpaceBase<Scalar> > FSpaceCPtr; \
77  typedef boost::shared_ptr< FeatureMatrixBase<Scalar> > FMatrixBasePtr; \
78  typedef boost::shared_ptr< const FeatureMatrixBase<Scalar> > FMatrixBaseCPtr;
79 
81 #define KQP_MATRIX_TYPEDEFS(Scalar) \
82  KQP_SCALAR_TYPEDEFS(Scalar) \
83  typedef boost::shared_ptr<Self> SelfPtr; \
84  typedef boost::shared_ptr<const Self> SelfCPtr;
85 
87 #define KQP_SPACE_TYPEDEFS(SpaceName, Scalar) \
88  KQP_SCALAR_TYPEDEFS(Scalar) \
89  typedef boost::shared_ptr<Self> SelfPtr; \
90  typedef boost::shared_ptr<const Self> SelfCPtr; \
91  static const std::string &NAME() { static std::string _NAME(SpaceName); return _NAME; } \
92  const std::string &name() const { return NAME(); }
93 
94 
96 template<typename Scalar> struct ScalarInfo;
97 
98 template<> struct ScalarInfo<double>
99 {
100  static std::string name()
101  {
102  return "double";
103  }
104 };
105 template<> struct ScalarInfo<float>
106 {
107  static std::string name()
108  {
109  return "float";
110  }
111 };
112 template<> struct ScalarInfo<std::complex<double> >
113 {
114  static std::string name()
115  {
116  return "complex/double";
117  }
118 };
119 template<> struct ScalarInfo<std::complex<float> >
120 {
121  static std::string name()
122  {
123  return "complex/float";
124  }
125 };
126 
127 
128 template<typename Scalar> class FeatureMatrixBase;
129 template<typename Scalar> class SpaceBase;
130 
131 
135 template<typename Scalar> class KernelValues {
136 public:
137  KernelValues(Scalar inner, Scalar innerX, Scalar innerY) :
138  _inner(inner), _innerX(innerX), _innerY(innerY) {}
139 
140  KernelValues() {}
141 
142  inline Scalar inner(int mode = 0) const { if (mode == -1) return _innerX; return mode == 1 ? _innerY : _inner; }
143  inline Scalar innerX(int mode = 0) const { if (mode == 1) return _innerY; return _innerX; }
144  inline Scalar innerY(int mode = 0) const { if (mode == -1) return _innerX; return _innerY; }
145 
146  Scalar _inner, _innerX, _innerY;
147 };
148 
149 
156 template <class Scalar>
158 {
159 public:
160  KQP_SCALAR_TYPEDEFS(Scalar);
161 
162  virtual ~FeatureMatrixBase() {}
163 
165  virtual Index size() const = 0;
166 
168  virtual void add(const FMatrixBase &f, const std::vector<bool> *which = NULL) = 0;
169 
177  virtual FMatrixPtr subset(const std::vector<bool>::const_iterator &begin, const std::vector<bool>::const_iterator &end) const = 0;
178 
183  inline FMatrixPtr subset(const std::vector<bool> &list) const
184  {
185  return subset(list.begin(), list.end());
186  }
187 
189  virtual FMatrixBase &operator=(const FeatureMatrixBase<Scalar> &other) = 0;
190 
192  virtual FMatrixBasePtr copy() const = 0;
193 
194 
195 #ifndef SWIG
196 
197  template<typename T> inline const T &as() const
198  {
199  return kqp::our_dynamic_cast<const T &>(*this);
200  }
201  template<typename T> inline T &as()
202  {
203  return kqp::our_dynamic_cast<T &>(*this);
204  }
205 
207  inline void add(const FMatrixCPtr &f, const std::vector<bool> *which = NULL)
208  {
209  return add(*f, which);
210  }
211 #endif
212 
213 
214 
215 
216 };
217 
218 
223 {
224 public:
225  virtual ~AbstractSpace() {}
230  virtual void load(const picojson::object &json) = 0;
231 
233  virtual picojson::object save() const = 0;
234 
236  virtual const std::string & name() const = 0;
237 };
238 
239 
244 template<typename Scalar>
245 class SpaceBase : public AbstractSpace
246 {
247 public:
248  KQP_SCALAR_TYPEDEFS(Scalar);
249 
250 
251  static int &counter() { static int counter = 0; return counter; }
252 
253  virtual ~SpaceBase() {
254  KQP_HLOG_DEBUG_F("Destroying %s (%d / counter = %d)", %KQP_DEMANGLE(*this) %this %--counter());
255  }
256 
257  SpaceBase() {
258  KQP_HLOG_DEBUG_F("Creating %s (%d / counter = %d)", %KQP_DEMANGLE(*this) %this %++counter());
259  }
260 
262  virtual Index dimension() const = 0;
263 
264 
265 
267  virtual boost::shared_ptr< SpaceBase<Scalar> > copy() const = 0;
268 
270  virtual const ScalarMatrix &k(const FMatrixBase &mX) const = 0;
271 
273  virtual ScalarMatrix k(const FMatrixBase &mX, const ScalarAltMatrix &mY, const RealAltVector &mD) const
274  {
275  return mD.asDiagonal() * mY.adjoint() * k(mX) * mY * mD.asDiagonal();
276  }
277 
279  virtual ScalarMatrix k(const FMatrixBase &mX1, const ScalarAltMatrix &mY1, const RealAltVector &mD1,
280  const FMatrixBase &mX2, const ScalarAltMatrix &mY2, const RealAltVector &mD2) const = 0;
281 
283  inline ScalarMatrix k(const FMatrixBase &mX1, const ScalarAltMatrix &mY1,
284  const FMatrixBase &mX2, const ScalarAltMatrix &mY2) const
285  {
286  return k(mX1, mY1, RealVector::Ones(mY1.cols()), mX2, mY2, RealVector::Ones(mY2.cols()));
287  }
288 
290  inline ScalarMatrix k(const FMatrixBase &mX1, const FMatrixBase &mX2) const
291  {
292  return k(mX1, Eigen::Identity<Scalar>(mX1.size(), mX1.size()), RealVector::Ones(mX1.size()),
293  mX2, Eigen::Identity<Scalar>(mX2.size(), mX2.size()), RealVector::Ones(mX2.size()));
294  }
295 
307  virtual void updatePartials(Real /*alpha*/, std::vector<Real> &/*partials*/, int /*offset*/,
308  const std::vector< KernelValues<Scalar> > &/* kernelValues */, int /* kOffset */, int /* mode */) const {}
309 
310  void updatePartials(Real alpha, std::vector<Real> &partials,
311  const std::vector< KernelValues<Scalar> > &kernelValues, int mode) const {
312  updatePartials(alpha, partials, 0, kernelValues, 0, mode);
313  }
314 
315  virtual void update(std::vector< KernelValues<Scalar> > &, int /* kOffset */ = 0) const {}
316 
320  virtual int numberOfParameters(bool /*freeParameters*/) const
321  {
322  return 0;
323  }
324 
325  virtual void getBounds(bool /*onlyFreeParameters*/, std::vector<Real> &/*lower*/, std::vector<Real> &/*upper*/, int /*offset*/ = 0) const
326  {
327  }
328 
329  virtual void getParameters(bool /*onlyFreeParameters*/, std::vector<Real> &/*parameters*/, int /*offset*/ = 0) const
330  {
331  }
332 
333  virtual void setParameters(bool /*onlyFreeParameters*/, const std::vector<Real> &/*parameters*/, int /*offset*/ = 0)
334  {
335  }
336 
337  virtual int getNumberOfConstraints(bool /*onlyFreeParameters*/) const
338  {
339  return 0;
340  }
341 
345  virtual void getConstraints(bool /*onlyFreeParameters*/, std::vector<Real> &/*constraintValues*/, int /*offset*/ = 0) const
346  {
347  }
348 
349 
350 
354  virtual int numberOfKernelValues() const
355  {
356  return 1;
357  }
358 
359  virtual picojson::object save() const override {
360  picojson::object json;
361  json["name"] = picojson::value(name());
362  json["scalar"] = picojson::value(ScalarInfo<Scalar>::name());
363  return json;
364  }
365 
367  virtual FMatrixBasePtr newMatrix() const = 0;
368 
370  virtual FMatrixBasePtr linearCombination(const FMatrixBase &, const ScalarAltMatrix &, Scalar , const FMatrixBase *, const ScalarAltMatrix *, Scalar) const
371  {
372  KQP_THROW_EXCEPTION_F(illegal_operation_exception, "Cannot compute the linear combination in feature space [%s]", % KQP_DEMANGLE(*this));
373  }
374 
375 #ifndef SWIG
376 
377  template<typename T> inline const T &as() const
378  {
379  return kqp::our_dynamic_cast<T &>(*this);
380  }
381  template<typename T> inline T &as()
382  {
383  return kqp::our_dynamic_cast<T &>(*this);
384  }
385 #endif
386 
387  template<typename T> inline bool castable_as() const
388  {
389  return kqp::our_dynamic_cast<T *>(this) != 0;
390  }
391 
392 #ifndef SWIG
393 
394  inline const ScalarMatrix &k(const FMatrixCPtr &mX) const
395  {
396  return k(*mX);
397  }
398 
400  inline ScalarMatrix k(const FMatrixCPtr &mX, const ScalarAltMatrix &mY, const RealAltVector &mD) const
401  {
402  return k(*mX, mY, mD);
403  }
404 
406  inline ScalarMatrix k(const FMatrixCPtr &mX, const ScalarAltMatrix &mY) const
407  {
408  return k(*mX, mY, RealVector::Ones(mY.cols()));
409  }
410 
412  inline ScalarMatrix k(const FMatrixBase &mX, const ScalarAltMatrix &mY) const
413  {
414  return k(mX, mY, RealVector::Ones(mY.cols()));
415  }
416 
418  inline ScalarMatrix k(const FMatrixCPtr &mX1, const ScalarAltMatrix &mY1, const RealAltVector &mD1,
419  const FMatrixCPtr &mX2, const ScalarAltMatrix &mY2, const RealAltVector &mD2) const
420  {
421  return k(*mX1, mY1, mD1, *mX2, mY2, mD2);
422  }
423 
424 
426  inline ScalarMatrix k(const FMatrixCPtr &mX1, const ScalarAltMatrix &mY1, const FMatrixCPtr &mX2, const ScalarAltMatrix &mY2) const
427  {
428  return k(mX1, mY1, RealVector::Ones(mY1.cols()), mX2, mY2, RealVector::Ones(mY2.cols()));
429  }
430 
432  inline ScalarMatrix k(const FMatrixCPtr &mX1, const FMatrixCPtr &mX2) const
433  {
434  return k(mX1, IdentityScalarMatrix(mX1->size(), mX1->size()), RealVector::Ones(mX1->size()),
435  mX2, IdentityScalarMatrix(mX2->size(), mX2->size()), RealVector::Ones(mX2->size()));
436  }
437 #endif
438 
439 
441  inline bool canLinearlyCombine() const
442  {
443  return m_useLinearCombination && _canLinearlyCombine();
444  };
445 
446 
452  inline FMatrixPtr linearCombination(const FMatrixBase &mX, const ScalarAltMatrix &mA, Scalar alpha, const FMatrixBase &mY, const ScalarAltMatrix &mB, Scalar beta) const
453  {
454  // Check for correctedness
455  if (mA.rows() != mX.size())
456  KQP_THROW_EXCEPTION_F(out_of_bound_exception, "Cannot linearly combine with a matrix with %d rows (we have %d pre-images)", % mA.rows() % mX.size());
457 
458  // If we have no columns in A, then return an empty feature matrix
459  if (mA.cols() == 0)
460  return newMatrix();
461 
462  // Call the derived
463  return FMatrixPtr(linearCombination(mX, mA, alpha, &mY, &mB, beta));
464  }
465 
466 
472  inline FMatrixPtr linearCombination(const FMatrixBase &mX, const ScalarAltMatrix &mA, Scalar alpha = (Scalar)1) const {
473  // Check for correctedness
474  if (mA.rows() != mX.size())
475  KQP_THROW_EXCEPTION_F(out_of_bound_exception, "Cannot linearly combine with a matrix with %d rows (we have %d pre-images)", % mA.rows() % mX.size());
476 
477  // If we have no columns in A, then return an empty feature matrix
478  if (mA.cols() == 0)
479  return newMatrix();
480 
481  // Call the derived
482  return FMatrixPtr(linearCombination(mX, mA, alpha, NULL, NULL, 0));
483  }
484 
485 #ifndef SWIG
486 
487  inline FMatrixPtr linearCombination(const FMatrixCPtr &mX, const ScalarAltMatrix &mA, Scalar alpha = (Scalar)1) const {
488  return linearCombination(*mX, mA, alpha);
489  }
490 #endif
491 
493  void setUseLinearCombination(bool flag)
494  {
495  m_useLinearCombination = flag;
496  }
497 
498 
499 
500  std::string demangle() const {
501  return KQP_DEMANGLE(*this);
502  }
503 
504 protected:
506  virtual bool _canLinearlyCombine() const
507  {
508  return false;
509  }
510 
511 private:
512 
514  bool m_useLinearCombination;
515 };
516 
517 
518 }
519 
520 #endif