[shark] 66/79: removed almost all uses of axpy_prod
Ghislain Vaillant
ghisvail-guest at moszumanska.debian.org
Thu Nov 26 15:41:23 UTC 2015
This is an automated email from the git hooks/post-receive script.
ghisvail-guest pushed a commit to branch master
in repository shark.
commit 62dec7c052c168c91af499a1fbd356514395f750
Author: Oswin Krause <oswin.krause at di.ku.dk>
Date: Wed Oct 28 14:40:34 2015 +0100
removed almost all uses of axpy_prod
---
.../shark/Algorithms/Trainers/KernelSGDTrainer.h | 3 +-
include/shark/LinAlg/Impl/Cholesky.inl | 6 +-
include/shark/LinAlg/Impl/pivotingRQ.inl | 5 +-
include/shark/LinAlg/Impl/solveSystem.inl | 32 +--
include/shark/LinAlg/Metrics.h | 4 +-
include/shark/Models/Autoencoder.h | 10 +-
include/shark/Models/ConvexCombination.h | 7 +-
include/shark/Models/FFNet.h | 10 +-
include/shark/Models/Impl/LinearModel.inl | 264 ---------------------
include/shark/Models/Kernels/CSvmDerivative.h | 7 +-
include/shark/Models/Kernels/GaussianRbfKernel.h | 2 +-
include/shark/Models/Kernels/KernelExpansion.h | 2 +-
include/shark/Models/Kernels/LinearKernel.h | 4 +-
include/shark/Models/Kernels/MonomialKernel.h | 8 +-
include/shark/Models/Kernels/PolynomialKernel.h | 9 +-
include/shark/Models/LinearModel.h | 6 +-
include/shark/Models/TiedAutoencoder.h | 12 +-
.../Benchmarks/RotatedErrorFunction.h | 3 +-
.../Impl/SparseAutoencoderError.inl | 2 +-
.../Unsupervised/RBM/Impl/AverageEnergyGradient.h | 4 +-
.../Unsupervised/RBM/Neuronlayers/GaussianLayer.h | 3 +-
include/shark/Unsupervised/RBM/RBM.h | 4 +-
src/Algorithms/GradientDescent/BFGS.cpp | 6 +-
src/Algorithms/LinearRegression.cpp | 2 +-
src/Algorithms/PCA.cpp | 8 +-
src/Models/OnlineRNNet.cpp | 10 +-
src/Models/RBFLayer.cpp | 3 +-
27 files changed, 68 insertions(+), 368 deletions(-)
diff --git a/include/shark/Algorithms/Trainers/KernelSGDTrainer.h b/include/shark/Algorithms/Trainers/KernelSGDTrainer.h
index ba0f75d..b880a1a 100644
--- a/include/shark/Algorithms/Trainers/KernelSGDTrainer.h
+++ b/include/shark/Algorithms/Trainers/KernelSGDTrainer.h
@@ -174,9 +174,8 @@ public:
const double eta = 1.0 / (lambda * (iter + ell));
// compute prediction
- f_b.clear();
K.row(b, kernelRow);
- axpy_prod(trans(alpha), kernelRow, f_b, false, alphaScale);
+ noalias(f_b) = alphaScale * prod(trans(alpha), kernelRow);
if(m_offset) noalias(f_b) += model.offset();
// stochastic gradient descent (SGD) step
diff --git a/include/shark/LinAlg/Impl/Cholesky.inl b/include/shark/LinAlg/Impl/Cholesky.inl
index c5ee286..54046b2 100644
--- a/include/shark/LinAlg/Impl/Cholesky.inl
+++ b/include/shark/LinAlg/Impl/Cholesky.inl
@@ -133,11 +133,9 @@ std::size_t shark::blas::pivotingCholeskyDecompositionInPlace(
//(L1L1T)^(j)+...+(Lj-1Lj-1)^(j)
//=L1*L1j+L2*L2j+L3*L3j...
//which is a matrix-vector product
- axpy_prod(
+ unitTriangularColumn(Lk,j) -= prod(
LkBlocked.lowerLeft(),
- subrange(row(Lk,j),0,j),
- unitTriangularColumn(Lk,j),
- false,-1.0
+ subrange(row(Lk,j),0,j)
);
}
unitTriangularColumn(Lk,j) /= Lk(j,j);
diff --git a/include/shark/LinAlg/Impl/pivotingRQ.inl b/include/shark/LinAlg/Impl/pivotingRQ.inl
index 4190574..52415c2 100644
--- a/include/shark/LinAlg/Impl/pivotingRQ.inl
+++ b/include/shark/LinAlg/Impl/pivotingRQ.inl
@@ -161,9 +161,8 @@ std::size_t shark::blas::pivotingRQ
solveTriangularSystemInPlace<SolveAXB,lower>(trans(T), InvTU);
matrixQ().resize(n,n);
//now Compute U^T temp = U^T T^-1 U
- axpy_prod(trans(rows(U,0,rank)),InvTU,matrixQ);
- matrixQ()*=-1;
- matrixQ()+=RealIdentityMatrix(n);
+ noalias(matrixQ) = -prod(trans(rows(U,0,rank)),InvTU);
+ diag(matrixQ) += n;
//testing algorithm
// matrixQ().resize(n,n);
diff --git a/include/shark/LinAlg/Impl/solveSystem.inl b/include/shark/LinAlg/Impl/solveSystem.inl
index 8bfb77c..822358d 100644
--- a/include/shark/LinAlg/Impl/solveSystem.inl
+++ b/include/shark/LinAlg/Impl/solveSystem.inl
@@ -276,9 +276,8 @@ void shark::blas::solveSymmSemiDefiniteSystemInPlace(
RealMatrix LTL(rank,rank);
symm_prod(trans(L),LTL);
- //compute z= L^TB
- RealMatrix Z(rank,B().size2());
- axpy_prod(trans(L),B,Z);
+ //compute Z= L^TB
+ RealMatrix Z = prod(trans(L),B,Z);
//compute cholesky factor of L^TL
RealMatrix LTLcholesky(rank,rank);
@@ -287,7 +286,7 @@ void shark::blas::solveSymmSemiDefiniteSystemInPlace(
//A'b = L(L^TL)^-1(L^TL)^-1z
solveTriangularCholeskyInPlace<System>(LTLcholesky,Z);
solveTriangularCholeskyInPlace<System>(LTLcholesky,Z);
- axpy_prod(L,Z,B);
+ noalias(B) = prod(L,Z);
}else{
//complex case.
//X=BL(L^TL)^-1(L^TL)^-1 L^T
@@ -295,8 +294,7 @@ void shark::blas::solveSymmSemiDefiniteSystemInPlace(
symm_prod(trans(L),LTL);
//compute z= L^TB
- RealMatrix Z(B().size1(),rank);
- axpy_prod(B,L,Z);
+ RealMatrix Z = prod(B,L);
//compute cholesky factor of L^TL
RealMatrix LTLcholesky(rank,rank);
@@ -305,7 +303,7 @@ void shark::blas::solveSymmSemiDefiniteSystemInPlace(
//A'b = L(L^TL)^-1(L^TL)^-1z
solveTriangularCholeskyInPlace<System>(LTLcholesky,Z);
solveTriangularCholeskyInPlace<System>(LTLcholesky,Z);
- axpy_prod(Z,trans(L),B);
+ noalias(B) = prod(Z,trans(L));
}
//finally swap back into the unpermuted coordinate system
if(System::left)
@@ -328,8 +326,7 @@ void shark::blas::generalSolveSystemInPlace(
//Ax=b => A^TAx=A'b => x= A'b = (A^TA)' Ab
// with z = Ab => (A^TA) x= z
//compute A^TA
- RealMatrix ATA(n,n);
- axpy_prod(trans(A),A,ATA);
+ RealMatrix ATA = prod(trans(A),A);
//compute z=Ab
RealVector z = prod(trans(A),b);
@@ -344,8 +341,7 @@ void shark::blas::generalSolveSystemInPlace(
//x^TA=b^T => x^TAA'=b^TA' => x^T= b^TA' = b^TA^T(AA^T)'
// with z = Ab => x^T(AA^T) = z^T
//compute AAT
- RealMatrix AAT(m,m);
- axpy_prod(A,trans(A),AAT);
+ RealMatrix AAT = prod(A,trans(A));
//compute z=Ab
RealVector z = prod(A,b);
@@ -370,12 +366,9 @@ void shark::blas::generalSolveSystemInPlace(
//AX=B => A'AX=A'B => X= A'B = (A^TA)' A^TB
// with Z = A^TB => (A^TA) X= Z
//compute A^TA
- RealMatrix ATA(n,n);
- axpy_prod(trans(A),A,ATA);
+ RealMatrix ATA = prod(trans(A),A);
- //compute Z=AB
- RealMatrix Z(n,B().size2());
- axpy_prod(trans(A),B,Z);
+ RealMatrix Z = prod(trans(A),B);
//call recursively for the quadratic case
solveSymmSemiDefiniteSystemInPlace<System>(ATA,Z);
@@ -387,12 +380,9 @@ void shark::blas::generalSolveSystemInPlace(
//~ //XA=B => XAA'=BA' => X = BA' = BA^T(AA^T)'
//~ // with Z = BA^T => X(AA^T) = Z
//~ //compute AAT
- RealMatrix AAT(m,m);
- axpy_prod(A,trans(A),AAT);
+ RealMatrix AAT = prod(A,trans(A));
- //compute z=Ab
- RealMatrix Z(B().size1(),m);
- axpy_prod(B,trans(A),Z);
+ RealMatrix Z= prod(B,trans(A));
//call recursively for the quadratic case
solveSymmSemiDefiniteSystemInPlace<System>(AAT,Z);
diff --git a/include/shark/LinAlg/Metrics.h b/include/shark/LinAlg/Metrics.h
index abe2f6e..a2885ba 100644
--- a/include/shark/LinAlg/Metrics.h
+++ b/include/shark/LinAlg/Metrics.h
@@ -239,14 +239,14 @@ namespace detail{
typedef typename Result::value_type value_type;
std::size_t sizeX=X.size1();
std::size_t sizeY=Y.size1();
+ ensureSize(distances,X.size1(),Y.size1());
if(sizeX < 10 || sizeY<10){
distanceSqrBlockBlockRowWise(X,Y,distances);
return;
}
//fast blockwise iteration
//uses: (a-b)^2 = a^2 -2ab +b^2
- axpy_prod(X,trans(Y),distances);
- distances*=-2;
+ noalias(distances) = -2*prod(X,trans(Y));
//first a^2+b^2
vector<value_type> ySqr(sizeY);
for(std::size_t i = 0; i != sizeY; ++i){
diff --git a/include/shark/Models/Autoencoder.h b/include/shark/Models/Autoencoder.h
index 8e15df5..0a2912b 100644
--- a/include/shark/Models/Autoencoder.h
+++ b/include/shark/Models/Autoencoder.h
@@ -170,8 +170,7 @@ public:
std::size_t numOutputs = encoderMatrix().size1();
outputs.resize(numPatterns,numOutputs);
outputs.clear();
- axpy_prod(patterns,trans(encoderMatrix()),outputs);
- noalias(outputs) += repeat(hiddenBias(),numPatterns);
+ noalias(outputs) = prod(patterns,trans(encoderMatrix())) + repeat(hiddenBias(),numPatterns);
noalias(outputs) = m_hiddenNeuron(outputs);
}
else{//hidden->output
@@ -179,8 +178,7 @@ public:
std::size_t numOutputs = decoderMatrix().size1();
outputs.resize(numPatterns,numOutputs);
outputs.clear();
- axpy_prod(patterns,trans(decoderMatrix()),outputs);
- noalias(outputs) += repeat(outputBias(),numPatterns);
+ noalias(outputs) = prod(patterns,trans(decoderMatrix()),outputs) + repeat(outputBias(),numPatterns);
noalias(outputs) = m_outputNeuron(outputs);
}
}
@@ -304,7 +302,7 @@ private:
noalias(outputDelta) *= m_outputNeuron.derivative(s.outputResponses);
hiddenDelta.resize(outputDelta.size1(),numberOfHiddenNeurons());
- axpy_prod(outputDelta,decoderMatrix(),hiddenDelta,true);
+ noalias(hiddenDelta) = prod(outputDelta,decoderMatrix());
noalias(hiddenDelta) *= m_hiddenNeuron.derivative(s.hiddenResponses);
}
@@ -313,7 +311,7 @@ private:
)const{
computeDelta(state,outputDelta,hiddenDelta);
inputDelta.resize(outputDelta.size1(),inputSize());
- axpy_prod(hiddenDelta,encoderMatrix(),inputDelta,true);
+ noalias(inputDelta) = prod(hiddenDelta,encoderMatrix());
}
void computeParameterDerivative(
diff --git a/include/shark/Models/ConvexCombination.h b/include/shark/Models/ConvexCombination.h
index c3a60cf..58c0aaa 100644
--- a/include/shark/Models/ConvexCombination.h
+++ b/include/shark/Models/ConvexCombination.h
@@ -142,7 +142,7 @@ public:
/// Evaluate the model: output = w * input
void eval(BatchInputType const& inputs, BatchOutputType& outputs)const{
outputs.resize(inputs.size1(),m_w.size1());
- axpy_prod(inputs,trans(m_w),outputs);
+ noalias(outputs) = prod(inputs,trans(m_w));
}
/// Evaluate the model: output = w *input
void eval(BatchInputType const& inputs, BatchOutputType& outputs, State& state)const{
@@ -162,8 +162,7 @@ public:
//derivative is
//sum_i sum_j c_ij sum_k x_ik grad_q w_jk= sum_k sum_j grad_q w_jk (sum_i c_ij x_ik)
//and we set d_jk=sum_i c_ij x_ik => d = C^TX
- RealMatrix d(outputSize(),inputSize());
- axpy_prod(trans(coefficients), patterns,d);
+ RealMatrix d = prod(trans(coefficients), patterns);
//use the same drivative as in the softmax model
for(std::size_t i = 0; i != outputSize(); ++i){
@@ -185,7 +184,7 @@ public:
SIZE_CHECK(coefficients.size1() == patterns.size1());
derivative.resize(patterns.size1(),inputSize());
- axpy_prod(coefficients,m_w,derivative);
+ noalias(derivative) = prod(coefficients,m_w);
}
/// From ISerializable
diff --git a/include/shark/Models/FFNet.h b/include/shark/Models/FFNet.h
index 11ff515..b77b114 100644
--- a/include/shark/Models/FFNet.h
+++ b/include/shark/Models/FFNet.h
@@ -276,7 +276,7 @@ public:
//calculate activation. first compute the linear part and the optional bias and then apply
// the non-linearity
- axpy_prod(patterns,trans(layerMatrix(layer)),outputs);
+ noalias(outputs) = prod(patterns,trans(layerMatrix(layer)));
if(!bias().empty()){
noalias(outputs) += repeat(bias(layer),numPatterns);
}
@@ -323,7 +323,7 @@ public:
//calculate activation. first compute the linear part and the optional bias and then apply
// the non-linearity
- axpy_prod(weights,input,responses);
+ noalias(responses) = prod(weights,input);
if(!bias().empty()){
//the bias of the layer is shifted as input units can not have bias.
ConstRealVectorRange bias = subrange(m_bias,beginNeuron-inputSize(),endNeuron-inputSize());
@@ -337,7 +337,7 @@ public:
else {
//add shortcuts if necessary
if(m_inputOutputShortcut.size1() != 0){
- axpy_prod(m_inputOutputShortcut,trans(patterns),responses,false);
+ noalias(responses) += prod(m_inputOutputShortcut,trans(patterns));
}
noalias(responses) = m_outputNeuron(responses);
}
@@ -621,7 +621,7 @@ private:
RealSubMatrix layerDeltaInput = rows(delta,endNeuron,endNeuron+weights.size2());
ConstRealSubMatrix layerResponse = rows(s.responses,beginNeuron,endNeuron);
- axpy_prod(weights,layerDeltaInput,layerDelta,false);//add the values to the maybe non-empty delta part
+ noalias(layerDelta) += prod(weights,layerDeltaInput);//add the values to the maybe non-empty delta part
if(layer != 0){
noalias(layerDelta) *= m_hiddenNeuron.derivative(layerResponse);
}
@@ -632,7 +632,7 @@ private:
//add the shortcut deltas if necessary
if(inputOutputShortcut().size1() != 0)
- axpy_prod(trans(inputOutputShortcut()),outputDelta,rows(delta,0,inputSize()),false);
+ noalias(rows(delta,0,inputSize())) += prod(trans(inputOutputShortcut()),outputDelta);
}
void computeParameterDerivative(RealMatrix const& delta, State const& state, RealVector& gradient)const{
diff --git a/include/shark/Models/Impl/LinearModel.inl b/include/shark/Models/Impl/LinearModel.inl
deleted file mode 100644
index e8f4c38..0000000
--- a/include/shark/Models/Impl/LinearModel.inl
+++ /dev/null
@@ -1,264 +0,0 @@
-/*!
-* \brief Implements a Model using a linear function.
-*
-* \author T. Glasmachers, O. Krause
-* \date 2010-2011
-*
-* \par Copyright (c) 1999-2011:
-* Institut für Neuroinformatik<BR>
-* Ruhr-Universität Bochum<BR>
-* D-44780 Bochum, Germany<BR>
-* Phone: +49-234-32-27974<BR>
-* Fax: +49-234-32-14209<BR>
-* eMail: Shark-admin at neuroinformatik.ruhr-uni-bochum.de<BR>
-* www: http://www.neuroinformatik.ruhr-uni-bochum.de<BR>
-*
-*
-*
-* <BR><HR>
-* This file is part of Shark. This library is free software;
-* you can redistribute it and/or modify it under the terms of the
-* GNU General Public License as published by the Free Software
-* Foundation; either version 3, or (at your option) any later version.
-*
-* This library is distributed in the hope that it will be useful,
-* but WITHOUT ANY WARRANTY; without even the implied warranty of
-* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
-* GNU General Public License for more details.
-*
-* You should have received a copy of the GNU General Public License
-* along with this library; if not, see <http://www.gnu.org/licenses/>.
-*
-*/
-#include <shark/Models/AbstractModel.h>
-#include <shark/LinAlg/BLAS/Initialize.h>
-namespace shark {
-
-namespace detail {
-
-/// \brief Baseclass of LinearModelWrapper.
-/// This class is needed for the type
-/// erasure in LinearModel.
-template <class InputType, class OutputType>
-class LinearModelWrapperBase : public AbstractModel<InputType, OutputType>
-{
-protected:
- OutputType m_offset;
-
- LinearModelWrapperBase()
- { }
-
- LinearModelWrapperBase(OutputType const& offset)
- : m_offset(offset)
- { }
-
- LinearModelWrapperBase(bool offset, std::size_t outputSize)
- :m_offset(offset? outputSize:0,0.0){}
-
-public:
-
- /// \brief From INameable: return the class name.
- std::string name() const
- { return "LinearModelWrapperBase"; }
-
- bool hasOffset() const
- { return (m_offset.size() > 0); }
-
- virtual std::size_t inputSize() const = 0;
- virtual std::size_t outputSize() const = 0;
- virtual void matrixRow(std::size_t index, RealVector& row) const = 0;
- virtual void matrixColumn(std::size_t index, RealVector& column) const = 0;
- virtual void matrix(RealMatrix& mat) const = 0;
- OutputType const& offset() const
- { return m_offset; }
-
- void read(InArchive& archive){
- archive >> m_offset;
- }
- void write(OutArchive& archive) const{
- archive << m_offset;
- }
-
- //workaround for base pointer serialization
- template<class Archive>
- void serialize(Archive& ar,unsigned int version){
- ar & m_offset;
- }
-
- virtual LinearModelWrapperBase<InputType, OutputType>* clone()const = 0;
-};
-
-///\brief Implementation of the Linear Model for specific matrix types
-///
-/// The Matrix template parameter is a real matrix.
-/// It is open to sparse and dense representations.
-template <class Matrix, class InputType, class OutputType>
-class LinearModelWrapper : public LinearModelWrapperBase<InputType, OutputType>
-{
- typedef LinearModelWrapperBase<InputType, OutputType> base_type;
- typedef typename Batch<InputType>::type BatchInputType;
- typedef typename Batch<OutputType>::type BatchOutputType;
-public:
- LinearModelWrapper(){}
- LinearModelWrapper(unsigned int inputs, unsigned int outputs, bool offset)
- :base_type(offset,outputs),m_matrix(outputs,inputs,0.0){
- if(!blas::traits::IsSparse<Matrix>::value){
- base_type::m_features |= base_type::HAS_FIRST_PARAMETER_DERIVATIVE;
- }
- }
-
- LinearModelWrapper(Matrix const& matrix)
- : m_matrix(matrix){
- if(!blas::traits::IsSparse<Matrix>::value){
- base_type::m_features |= base_type::HAS_FIRST_PARAMETER_DERIVATIVE;
- }
- }
-
- LinearModelWrapper(Matrix const& matrix, OutputType const& offset)
- : base_type(offset)
- , m_matrix(matrix){
- if(!blas::traits::IsSparse<Matrix>::value){
- base_type::m_features |= base_type::HAS_FIRST_PARAMETER_DERIVATIVE;
- }
- }
-
- /// \brief From INameable: return the class name.
- std::string name() const
- { return "LinearModelWrapper"; }
-
- std::size_t inputSize() const{
- return m_matrix.size2();
- }
- std::size_t outputSize() const{
- return m_matrix.size1();
- }
- void matrix(RealMatrix& mat) const{
- mat = m_matrix;
- }
- void matrixRow(std::size_t index, RealVector& rowStorage) const{
- rowStorage = row(m_matrix, index);
- }
- void matrixColumn(std::size_t index, RealVector& columnStorage) const{
- columnStorage = column(m_matrix, index);
- }
-
- RealVector parameterVector() const{
- RealVector ret(numberOfParameters());
- if(base_type::hasOffset())
- init(ret) << toVector(m_matrix),base_type::m_offset;
- else
- init(ret) << toVector(m_matrix);
-
- return ret;
- }
-
- void setParameterVector(RealVector const& newParameters){
- SHARK_CHECK(newParameters.size() == numberOfParameters(), "[LinearModelWrapper::setParameterVector] invalid size of the parameter vector");
- if(base_type::hasOffset())
- init(newParameters) >> toVector(m_matrix),base_type::m_offset;
- else
- init(newParameters) >> toVector(m_matrix);
- }
-
- std::size_t numberOfParameters() const{
- if(blas::traits::IsSparse<Matrix>::value){
- std::size_t elements = base_type::m_offset.size();
- //works for sparse, but not very efficient for big dense matrices
- for(typename Matrix::const_iterator1 pos = m_matrix.begin1(); pos != m_matrix.end1(); ++pos){
- typename Matrix::const_iterator2 pos1Iter=pos.begin();
- typename Matrix::const_iterator2 pos2Iter=pos.end();
- elements += std::distance(pos1Iter,pos2Iter);
- }
- return elements;
- }
- else{
- return m_matrix.size1()*m_matrix.size2()+base_type::m_offset.size();
- }
- }
-
- boost::shared_ptr<State> createState()const{
- return boost::shared_ptr<State>(new EmptyState());
- }
-
- void eval(BatchInputType const& inputs, BatchOutputType& outputs)const{
- outputs.resize(inputs.size1(),m_matrix.size1()
- ,false // Oh thanks ublas for making this necessary!!! damn!
- );
- //we multiply with a set of row vectors from the left
- axpy_prod(inputs,trans(m_matrix),outputs);
- if (base_type::hasOffset()){
- noalias(outputs)+=repeat(base_type::m_offset,inputs.size1());
- }
- }
-
-
- /// \param inputs Batch of input patterns
- /// \param outputs Outputs corresponding to input batch
- /// \param state Internal state (e.g., for passing to derivative computation) of the model
- void eval(BatchInputType const& inputs, BatchOutputType& outputs, State& state)const{
- eval(inputs,outputs);
- }
-
- void weightedParameterDerivative(
- BatchInputType const& patterns, RealMatrix const& coefficients, State const& state, RealVector& gradient
- )const{
- //todo: doesn't work for sparse.
- SIZE_CHECK(coefficients.size2()==outputSize());
- SIZE_CHECK(coefficients.size1()==patterns.size1());
- SHARK_CHECK(!blas::traits::IsSparse<Matrix>::value,"Derivative for sparse matrices is not implemented.");
- gradient.resize(numberOfParameters());
- std::size_t inputs = inputSize();
- std::size_t outputs = outputSize();
- gradient.clear();
- std::size_t first = 0;
- for (std::size_t output = 0; output < outputs; output++){
- //sum_i coefficients(output,i)*pattern(i))
- axpy_prod(trans(patterns),column(coefficients,output),subrange(gradient, first, first + inputs),0.0);
- first += inputs;
- }
- if (base_type::hasOffset()){
- for(std::size_t i = 0; i != patterns.size1();++i){
- noalias(subrange(gradient, first, first + outputs))+= row(coefficients,i);
- }
- }
- }
-
- void read(InArchive& archive){
- archive >> boost::serialization::base_object<base_type>(*this);
- archive >> m_matrix;
- }
- void write(OutArchive& archive) const{
- //VS2010 workaround, produces warnings in gcc :(
- archive << boost::serialization::base_object<base_type>(
- const_cast<LinearModelWrapper<Matrix,InputType, OutputType>& >(*this)
- );
-// archive & boost::serialization::base_object<base_type >(*this);
- archive << m_matrix;
- }
-
- /**
- * \brief Versioned loading of components, calls read(...).
- */
- void load(InArchive & archive,unsigned int version){
- read(archive);
- }
-
- /**
- * \brief Versioned storing of components, calls write(...).
- */
- void save(OutArchive & archive,unsigned int version)const{
- write(archive);
- }
- BOOST_SERIALIZATION_SPLIT_MEMBER();
-
- LinearModelWrapperBase<InputType, OutputType>* clone()const{
- return new LinearModelWrapper<Matrix,InputType, OutputType>(*this);
- }
-
-protected:
- Matrix m_matrix;
-};
-
-
-}
-}
diff --git a/include/shark/Models/Kernels/CSvmDerivative.h b/include/shark/Models/Kernels/CSvmDerivative.h
index c512a79..435f546 100644
--- a/include/shark/Models/Kernels/CSvmDerivative.h
+++ b/include/shark/Models/Kernels/CSvmDerivative.h
@@ -356,14 +356,13 @@ private:
// compute the derivative of (\alpha, b) w.r.t. C
if ( m_noofBoundedSVs > 0 ) {
- axpy_prod( R, m_boundedLabels, column(m_d_alphab_d_theta,m_nkp));
+ noalias(column(m_d_alphab_d_theta,m_nkp)) = prod( R, m_boundedLabels);
}
// compute the derivative of (\alpha, b) w.r.t. the kernel parameters
for ( std::size_t k=0; k<m_nkp; k++ ) {
- RealVector sum( m_noofFreeSVs+1);
- axpy_prod( dH[k], m_freeAlphas, sum ); //sum = dH * \alpha_f
+ RealVector sum = prod( dH[k], m_freeAlphas ); //sum = dH * \alpha_f
if(m_noofBoundedSVs > 0)
- axpy_prod( dR[k], m_boundedAlphas, sum, false ); // sum += dR * \alpha_r , i.e., the C*y_g is expressed as alpha_g
+ noalias(sum) += prod( dR[k], m_boundedAlphas ); // sum += dR * \alpha_r , i.e., the C*y_g is expressed as alpha_g
//fill the remaining columns of the derivative matrix (except the last, which is for C)
noalias(column(m_d_alphab_d_theta,k)) = sum;
}
diff --git a/include/shark/Models/Kernels/GaussianRbfKernel.h b/include/shark/Models/Kernels/GaussianRbfKernel.h
index 9f47d90..eb85497 100644
--- a/include/shark/Models/Kernels/GaussianRbfKernel.h
+++ b/include/shark/Models/Kernels/GaussianRbfKernel.h
@@ -222,7 +222,7 @@ public:
gradient.resize(sizeX1,batchX1.size2());
RealMatrix W = coefficientsX2*s.expNorm;
- axpy_prod(W,batchX2,gradient);
+ noalias(gradient) = prod(W,batchX2);
RealVector columnSum = sum_columns(coefficientsX2*s.expNorm);
for(std::size_t i = 0; i != sizeX1; ++i){
diff --git a/include/shark/Models/Kernels/KernelExpansion.h b/include/shark/Models/Kernels/KernelExpansion.h
index b5a2871..e2405db 100644
--- a/include/shark/Models/Kernels/KernelExpansion.h
+++ b/include/shark/Models/Kernels/KernelExpansion.h
@@ -252,7 +252,7 @@ public:
//get the part of the alpha matrix which is suitable for this batch
ConstRealSubMatrix batchAlpha = subrange(m_alpha,batchStart,batchEnd,0,outputSize());
- axpy_prod(trans(kernelEvaluations),batchAlpha,output,false);
+ noalias(output) += prod(trans(kernelEvaluations),batchAlpha);
batchStart = batchEnd;
}
}
diff --git a/include/shark/Models/Kernels/LinearKernel.h b/include/shark/Models/Kernels/LinearKernel.h
index 02b73bd..a5ae754 100644
--- a/include/shark/Models/Kernels/LinearKernel.h
+++ b/include/shark/Models/Kernels/LinearKernel.h
@@ -85,7 +85,7 @@ public:
void eval(ConstBatchInputReference x1, ConstBatchInputReference x2, RealMatrix& result) const{
SIZE_CHECK(x1.size2() == x2.size2());
result.resize(x1.size1(),x2.size1());
- axpy_prod(x1,trans(x2),result);
+ noalias(result) = prod(x1,trans(x2));
}
void weightedParameterDerivative(
@@ -110,7 +110,7 @@ public:
//~ SIZE_CHECK(cofficientsX2.size2() == batchX2.size1());
gradient.resize(batchX1.size1(),batchX1.size2());
- axpy_prod(coefficientsX2,batchX2,gradient);
+ noalias(gradient) = prod(coefficientsX2,batchX2);
}
virtual double featureDistanceSqr(ConstInputReference x1, ConstInputReference x2) const{
diff --git a/include/shark/Models/Kernels/MonomialKernel.h b/include/shark/Models/Kernels/MonomialKernel.h
index 2558adf..1c2f80d 100644
--- a/include/shark/Models/Kernels/MonomialKernel.h
+++ b/include/shark/Models/Kernels/MonomialKernel.h
@@ -107,8 +107,7 @@ public:
std::size_t sizeX1 = batchX1.size1();
std::size_t sizeX2 = batchX2.size1();
result.resize(sizeX1,sizeX2);
- //calculate the inner product
- axpy_prod(batchX1,trans(batchX2),result);
+ noalias(result) = prod(batchX1,trans(batchX2));
if(m_exponent != 1)
noalias(result) = pow(result,m_exponent);
}
@@ -124,7 +123,7 @@ public:
s.resize(sizeX1,sizeX2);
//calculate the inner product
- axpy_prod(batchX1,trans(batchX2),s.base);
+ noalias(s.base) = prod(batchX1,trans(batchX2));
//now do exponentiation
if(m_exponent != 1)
noalias(result) = pow(s.base,m_exponent);
@@ -174,8 +173,7 @@ public:
//The derivative of input i of batch x1 is
//g = sum_j m_exponent*weights(i,j)*x2_j
//we now sum over j which is a matrix-matrix product
- axpy_prod(weights,batchX2,gradient);
- gradient*= m_exponent;
+ noalias(gradient) = m_exponent * prod(weights,batchX2);
}
void read(InArchive& ar){
diff --git a/include/shark/Models/Kernels/PolynomialKernel.h b/include/shark/Models/Kernels/PolynomialKernel.h
index 640dfad..ba5d5d3 100644
--- a/include/shark/Models/Kernels/PolynomialKernel.h
+++ b/include/shark/Models/Kernels/PolynomialKernel.h
@@ -170,7 +170,7 @@ public:
result.resize(sizeX1,sizeX2);
//calculate the inner product
- axpy_prod(batchX1,trans(batchX2),result);
+ noalias(result) = prod(batchX1,trans(batchX2));
result += m_offset;
//now do exponentiation
if(m_degree != 1)
@@ -188,7 +188,7 @@ public:
result.resize(sizeX1,sizeX2);
//calculate the inner product
- axpy_prod(batchX1,trans(batchX2),s.base);
+ noalias(s.base) = prod(batchX1,trans(batchX2));
s.base += m_offset;
//now do exponentiation
@@ -265,7 +265,7 @@ public:
//again m_degree == 1 is easy, as it is for the i-th row
//just c_i X2;
if(m_degree == 1){
- axpy_prod(coefficientsX2,batchX2,gradient);
+ noalias(gradient) = prod(coefficientsX2,batchX2);
return;
}
@@ -275,8 +275,7 @@ public:
//and the derivative of input i of batch x1 is
//g = sum_j m_n*weights(i,j)*x2_j
//we now sum over j which is a matrix-matrix product
- axpy_prod(weights,batchX2,gradient);
- gradient*= m_degree;
+ noalias(gradient) = m_degree * prod(weights,batchX2);
}
void read(InArchive& ar){
diff --git a/include/shark/Models/LinearModel.h b/include/shark/Models/LinearModel.h
index c4438fc..5042768 100644
--- a/include/shark/Models/LinearModel.h
+++ b/include/shark/Models/LinearModel.h
@@ -177,7 +177,7 @@ public:
void eval(BatchInputType const& inputs, BatchOutputType& outputs)const{
outputs.resize(inputs.size1(),m_matrix.size1());
//we multiply with a set of row vectors from the left
- axpy_prod(inputs,trans(m_matrix),outputs);
+ noalias(gradient) = prod(inputs,trans(m_matrix));
if (hasOffset()){
noalias(outputs)+=repeat(m_offset,inputs.size1());
}
@@ -201,7 +201,7 @@ public:
blas::dense_matrix_adaptor<double> weightGradient = blas::adapt_matrix(outputs,inputs,gradient.storage());
//sum_i coefficients(output,i)*pattern(i))
- axpy_prod(trans(coefficients),patterns,weightGradient,false);
+ noalias(weightGradient) = prod(trans(coefficients),patterns);
if (hasOffset()){
std::size_t start = inputs*outputs;
@@ -219,7 +219,7 @@ public:
SIZE_CHECK(coefficients.size1() == patterns.size1());
derivative.resize(patterns.size1(),inputSize());
- axpy_prod(coefficients,m_matrix,derivative);
+ noalias(derivative) = prod(coefficients,m_matrix);
}
/// From ISerializable
diff --git a/include/shark/Models/TiedAutoencoder.h b/include/shark/Models/TiedAutoencoder.h
index 56815ab..45d3b47 100644
--- a/include/shark/Models/TiedAutoencoder.h
+++ b/include/shark/Models/TiedAutoencoder.h
@@ -171,18 +171,14 @@ public:
SIZE_CHECK(patterns.size2() == encoderMatrix().size2());
std::size_t numOutputs = encoderMatrix().size1();
outputs.resize(numPatterns,numOutputs);
- outputs.clear();
- axpy_prod(patterns,trans(encoderMatrix()),outputs);
- noalias(outputs) += repeat(hiddenBias(),numPatterns);
+ noalias(outputs) = prod(patterns,trans(encoderMatrix())) + repeat(hiddenBias(),numPatterns);
noalias(outputs) = m_hiddenNeuron(outputs);
}
else{//hidden->output
SIZE_CHECK(patterns.size2() == decoderMatrix().size2());
std::size_t numOutputs = decoderMatrix().size1();
outputs.resize(numPatterns,numOutputs);
- outputs.clear();
- axpy_prod(patterns,trans(decoderMatrix()),outputs);
- noalias(outputs) += repeat(outputBias(),numPatterns);
+ noalias(outputs) = prod(patterns,trans(decoderMatrix()),outputs) + repeat(outputBias(),numPatterns);
noalias(outputs) = m_outputNeuron(outputs);
}
}
@@ -302,7 +298,7 @@ private:
noalias(outputDelta) *= m_outputNeuron.derivative(s.outputResponses);
hiddenDelta.resize(outputDelta.size1(),numberOfHiddenNeurons());
- axpy_prod(outputDelta,decoderMatrix(),hiddenDelta,true);
+ noalias(hiddenDelta) = prod(outputDelta,decoderMatrix());
noalias(hiddenDelta) *= m_hiddenNeuron.derivative(s.hiddenResponses);
}
@@ -311,7 +307,7 @@ private:
)const{
computeDelta(state,outputDelta,hiddenDelta);
inputDelta.resize(outputDelta.size1(),inputSize());
- axpy_prod(hiddenDelta,encoderMatrix(),inputDelta,true);
+ noalias(inputDelta) = prod(hiddenDelta,encoderMatrix());
}
void computeParameterDerivative(
diff --git a/include/shark/ObjectiveFunctions/Benchmarks/RotatedErrorFunction.h b/include/shark/ObjectiveFunctions/Benchmarks/RotatedErrorFunction.h
index bf5c2d5..a92a2ae 100644
--- a/include/shark/ObjectiveFunctions/Benchmarks/RotatedErrorFunction.h
+++ b/include/shark/ObjectiveFunctions/Benchmarks/RotatedErrorFunction.h
@@ -78,8 +78,7 @@ struct RotatedObjectiveFunction : public SingleObjectiveFunction {
double eval( const SearchPointType & p ) const {
m_evaluationCounter++;
- RealVector x(numberOfVariables());
- axpy_prod(m_rotation,p,x);
+ RealVector x = prod(m_rotation,p);
return m_objective->eval(x);
}
private:
diff --git a/include/shark/ObjectiveFunctions/Impl/SparseAutoencoderError.inl b/include/shark/ObjectiveFunctions/Impl/SparseAutoencoderError.inl
index 0814973..9feed47 100644
--- a/include/shark/ObjectiveFunctions/Impl/SparseAutoencoderError.inl
+++ b/include/shark/ObjectiveFunctions/Impl/SparseAutoencoderError.inl
@@ -161,7 +161,7 @@ public:
noalias(hiddenDerivativeSum) += sum_rows(hiddenDerivative);
// Calculate the gradient with respect to the lower weight matrix
- axpy_prod(trans(hiddenDerivative),batch.input,W1Derivatives,false);
+ noalias(W1Derivatives) +=prod(trans(hiddenDerivative),batch.input);
}
error /= dataSize;
meanActivation /= dataSize;
diff --git a/include/shark/Unsupervised/RBM/Impl/AverageEnergyGradient.h b/include/shark/Unsupervised/RBM/Impl/AverageEnergyGradient.h
index 31f4901..0a9fdca 100644
--- a/include/shark/Unsupervised/RBM/Impl/AverageEnergyGradient.h
+++ b/include/shark/Unsupervised/RBM/Impl/AverageEnergyGradient.h
@@ -50,7 +50,7 @@ public:
for(std::size_t i = 0; i != batchSize; ++i){
row(weightedFeatures,i)*= weights(i);
}
- axpy_prod(trans(mpe_rbm->hiddenNeurons().expectedPhiValue(hiddens.statistics)),weightedFeatures,m_deltaWeights,false);
+ noalias(m_deltaWeights) += prod(trans(mpe_rbm->hiddenNeurons().expectedPhiValue(hiddens.statistics)),weightedFeatures);
mpe_rbm->visibleNeurons().parameterDerivative(m_deltaBiasVisible,visibles,weights);
mpe_rbm->hiddenNeurons().expectedParameterDerivative(m_deltaBiasHidden,hiddens,weights);
}
@@ -80,7 +80,7 @@ public:
row(weightedFeatures,i)*= weights(i);
}
- axpy_prod(trans(weightedFeatures),mpe_rbm->visibleNeurons().expectedPhiValue(visibles.statistics),m_deltaWeights,false);
+ noalias(m_deltaWeights) += prod(trans(weightedFeatures),mpe_rbm->visibleNeurons().expectedPhiValue(visibles.statistics));
mpe_rbm->hiddenNeurons().parameterDerivative(m_deltaBiasHidden,hiddens,weights);
mpe_rbm->visibleNeurons().expectedParameterDerivative(m_deltaBiasVisible,visibles,weights);
}
diff --git a/include/shark/Unsupervised/RBM/Neuronlayers/GaussianLayer.h b/include/shark/Unsupervised/RBM/Neuronlayers/GaussianLayer.h
index f7f1aa7..91a802f 100644
--- a/include/shark/Unsupervised/RBM/Neuronlayers/GaussianLayer.h
+++ b/include/shark/Unsupervised/RBM/Neuronlayers/GaussianLayer.h
@@ -184,8 +184,7 @@ public:
//return beta * inner_prod(m_bias,state) - norm_sqr(state)/2.0;
std::size_t batchSize = state.size1();
- RealVector energies(batchSize);
- axpy_prod(state,m_bias,energies);
+ RealVector energies = prod(state,m_bias);
noalias(energies) *= beta;
for(std::size_t i = 0; i != batchSize; ++i){
energies(i) -= norm_sqr(row(state,i))/2.0;
diff --git a/include/shark/Unsupervised/RBM/RBM.h b/include/shark/Unsupervised/RBM/RBM.h
index eefbe30..e6d0a84 100644
--- a/include/shark/Unsupervised/RBM/RBM.h
+++ b/include/shark/Unsupervised/RBM/RBM.h
@@ -241,7 +241,7 @@ public:
SIZE_CHECK(inputs.size2() == m_hiddenNeurons.size());
SIZE_CHECK( visibleStates.size2() == m_visibleNeurons.size());
- axpy_prod(m_visibleNeurons.phi(visibleStates),trans(m_weightMatrix),inputs);
+ noalias(inputs) = prod(m_visibleNeurons.phi(visibleStates),trans(m_weightMatrix));
}
@@ -253,7 +253,7 @@ public:
SIZE_CHECK(hiddenStates.size1() == inputs.size1());
SIZE_CHECK(inputs.size2() == m_visibleNeurons.size());
- axpy_prod(m_hiddenNeurons.phi(hiddenStates),m_weightMatrix,inputs);
+ noalias(inputs) = prod(m_hiddenNeurons.phi(hiddenStates),m_weightMatrix);
}
using base_type::eval;
diff --git a/src/Algorithms/GradientDescent/BFGS.cpp b/src/Algorithms/GradientDescent/BFGS.cpp
index f9ed027..56b7675 100644
--- a/src/Algorithms/GradientDescent/BFGS.cpp
+++ b/src/Algorithms/GradientDescent/BFGS.cpp
@@ -50,8 +50,7 @@ void BFGS::computeSearchDirection(){
RealVector delta = m_best.point - m_lastPoint;
double d = inner_prod(gamma,delta);
- RealVector Hg(m_dimension,0.0);
- axpy_prod(m_hessian,gamma,Hg);
+ RealVector Hg = prod(m_hessian,gamma);
//update hessian
if (d < 1e-20)
@@ -69,8 +68,7 @@ void BFGS::computeSearchDirection(){
}
//compute search direction
- axpy_prod(m_hessian,m_derivative,m_searchDirection);
- m_searchDirection *= -1;
+ noalias(m_searchDirection) = -prod(m_hessian,m_derivative);
}
//from ISerializable
diff --git a/src/Algorithms/LinearRegression.cpp b/src/Algorithms/LinearRegression.cpp
index 1df40da..715f080 100644
--- a/src/Algorithms/LinearRegression.cpp
+++ b/src/Algorithms/LinearRegression.cpp
@@ -73,7 +73,7 @@ void LinearRegression::train(LinearModel<>& model, LabeledData<RealVector, RealV
for (std::size_t b=0; b != numBatches; b++){
BatchRef batch = dataset.batch(b);
RealSubMatrix PTL = subrange(XTL,0,inputDim,0,outputDim);
- axpy_prod(trans(batch.input),batch.label,PTL,false);
+ noalias(PTL) += prod(trans(batch.input),batch.label);
noalias(row(XTL,inputDim))+=sum_rows(batch.label);
}
diff --git a/src/Algorithms/PCA.cpp b/src/Algorithms/PCA.cpp
index 81946e2..b77d9a5 100644
--- a/src/Algorithms/PCA.cpp
+++ b/src/Algorithms/PCA.cpp
@@ -38,7 +38,7 @@
using namespace shark;
- /// Set the input data, which is stored in the PCA object.
+/// Set the input data, which is stored in the PCA object.
void PCA::setData(UnlabeledData<RealVector> const& inputs) {
SHARK_CHECK(inputs.numberOfElements() >= 2, "[PCA::train] input needs to contain at least two points");
m_l = inputs.numberOfElements(); ///< number of data points
@@ -79,7 +79,7 @@ void PCA::setData(UnlabeledData<RealVector> const& inputs) {
RealMatrix X2 = inputs.batch(b2)-repeat(m_mean,batchSize2);
RealSubMatrix X1X2T= subrange(S,start1,start1+batchSize1,start2,start2+batchSize2);
RealSubMatrix X2X1T= subrange(S,start2,start2+batchSize2,start1,start1+batchSize1);
- axpy_prod(X1,trans(X2),X1X2T);// X1 X2^T
+ noalias(X1X2T) = prod(X1,trans(X2));// X1 X2^T
noalias(X2X1T) = trans(X1X2T);// X2 X1^T
start2+=batchSize2;
}
@@ -90,8 +90,6 @@ void PCA::setData(UnlabeledData<RealVector> const& inputs) {
}
S /= m_l;
m_eigenvalues.resize(m_l);
- m_eigenvectors.resize(m_n, m_l);
- m_eigenvectors.clear();
RealMatrix U(m_l, m_l);
eigensymm(S, U, m_eigenvalues);
// compute true eigenvectors
@@ -101,7 +99,7 @@ void PCA::setData(UnlabeledData<RealVector> const& inputs) {
std::size_t batchSize = inputs.batch(b).size1();
std::size_t batchEnd = batchStart+batchSize;
RealMatrix X = inputs.batch(b)-repeat(m_mean,batchSize);
- axpy_prod(trans(X),rows(U,batchStart,batchEnd),m_eigenvectors,false);
+ m_eigenvectors = prod(trans(X),rows(U,batchStart,batchEnd));
batchStart = batchEnd;
}
//normalize
diff --git a/src/Models/OnlineRNNet.cpp b/src/Models/OnlineRNNet.cpp
index 240ca4d..eb80112 100644
--- a/src/Models/OnlineRNNet.cpp
+++ b/src/Models/OnlineRNNet.cpp
@@ -112,18 +112,14 @@ void OnlineRNNet::weightedParameterDerivative(RealMatrix const& pattern, const R
);
//update the new gradient with the effect of last timestep
- RealMatrix newUnitGradient(mpe_structure->parameters(),numNeurons);
- axpy_prod(m_unitGradient,trans(hiddenWeights),newUnitGradient);
- swap(m_unitGradient,newUnitGradient);
- newUnitGradient = RealMatrix();//empty
-
+ noalias(m_unitGradient) = prod(m_unitGradient,trans(hiddenWeights));
//add the effect of the current time step
std::size_t param = 0;
for(std::size_t i = 0; i != numNeurons; ++i){
for(std::size_t j = 0; j != numUnits; ++j){
if(mpe_structure->connection(i,j)){
- m_unitGradient(param,i)+=m_lastActivation(j);
+ m_unitGradient(param,i) += m_lastActivation(j);
++param;
}
}
@@ -131,7 +127,7 @@ void OnlineRNNet::weightedParameterDerivative(RealMatrix const& pattern, const R
//multiply with outer derivative of the neurons
for(std::size_t i = 0; i != m_unitGradient.size1();++i){
- noalias(row(m_unitGradient,i))= element_prod(row(m_unitGradient,i),neuronDerivatives);
+ noalias(row(m_unitGradient,i)) = element_prod(row(m_unitGradient,i),neuronDerivatives);
}
//and formula 4 (the gradient itself)
noalias(gradient) = prod(
diff --git a/src/Models/RBFLayer.cpp b/src/Models/RBFLayer.cpp
index 69ec87d..3ced18f 100644
--- a/src/Models/RBFLayer.cpp
+++ b/src/Models/RBFLayer.cpp
@@ -142,8 +142,7 @@ void RBFLayer::weightedParameterDerivative(
//the second part is than just a matrix-diagonal multiplication
blas::dense_matrix_adaptor<double> centerDerivative = blas::adapt_matrix(outputSize(),inputSize(),&gradient(currentParameter));
- //compute first part
- axpy_prod(trans(delta),patterns,centerDerivative);
+ noalias(centerDerivative) = prod(trans(delta),patterns);
//compute second part
for(std::size_t i = 0; i != outputSize(); ++i){
noalias(row(centerDerivative,i)) -= deltaSum(i)*row(m_centers,i);
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/shark.git
More information about the debian-science-commits
mailing list