kronecker_product.cpp
Go to the documentation of this file.
00001 // This file is part of Eigen, a lightweight C++ template library
00002 // for linear algebra.
00003 //
00004 // Copyright (C) 2011 Kolja Brix <brix@igpm.rwth-aachen.de>
00005 // Copyright (C) 2011 Andreas Platen <andiplaten@gmx.de>
00006 // Copyright (C) 2012 Chen-Pang He <jdh8@ms63.hinet.net>
00007 //
00008 // This Source Code Form is subject to the terms of the Mozilla
00009 // Public License v. 2.0. If a copy of the MPL was not distributed
00010 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
00011 
00012 
00013 #include "sparse.h"
00014 #include <Eigen/SparseExtra>
00015 #include <Eigen/KroneckerProduct>
00016 
00017 
00018 template<typename MatrixType>
00019 void check_dimension(const MatrixType& ab, const int rows,  const int cols)
00020 {
00021   VERIFY_IS_EQUAL(ab.rows(), rows);
00022   VERIFY_IS_EQUAL(ab.cols(), cols);
00023 }
00024 
00025 
00026 template<typename MatrixType>
00027 void check_kronecker_product(const MatrixType& ab)
00028 {
00029   VERIFY_IS_EQUAL(ab.rows(), 6);
00030   VERIFY_IS_EQUAL(ab.cols(), 6);
00031   VERIFY_IS_EQUAL(ab.nonZeros(),  36);
00032   VERIFY_IS_APPROX(ab.coeff(0,0), -0.4017367630386106);
00033   VERIFY_IS_APPROX(ab.coeff(0,1),  0.1056863433932735);
00034   VERIFY_IS_APPROX(ab.coeff(0,2), -0.7255206194554212);
00035   VERIFY_IS_APPROX(ab.coeff(0,3),  0.1908653336744706);
00036   VERIFY_IS_APPROX(ab.coeff(0,4),  0.350864567234111);
00037   VERIFY_IS_APPROX(ab.coeff(0,5), -0.0923032108308013);
00038   VERIFY_IS_APPROX(ab.coeff(1,0),  0.415417514804677);
00039   VERIFY_IS_APPROX(ab.coeff(1,1), -0.2369227701722048);
00040   VERIFY_IS_APPROX(ab.coeff(1,2),  0.7502275131458511);
00041   VERIFY_IS_APPROX(ab.coeff(1,3), -0.4278731019742696);
00042   VERIFY_IS_APPROX(ab.coeff(1,4), -0.3628129162264507);
00043   VERIFY_IS_APPROX(ab.coeff(1,5),  0.2069210808481275);
00044   VERIFY_IS_APPROX(ab.coeff(2,0),  0.05465890160863986);
00045   VERIFY_IS_APPROX(ab.coeff(2,1), -0.2634092511419858);
00046   VERIFY_IS_APPROX(ab.coeff(2,2),  0.09871180285793758);
00047   VERIFY_IS_APPROX(ab.coeff(2,3), -0.4757066334017702);
00048   VERIFY_IS_APPROX(ab.coeff(2,4), -0.04773740823058334);
00049   VERIFY_IS_APPROX(ab.coeff(2,5),  0.2300535609645254);
00050   VERIFY_IS_APPROX(ab.coeff(3,0), -0.8172945853260133);
00051   VERIFY_IS_APPROX(ab.coeff(3,1),  0.2150086428359221);
00052   VERIFY_IS_APPROX(ab.coeff(3,2),  0.5825113847292743);
00053   VERIFY_IS_APPROX(ab.coeff(3,3), -0.1532433770097174);
00054   VERIFY_IS_APPROX(ab.coeff(3,4), -0.329383387282399);
00055   VERIFY_IS_APPROX(ab.coeff(3,5),  0.08665207912033064);
00056   VERIFY_IS_APPROX(ab.coeff(4,0),  0.8451267514863225);
00057   VERIFY_IS_APPROX(ab.coeff(4,1), -0.481996458918977);
00058   VERIFY_IS_APPROX(ab.coeff(4,2), -0.6023482390791535);
00059   VERIFY_IS_APPROX(ab.coeff(4,3),  0.3435339347164565);
00060   VERIFY_IS_APPROX(ab.coeff(4,4),  0.3406002157428891);
00061   VERIFY_IS_APPROX(ab.coeff(4,5), -0.1942526344200915);
00062   VERIFY_IS_APPROX(ab.coeff(5,0),  0.1111982482925399);
00063   VERIFY_IS_APPROX(ab.coeff(5,1), -0.5358806424754169);
00064   VERIFY_IS_APPROX(ab.coeff(5,2), -0.07925446559335647);
00065   VERIFY_IS_APPROX(ab.coeff(5,3),  0.3819388757769038);
00066   VERIFY_IS_APPROX(ab.coeff(5,4),  0.04481475387219876);
00067   VERIFY_IS_APPROX(ab.coeff(5,5), -0.2159688616158057);
00068 }
00069 
00070 
00071 template<typename MatrixType>
00072 void check_sparse_kronecker_product(const MatrixType& ab)
00073 {
00074   VERIFY_IS_EQUAL(ab.rows(), 12);
00075   VERIFY_IS_EQUAL(ab.cols(), 10);
00076   VERIFY_IS_EQUAL(ab.nonZeros(), 3*2);
00077   VERIFY_IS_APPROX(ab.coeff(3,0), -0.04);
00078   VERIFY_IS_APPROX(ab.coeff(5,1),  0.05);
00079   VERIFY_IS_APPROX(ab.coeff(0,6), -0.08);
00080   VERIFY_IS_APPROX(ab.coeff(2,7),  0.10);
00081   VERIFY_IS_APPROX(ab.coeff(6,8),  0.12);
00082   VERIFY_IS_APPROX(ab.coeff(8,9), -0.15);
00083 }
00084 
00085 
00086 void test_kronecker_product()
00087 {
00088   // DM = dense matrix; SM = sparse matrix
00089 
00090   Matrix<double, 2, 3> DM_a;
00091   SparseMatrix<double> SM_a(2,3);
00092   SM_a.insert(0,0) = DM_a.coeffRef(0,0) = -0.4461540300782201;
00093   SM_a.insert(0,1) = DM_a.coeffRef(0,1) = -0.8057364375283049;
00094   SM_a.insert(0,2) = DM_a.coeffRef(0,2) =  0.3896572459516341;
00095   SM_a.insert(1,0) = DM_a.coeffRef(1,0) = -0.9076572187376921;
00096   SM_a.insert(1,1) = DM_a.coeffRef(1,1) =  0.6469156566545853;
00097   SM_a.insert(1,2) = DM_a.coeffRef(1,2) = -0.3658010398782789;
00098  
00099   MatrixXd             DM_b(3,2);
00100   SparseMatrix<double> SM_b(3,2);
00101   SM_b.insert(0,0) = DM_b.coeffRef(0,0) =  0.9004440976767099;
00102   SM_b.insert(0,1) = DM_b.coeffRef(0,1) = -0.2368830858139832;
00103   SM_b.insert(1,0) = DM_b.coeffRef(1,0) = -0.9311078389941825;
00104   SM_b.insert(1,1) = DM_b.coeffRef(1,1) =  0.5310335762980047;
00105   SM_b.insert(2,0) = DM_b.coeffRef(2,0) = -0.1225112806872035;
00106   SM_b.insert(2,1) = DM_b.coeffRef(2,1) =  0.5903998022741264;
00107 
00108   SparseMatrix<double,RowMajor> SM_row_a(SM_a), SM_row_b(SM_b);
00109 
00110   // test kroneckerProduct(DM_block,DM,DM_fixedSize)
00111   Matrix<double, 6, 6> DM_fix_ab = kroneckerProduct(DM_a.topLeftCorner<2,3>(),DM_b);
00112 
00113   CALL_SUBTEST(check_kronecker_product(DM_fix_ab));
00114 
00115   for(int i=0;i<DM_fix_ab.rows();++i)
00116     for(int j=0;j<DM_fix_ab.cols();++j)
00117        VERIFY_IS_APPROX(kroneckerProduct(DM_a,DM_b).coeff(i,j), DM_fix_ab(i,j));
00118 
00119   // test kroneckerProduct(DM,DM,DM_block)
00120   MatrixXd DM_block_ab(10,15);
00121   DM_block_ab.block<6,6>(2,5) = kroneckerProduct(DM_a,DM_b);
00122   CALL_SUBTEST(check_kronecker_product(DM_block_ab.block<6,6>(2,5)));
00123 
00124   // test kroneckerProduct(DM,DM,DM)
00125   MatrixXd DM_ab = kroneckerProduct(DM_a,DM_b);
00126   CALL_SUBTEST(check_kronecker_product(DM_ab));
00127 
00128   // test kroneckerProduct(SM,DM,SM)
00129   SparseMatrix<double> SM_ab = kroneckerProduct(SM_a,DM_b);
00130   CALL_SUBTEST(check_kronecker_product(SM_ab));
00131   SparseMatrix<double,RowMajor> SM_ab2 = kroneckerProduct(SM_a,DM_b);
00132   CALL_SUBTEST(check_kronecker_product(SM_ab2));
00133 
00134   // test kroneckerProduct(DM,SM,SM)
00135   SM_ab.setZero();
00136   SM_ab.insert(0,0)=37.0;
00137   SM_ab = kroneckerProduct(DM_a,SM_b);
00138   CALL_SUBTEST(check_kronecker_product(SM_ab));
00139   SM_ab2.setZero();
00140   SM_ab2.insert(0,0)=37.0;
00141   SM_ab2 = kroneckerProduct(DM_a,SM_b);
00142   CALL_SUBTEST(check_kronecker_product(SM_ab2));
00143 
00144   // test kroneckerProduct(SM,SM,SM)
00145   SM_ab.resize(2,33);
00146   SM_ab.insert(0,0)=37.0;
00147   SM_ab = kroneckerProduct(SM_a,SM_b);
00148   CALL_SUBTEST(check_kronecker_product(SM_ab));
00149   SM_ab2.resize(5,11);
00150   SM_ab2.insert(0,0)=37.0;
00151   SM_ab2 = kroneckerProduct(SM_a,SM_b);
00152   CALL_SUBTEST(check_kronecker_product(SM_ab2));
00153 
00154   // test kroneckerProduct(SM,SM,SM) with sparse pattern
00155   SM_a.resize(4,5);
00156   SM_b.resize(3,2);
00157   SM_a.resizeNonZeros(0);
00158   SM_b.resizeNonZeros(0);
00159   SM_a.insert(1,0) = -0.1;
00160   SM_a.insert(0,3) = -0.2;
00161   SM_a.insert(2,4) =  0.3;
00162   SM_a.finalize();
00163   
00164   SM_b.insert(0,0) =  0.4;
00165   SM_b.insert(2,1) = -0.5;
00166   SM_b.finalize();
00167   SM_ab.resize(1,1);
00168   SM_ab.insert(0,0)=37.0;
00169   SM_ab = kroneckerProduct(SM_a,SM_b);
00170   CALL_SUBTEST(check_sparse_kronecker_product(SM_ab));
00171 
00172   // test dimension of result of kroneckerProduct(DM,DM,DM)
00173   MatrixXd DM_a2(2,1);
00174   MatrixXd DM_b2(5,4);
00175   MatrixXd DM_ab2 = kroneckerProduct(DM_a2,DM_b2);
00176   CALL_SUBTEST(check_dimension(DM_ab2,2*5,1*4));
00177   DM_a2.resize(10,9);
00178   DM_b2.resize(4,8);
00179   DM_ab2 = kroneckerProduct(DM_a2,DM_b2);
00180   CALL_SUBTEST(check_dimension(DM_ab2,10*4,9*8));
00181 }

shape_reconstruction
Author(s): Roberto Martín-Martín
autogenerated on Sat Jun 8 2019 18:32:35