// Copyright 2015 Conrad Sanderson (http://conradsanderson.id.au) // Copyright 2015 National ICT Australia (NICTA) // // Licensed under the Apache License, Version 2.0 (the "License"); // you may not use this file except in compliance with the License. // You may obtain a copy of the License at // http://www.apache.org/licenses/LICENSE-2.0 // // Unless required by applicable law or agreed to in writing, software // distributed under the License is distributed on an "AS IS" BASIS, // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. // See the License for the specific language governing permissions and // limitations under the License. // ------------------------------------------------------------------------ #include #include "catch.hpp" using namespace arma; TEST_CASE("decomp_eig_pair_1") { mat A1 = "\ 0.840375529753905 -0.600326562133734 -2.138355269439939 0.124049800003193 2.908008030729362;\ -0.888032082329010 0.489965321173948 -0.839588747336614 1.436696622718939 0.825218894228491;\ 0.100092833139322 0.739363123604474 1.354594328004644 -1.960899999365033 1.378971977916614;\ -0.544528929990548 1.711887782981555 -1.072155288384252 -0.197698225974150 -1.058180257987362;\ 0.303520794649354 -0.194123535758265 0.960953869740567 -1.207845485259799 -0.468615581100624;\ "; mat A2 = "\ -0.272469409250187 -0.353849997774433 0.033479882244451 0.022889792751630 -0.979206305167302;\ 1.098424617888623 -0.823586525156853 -1.333677943428106 -0.261995434966092 -1.156401655664002;\ -0.277871932787639 -1.577057022799202 1.127492278341590 -1.750212368446790 -0.533557109315987;\ 0.701541458163284 0.507974650905946 0.350179410603312 -0.285650971595330 -2.002635735883060;\ -2.051816299911149 0.281984063670556 -0.299066030332982 -0.831366511567624 0.964229422631627;\ "; cx_vec eigvals1 = { cx_double(-2.467066249890401, +0.000000000000000), cx_double( 1.483137782196390, +0.595028644066690), cx_double( 1.483137782196390, -0.595028644066690), cx_double(-0.646831879916377, +0.000000000000000), cx_double( 0.099992295916005, +0.000000000000000) }; cx_vec eigvals2 = eig_pair(A1, A2); cx_vec eigvals3; bool status = eig_pair(eigvals3, A1, A2); cx_vec eigvals4; cx_mat eigvecs4; eig_pair(eigvals4, eigvecs4, A1, A2); cx_vec eigvals5; cx_mat leigvecs5; cx_mat reigvecs5; eig_pair(eigvals5, leigvecs5, reigvecs5, A1, A2); cx_mat B = A2 * eigvecs4 * diagmat(eigvals4) * inv( eigvecs4) ; cx_mat Cl = inv(trans(leigvecs5)) * diagmat(eigvals5) * trans(leigvecs5) * A2; cx_mat Cr = A2 * reigvecs5 * diagmat(eigvals5) * inv(reigvecs5) ; REQUIRE( status == true ); REQUIRE( accu(abs(eigvals2 - eigvals1)) == Approx(0.0).epsilon(0.0001) ); REQUIRE( accu(abs(eigvals3 - eigvals1)) == Approx(0.0).epsilon(0.0001) ); REQUIRE( accu(abs(eigvals4 - eigvals1)) == Approx(0.0).epsilon(0.0001) ); REQUIRE( accu(abs(eigvals5 - eigvals1)) == Approx(0.0).epsilon(0.0001) ); REQUIRE( accu(abs(A1 - B )) == Approx(0.0).epsilon(0.0001) ); REQUIRE( accu(abs(A1 - Cl )) == Approx(0.0).epsilon(0.0001) ); REQUIRE( accu(abs(A1 - Cr )) == Approx(0.0).epsilon(0.0001) ); } TEST_CASE("decomp_eig_pair_2") { cx_mat A1 = { { cx_double( 0.520060101455458, +0.451679418928238), cx_double(-0.133217479507735, -1.361694470870754), cx_double(-0.293753597735416, +1.039090653504956), cx_double( 0.307535159238252, -0.195221197898754), cx_double(-1.332004421315247, +0.826062790211595) }, { cx_double(-0.020027851642538, -0.130284653145721), cx_double(-0.714530163787158, +0.455029556444334), cx_double(-0.847926243637934, -1.117638683265208), cx_double(-1.257118359352053, -0.217606350143192), cx_double(-2.329867155805076, +1.526976686733373) }, { cx_double(-0.034771086028483, +0.183689095861942), cx_double( 1.351385768426657, -0.848709379933659), cx_double(-1.120128301243728, +1.260658709120896), cx_double(-0.865468030554804, -0.303107621351741), cx_double(-1.449097292838739, +0.466914435684700) }, { cx_double(-0.798163584564142, -0.476153016619074), cx_double(-0.224771056052584, -0.334886938964048), cx_double( 2.525999692118309, +0.660143141046978), cx_double(-0.176534114231451, +0.023045624425105), cx_double( 0.333510833065806, -0.209713338388737) }, { cx_double( 1.018685282128575, +0.862021611556922), cx_double(-0.589029030720801, +0.552783345944550), cx_double( 1.655497592887346, -0.067865553542687), cx_double( 0.791416061628634, +0.051290355848775), cx_double( 0.391353604432901, +0.625190357087626) } }; cx_mat A2 = { { cx_double( 0.183227263001437, -0.444627816446985), cx_double( 0.515246335524849, +0.391894209432449), cx_double(-0.532011376808821, -0.320575506600239), cx_double(-1.174212331456816, -1.066701398984750), cx_double(-1.064213412889327, -1.565056014150725) }, { cx_double(-1.029767543566621, -0.155941035724769), cx_double( 0.261406324055383, -1.250678906826407), cx_double( 1.682103594663179, +0.012469041361618), cx_double(-0.192239517539275, +0.933728162671238), cx_double( 1.603457298120044, -0.084539479817724) }, { cx_double( 0.949221831131023, +0.276068253931536), cx_double(-0.941485770955434, -0.947960922331432), cx_double(-0.875729346160017, -3.029177341404146), cx_double(-0.274070229932602, +0.350321001356112), cx_double( 1.234679146890778, +1.603946350602880) }, { cx_double( 0.307061919146703, -0.261163645776479), cx_double(-0.162337672803828, -0.741106093940411), cx_double(-0.483815050110121, -0.457014640871583), cx_double( 1.530072514424096, -0.029005763708726), cx_double(-0.229626450963180, +0.098347774640108) }, { cx_double( 0.135174942099456, +0.443421912904091), cx_double(-0.146054634331526, -0.507817550278174), cx_double(-0.712004549027422, +1.242448406390738), cx_double(-0.249024742513714, +0.182452167505983), cx_double(-1.506159703979719, +0.041373613489615) } }; cx_vec eigvals1 = { cx_double(-0.567948485992280, -1.314594536444777), cx_double( 1.051873153748018, -0.162676262480913), cx_double( 0.610087089288344, +0.562148335263468), cx_double(-0.890023643973463, -0.352930772605452), cx_double(-0.365476750045154, +0.305826583179225) }; cx_vec eigvals2 = eig_pair(A1, A2); cx_vec eigvals3; bool status = eig_pair(eigvals3, A1, A2); cx_vec eigvals4; cx_mat eigvecs4; eig_pair(eigvals4, eigvecs4, A1, A2); cx_vec eigvals5; cx_mat leigvecs5; cx_mat reigvecs5; eig_pair(eigvals5, leigvecs5, reigvecs5, A1, A2); cx_mat B = A2 * eigvecs4 * diagmat(eigvals4) * inv( eigvecs4) ; cx_mat Cl = inv(trans(leigvecs5)) * diagmat(eigvals5) * trans(leigvecs5) * A2; cx_mat Cr = A2 * reigvecs5 * diagmat(eigvals5) * inv(reigvecs5) ; REQUIRE( status == true ); REQUIRE( accu(abs(eigvals2 - eigvals1)) == Approx(0.0).epsilon(0.0001) ); REQUIRE( accu(abs(eigvals3 - eigvals1)) == Approx(0.0).epsilon(0.0001) ); REQUIRE( accu(abs(eigvals4 - eigvals1)) == Approx(0.0).epsilon(0.0001) ); REQUIRE( accu(abs(eigvals5 - eigvals1)) == Approx(0.0).epsilon(0.0001) ); REQUIRE( accu(abs(A1 - B )) == Approx(0.0).epsilon(0.0001) ); REQUIRE( accu(abs(A1 - Cl )) == Approx(0.0).epsilon(0.0001) ); REQUIRE( accu(abs(A1 - Cr )) == Approx(0.0).epsilon(0.0001) ); }