Skip to content

[src] Add Vector strides, beginning draft of tensor stuff #3120

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 5 commits into from
Mar 15, 2019
Merged

[src] Add Vector strides, beginning draft of tensor stuff #3120

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
fc2423d
[src] starting to refactor matrix/vector; drafting tensor stuff
danpovey Mar 15, 2019
0688a87
[src] removing matrix column stride (not supported by BLAS); various …
danpovey Mar 15, 2019
dda4479
[src] some updates to kaldi-vector.cc RE column stride, etc.
danpovey Mar 15, 2019
4ded1fc
[src] Support strides for type Vector
danpovey Mar 15, 2019
d32f22e
[egs] Revert unwanted change to cmd.sh
danpovey Mar 15, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
56 changes: 32 additions & 24 deletions src/matrix/kaldi-matrix.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -28,7 +28,7 @@
#include "matrix/compressed-matrix.h"
#include "matrix/sparse-matrix.h"

static_assert(int(kaldi::kNoTrans) == int(CblasNoTrans) && int(kaldi::kTrans) == int(CblasTrans),
static_assert(int(kaldi::kNoTrans) == int(CblasNoTrans) && int(kaldi::kTrans) == int(CblasTrans),
"kaldi::kNoTrans and kaldi::kTrans must be equal to the appropriate CBLAS library constants!");

namespace kaldi {
Expand Down Expand Up @@ -117,10 +117,10 @@ template<>
template<>
void MatrixBase<float>::AddVecVec(const float alpha,
const VectorBase<float> &a,
const VectorBase<float> &rb) {
KALDI_ASSERT(a.Dim() == num_rows_ && rb.Dim() == num_cols_);
cblas_Xger(a.Dim(), rb.Dim(), alpha, a.Data(), 1, rb.Data(),
1, data_, stride_);
const VectorBase<float> &b) {
KALDI_ASSERT(a.Dim() == num_rows_ && b.Dim() == num_cols_);
cblas_Xger(a.Dim(), b.Dim(), alpha, a.Data(), a.Stride(),
b.Data(), b.Stride(), data_, stride_);
}

template<typename Real>
Expand All @@ -132,15 +132,18 @@ void MatrixBase<Real>::AddVecVec(const Real alpha,
if (num_rows_ * num_cols_ > 100) { // It's probably worth it to allocate
// temporary vectors of the right type and use BLAS.
Vector<Real> temp_a(a), temp_b(b);
cblas_Xger(num_rows_, num_cols_, alpha, temp_a.Data(), 1,
temp_b.Data(), 1, data_, stride_);
cblas_Xger(num_rows_, num_cols_, alpha,
temp_a.Data(), temp_a.Stride(),
temp_b.Data(), temp_b.Stride(),
data_, stride_);
} else {
const OtherReal *a_data = a.Data(), *b_data = b.Data();
MatrixIndexT a_stride = a.Stride(), b_stride = b.Stride();
Real *row_data = data_;
for (MatrixIndexT i = 0; i < num_rows_; i++, row_data += stride_) {
BaseFloat alpha_ai = alpha * a_data[i];
BaseFloat alpha_ai = alpha * a_data[i * a_stride];
for (MatrixIndexT j = 0; j < num_cols_; j++)
row_data[j] += alpha_ai * b_data[j];
row_data[j] += alpha_ai * b_data[j * b_stride];
}
}
}
Expand All @@ -159,11 +162,11 @@ template<>
template<>
void MatrixBase<double>::AddVecVec(const double alpha,
const VectorBase<double> &a,
const VectorBase<double> &rb) {
KALDI_ASSERT(a.Dim() == num_rows_ && rb.Dim() == num_cols_);
const VectorBase<double> &b) {
KALDI_ASSERT(a.Dim() == num_rows_ && b.Dim() == num_cols_);
if (num_rows_ == 0) return;
cblas_Xger(a.Dim(), rb.Dim(), alpha, a.Data(), 1, rb.Data(),
1, data_, stride_);
cblas_Xger(a.Dim(), b.Dim(), alpha, a.Data(), a.Stride(),
b.Data(), b.Stride(), data_, stride_);
}

template<typename Real>
Expand Down Expand Up @@ -538,7 +541,7 @@ void MatrixBase<Real>::AddMatSmat(Real alpha, const MatrixBase<Real> &A,
// pass stride to write a column as matrices are stored in row major order.
cblas_Xaxpy(this_num_rows, alpha_B_jk, a_col_k, A.stride_,
this_col_j, this->stride_);
//for (MatrixIndexT i = 0; i < this_num_rows; ++i)
//for (MatrixIndexT i = 0; i < this_num_rows; ++i)
// this_col_j[i*this->stride_] += alpha_B_jk * a_col_k[i*A.stride_];
}
}
Expand Down Expand Up @@ -591,8 +594,10 @@ void MatrixBase<Real>::AddDiagVecMat(
if (transM == kTrans) std::swap(M_row_stride, M_col_stride);
Real *data = data_;
const Real *Mdata = M.Data(), *vdata = v.Data();
MatrixIndexT v_stride = v.Stride();
if (num_rows_ == 0) return;
for (MatrixIndexT i = 0; i < num_rows; i++, data += stride, Mdata += M_row_stride, vdata++)
for (MatrixIndexT i = 0; i < num_rows;
i++, data += stride, Mdata += M_row_stride, vdata += v_stride)
cblas_Xaxpy(num_cols, alpha * *vdata, Mdata, M_col_stride, data, 1);
}

Expand Down Expand Up @@ -623,10 +628,11 @@ void MatrixBase<Real>::AddMatDiagVec(

Real *data = data_;
const Real *Mdata = M.Data(), *vdata = v.Data();
MatrixIndexT v_stride = v.Stride();
if (num_rows_ == 0) return;
for (MatrixIndexT i = 0; i < num_rows; i++){
for(MatrixIndexT j = 0; j < num_cols; j ++ ){
data[i*stride + j] += alpha * vdata[j] * Mdata[i*M_row_stride + j*M_col_stride];
data[i*stride + j] += alpha * vdata[j * v_stride] * Mdata[i*M_row_stride + j*M_col_stride];
}
}
}
Expand Down Expand Up @@ -658,7 +664,8 @@ void MatrixBase<Real>::AddMatMatElements(const Real alpha,
template<typename Real>
void MatrixBase<Real>::LapackGesvd(VectorBase<Real> *s, MatrixBase<Real> *U_in,
MatrixBase<Real> *V_in) {
KALDI_ASSERT(s != NULL && U_in != this && V_in != this);
KALDI_ASSERT(s != NULL && U_in != this && V_in != this &&
s->Stride() == 1);

Matrix<Real> tmpU, tmpV;
if (U_in == NULL) tmpU.Resize(this->num_rows_, 1); // work-space if U_in empty.
Expand Down Expand Up @@ -1786,7 +1793,7 @@ void MatrixBase<Real>::DestructiveSvd(VectorBase<Real> *s, MatrixBase<Real> *U,
// Throws exception on error.

KALDI_ASSERT(num_rows_>=num_cols_ && "Svd requires that #rows by >= #cols."); // For compatibility with JAMA code.
KALDI_ASSERT(s->Dim() == num_cols_); // s should be the smaller dim.
KALDI_ASSERT(s->Dim() == num_cols_ && s->Stride() == 1); // s should be the smaller dim.
KALDI_ASSERT(U == NULL || (U->num_rows_ == num_rows_&&U->num_cols_ == num_cols_));
KALDI_ASSERT(Vt == NULL || (Vt->num_rows_ == num_cols_&&Vt->num_cols_ == num_cols_));

Expand Down Expand Up @@ -1992,27 +1999,28 @@ void MatrixBase<Real>::OrthogonalizeRows() {
// symmetric positive definite).

template<typename Real>
void MatrixBase<Real>::SymPosSemiDefEig(VectorBase<Real> *rs, MatrixBase<Real> *rU, Real check_thresh) // e.g. check_thresh = 0.001
void MatrixBase<Real>::SymPosSemiDefEig(VectorBase<Real> *s, MatrixBase<Real> *U, Real check_thresh) // e.g. check_thresh = 0.001
{
const MatrixIndexT D = num_rows_;

KALDI_ASSERT(num_rows_ == num_cols_);
KALDI_ASSERT(IsSymmetric() && "SymPosSemiDefEig: expecting input to be symmetrical.");
KALDI_ASSERT(rU->num_rows_ == D && rU->num_cols_ == D && rs->Dim() == D);
KALDI_ASSERT(U->num_rows_ == D && U->num_cols_ == D && s->Dim() == D &&
s->Stride() == 1);

Matrix<Real> Vt(D, D);
Svd(rs, rU, &Vt);
Svd(s, U, &Vt);

// First just zero any singular values if the column of U and V do not have +ve dot product--
// this may mean we have small negative eigenvalues, and if we zero them the result will be closer to correct.
for (MatrixIndexT i = 0;i < D;i++) {
Real sum = 0.0;
for (MatrixIndexT j = 0;j < D;j++) sum += (*rU)(j, i) * Vt(i, j);
if (sum < 0.0) (*rs)(i) = 0.0;
for (MatrixIndexT j = 0;j < D;j++) sum += (*U)(j, i) * Vt(i, j);
if (sum < 0.0) (*s)(i) = 0.0;
}

{
Matrix<Real> tmpU(*rU); Vector<Real> tmps(*rs); tmps.ApplyPow(0.5);
Matrix<Real> tmpU(*U); Vector<Real> tmps(*s); tmps.ApplyPow(0.5);
tmpU.MulColsVec(tmps);
SpMatrix<Real> tmpThis(D);
tmpThis.AddMat2(1.0, tmpU, kNoTrans, 0.0);
Expand Down
26 changes: 16 additions & 10 deletions src/matrix/kaldi-matrix.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -59,13 +59,13 @@ class MatrixBase {
friend class SparseMatrix<float>;
friend class SparseMatrix<double>;

/// Returns number of rows (or zero for emtpy matrix).
/// Returns number of rows (or zero for empty matrix).
inline MatrixIndexT NumRows() const { return num_rows_; }

/// Returns number of columns (or zero for emtpy matrix).
inline MatrixIndexT NumCols() const { return num_cols_; }

/// Stride (distance in memory between each row). Will be >= NumCols.
/// Stride (distance in memory between each row). Must be >= NumCols().
inline MatrixIndexT Stride() const { return stride_; }

/// Returns size in bytes of the data held by the matrix.
Expand Down Expand Up @@ -183,18 +183,20 @@ class MatrixBase {

/* Accessing of sub-parts of the matrix. */

/// Return specific row of matrix [const].
inline const SubVector<Real> Row(MatrixIndexT i) const {
/// Return specific row of matrix. Warning: this can get
/// around const constraints.
inline SubVector<Real> Row(MatrixIndexT i) const {
KALDI_ASSERT(static_cast<UnsignedMatrixIndexT>(i) <
static_cast<UnsignedMatrixIndexT>(num_rows_));
return SubVector<Real>(data_ + (i * stride_), NumCols());
return SubVector<Real>(data_ + (i * stride_), num_cols_);
}

/// Return specific row of matrix.
inline SubVector<Real> Row(MatrixIndexT i) {
/// Return specific column of matrix. Warning: this can get
/// around const constraints.
inline const SubVector<Real> Col(MatrixIndexT i) const {
KALDI_ASSERT(static_cast<UnsignedMatrixIndexT>(i) <
static_cast<UnsignedMatrixIndexT>(num_rows_));
return SubVector<Real>(data_ + (i * stride_), NumCols());
static_cast<UnsignedMatrixIndexT>(num_cols_));
return SubVector<Real>(data_ + i, num_rows_, stride_);
}

/// Return a sub-part of matrix.
Expand Down Expand Up @@ -406,14 +408,17 @@ class MatrixBase {
Null pointers for U and/or Vt at input mean we do not want that output. We
expect that S.Dim() == m, U is either NULL or m by n,
and v is either NULL or n by n.
The singular values are not sorted (use SortSvd for that). */
The singular values are not sorted (use SortSvd for that).
Requires that s->Stride() == 1.
*/
void DestructiveSvd(VectorBase<Real> *s, MatrixBase<Real> *U,
MatrixBase<Real> *Vt); // Destroys calling matrix.

/// Compute SVD (*this) = U diag(s) Vt. Note that the V in the call is already
/// transposed; the normal formulation is U diag(s) V^T.
/// Null pointers for U or V mean we don't want that output (this saves
/// compute). The singular values are not sorted (use SortSvd for that).
/// Requires that s->Stride() == 1.
void Svd(VectorBase<Real> *s, MatrixBase<Real> *U,
MatrixBase<Real> *Vt) const;
/// Compute SVD but only retain the singular values.
Expand Down Expand Up @@ -531,6 +536,7 @@ class MatrixBase {
* positive semi-definite (check_thresh controls how stringent the check is;
* set it to 2 to ensure it won't ever complain, but it will zero out negative
* dimensions in your matrix.
* Requires s->Stride() == 1.
*/
void SymPosSemiDefEig(VectorBase<Real> *s, MatrixBase<Real> *P,
Real check_thresh = 0.001);
Expand Down
Loading