From 93ac51b73517c81bb563b0ef3effb73a04457e4f Mon Sep 17 00:00:00 2001 From: Kevin Date: Wed, 24 Jul 2013 23:17:42 -0400 Subject: [PATCH 1/3] changed ND Array input to const for the AJD class --- .../mathematics/ajd/ApproxJointDiagonalizer.h | 2 +- src/shogun/mathematics/ajd/FFDiag.cpp | 2 +- src/shogun/mathematics/ajd/FFDiag.h | 4 +- src/shogun/mathematics/ajd/JADiag.cpp | 2 +- src/shogun/mathematics/ajd/JADiag.h | 4 +- src/shogun/mathematics/ajd/JADiagOrth.cpp | 193 +++++++++--------- src/shogun/mathematics/ajd/JADiagOrth.h | 8 +- src/shogun/mathematics/ajd/UWedge.cpp | 2 +- src/shogun/mathematics/ajd/UWedge.h | 4 +- 9 files changed, 114 insertions(+), 107 deletions(-) diff --git a/src/shogun/mathematics/ajd/ApproxJointDiagonalizer.h b/src/shogun/mathematics/ajd/ApproxJointDiagonalizer.h index 689ccad61aa..2f6fc115c09 100644 --- a/src/shogun/mathematics/ajd/ApproxJointDiagonalizer.h +++ b/src/shogun/mathematics/ajd/ApproxJointDiagonalizer.h @@ -54,7 +54,7 @@ class CApproxJointDiagonalizer : public CSGObject * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - virtual SGMatrix compute(SGNDArray &C, + virtual SGMatrix compute(const SGNDArray &C, SGMatrix V0= SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200) = 0; diff --git a/src/shogun/mathematics/ajd/FFDiag.cpp b/src/shogun/mathematics/ajd/FFDiag.cpp index 4b0a74fc221..1adc82702a4 100644 --- a/src/shogun/mathematics/ajd/FFDiag.cpp +++ b/src/shogun/mathematics/ajd/FFDiag.cpp @@ -16,7 +16,7 @@ using namespace shogun; void getW(double *C, int *ptN, int *ptK, double *W); -SGMatrix CFFDiag::diagonalize(SGNDArray &C0, SGMatrix V0, +SGMatrix CFFDiag::diagonalize(const SGNDArray &C0, SGMatrix V0, double eps, int itermax) { int n = C0.dims[0]; diff --git a/src/shogun/mathematics/ajd/FFDiag.h b/src/shogun/mathematics/ajd/FFDiag.h index 39b1042aab4..ddb1e392e94 100644 --- a/src/shogun/mathematics/ajd/FFDiag.h +++ b/src/shogun/mathematics/ajd/FFDiag.h @@ -53,7 +53,7 @@ class CFFDiag : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - static SGMatrix diagonalize(SGNDArray &C, + static SGMatrix diagonalize(const SGNDArray &C, SGMatrix V0= SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200); @@ -65,7 +65,7 @@ class CFFDiag : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - virtual SGMatrix compute(SGNDArray &C, + virtual SGMatrix compute(const SGNDArray &C, SGMatrix V0= SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200) diff --git a/src/shogun/mathematics/ajd/JADiag.cpp b/src/shogun/mathematics/ajd/JADiag.cpp index a6641ebdc84..47268b3af05 100644 --- a/src/shogun/mathematics/ajd/JADiag.cpp +++ b/src/shogun/mathematics/ajd/JADiag.cpp @@ -17,7 +17,7 @@ using namespace shogun; void jadiagw(double c[], double w[], int *ptn, int *ptm, double a[], double *logdet, double *decr, double *result); -SGMatrix CJADiag::diagonalize(SGNDArray &C, SGMatrix V0, +SGMatrix CJADiag::diagonalize(const SGNDArray &C, SGMatrix V0, double eps, int itermax) { int d = C.dims[0]; diff --git a/src/shogun/mathematics/ajd/JADiag.h b/src/shogun/mathematics/ajd/JADiag.h index 5783ca9c02d..73dc1159850 100644 --- a/src/shogun/mathematics/ajd/JADiag.h +++ b/src/shogun/mathematics/ajd/JADiag.h @@ -53,7 +53,7 @@ class CJADiag : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - static SGMatrix diagonalize(SGNDArray &C, + static SGMatrix diagonalize(const SGNDArray &C, SGMatrix V0= SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200); @@ -65,7 +65,7 @@ class CJADiag : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - virtual SGMatrix compute(SGNDArray &C, + virtual SGMatrix compute(const SGNDArray &C, SGMatrix V0= SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200) diff --git a/src/shogun/mathematics/ajd/JADiagOrth.cpp b/src/shogun/mathematics/ajd/JADiagOrth.cpp index aed3c8ca8a8..565d27e7556 100644 --- a/src/shogun/mathematics/ajd/JADiagOrth.cpp +++ b/src/shogun/mathematics/ajd/JADiagOrth.cpp @@ -19,11 +19,18 @@ void left_rot_stack(double *A, int M, int N, int K, int p, int q, double c, doub void right_rot_stack(double *A, int M, int N, int K, int p, int q, double c, double s); void left_rot_simple(double *A, int m, int n, int p, int q, double c, double s); -SGMatrix CJADiagOrth::diagonalize(SGNDArray &C, SGMatrix V0, +SGMatrix CJADiagOrth::diagonalize(const SGNDArray &C0, SGMatrix V0, double eps, int itermax) -{ - int m = C.dims[0]; - int L = C.dims[2]; +{ + int m = C0.dims[0]; + int L = C0.dims[2]; + + index_t * C_dims = SG_MALLOC(index_t, 3); + C_dims[0] = C0.dims[0]; + C_dims[1] = C0.dims[1]; + C_dims[2] = C0.dims[2]; + SGNDArray C(C_dims,3); + memcpy(C.array, C0.array, C0.dims[0]*C0.dims[1]*C0.dims[2]*sizeof(float64_t)); SGMatrix V; if (V0.num_rows != 0) @@ -40,28 +47,28 @@ SGMatrix CJADiagOrth::diagonalize(SGNDArray &C, SGMatrix eps) - { - double c = cos(theta); - double s = sin(theta); - left_rot_stack (C.array, m, m, L, p, q, c, s); - right_rot_stack(C.array, m, m, L, p, q, c, s); - left_rot_simple(V.matrix, m, m, p, q, c, s); - rots++; - more = true; - } - } - } + more = false; + + for (int p = 0; p < m; p++) + { + for (int q = p+1; q < m; q++) + { + // computation of Givens angle + double theta = givens_stack(C.array, m, L, p, q); + + // Givens update + if (fabs(theta) > eps) + { + double c = cos(theta); + double s = sin(theta); + left_rot_stack (C.array, m, m, L, p, q, c, s); + right_rot_stack(C.array, m, m, L, p, q, c, s); + left_rot_simple(V.matrix, m, m, p, q, c, s); + rots++; + more = true; + } + } + } } return V; @@ -70,33 +77,33 @@ SGMatrix CJADiagOrth::diagonalize(SGNDArray &C, SGMatrix Ak(mxn) x R where R rotates the (p,q) columns R =[ c s ; -s c ] and Ak is the k-th M*N matrix in the stack */ void right_rot_stack(double *A, int M, int N, int K, int p, int q, double c, double s ) { - int k, ix, iy, cpt, kMN; - int pM = p*M; - int qM = q*M; - double nx, ny; - - for (k=0, kMN=0; k diagonalize(SGNDArray &C, - SGMatrix V0= SGMatrix(NULL,0,0), + static SGMatrix diagonalize(const SGNDArray &C, + SGMatrix V0 = SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200); @@ -65,8 +65,8 @@ class CJADiagOrth : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - virtual SGMatrix compute(SGNDArray &C, - SGMatrix V0= SGMatrix(NULL,0,0), + virtual SGMatrix compute(const SGNDArray &C, + SGMatrix V0 = SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200) { diff --git a/src/shogun/mathematics/ajd/UWedge.cpp b/src/shogun/mathematics/ajd/UWedge.cpp index e3702382166..ee2f69e6c70 100644 --- a/src/shogun/mathematics/ajd/UWedge.cpp +++ b/src/shogun/mathematics/ajd/UWedge.cpp @@ -14,7 +14,7 @@ typedef Matrix< float64_t, Dynamic, Dynamic, ColMajor > EMatrix; using namespace shogun; -SGMatrix CUWedge::diagonalize(SGNDArray &C, SGMatrix V0, +SGMatrix CUWedge::diagonalize(const SGNDArray &C, SGMatrix V0, double eps, int itermax) { int d = C.dims[0]; diff --git a/src/shogun/mathematics/ajd/UWedge.h b/src/shogun/mathematics/ajd/UWedge.h index a46edc44a04..5f4fc9799c2 100644 --- a/src/shogun/mathematics/ajd/UWedge.h +++ b/src/shogun/mathematics/ajd/UWedge.h @@ -52,7 +52,7 @@ class CUWedge : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - static SGMatrix diagonalize(SGNDArray &C, + static SGMatrix diagonalize(const SGNDArray &C, SGMatrix V0= SGMatrix(NULL,0,0), double eps=1e-12, int itermax=200); @@ -64,7 +64,7 @@ class CUWedge : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - virtual SGMatrix compute(SGNDArray &C, + virtual SGMatrix compute(const SGNDArray &C, SGMatrix V0= SGMatrix(NULL,0,0), double eps=1e-12, int itermax=200) From d2dea4311592a27919437ae0750393b32d32bf71 Mon Sep 17 00:00:00 2001 From: Kevin Date: Wed, 24 Jul 2013 23:19:02 -0400 Subject: [PATCH 2/3] fixed a memory leak in the ICA unit tests and example --- examples/undocumented/libshogun/converter_sobi_bss.cpp | 2 ++ tests/unit/converter/ica/Jade_unittest.cc | 1 + tests/unit/converter/ica/SOBI_unittest.cc | 1 + 3 files changed, 4 insertions(+) diff --git a/examples/undocumented/libshogun/converter_sobi_bss.cpp b/examples/undocumented/libshogun/converter_sobi_bss.cpp index 14f663d3f4a..fc021817dc7 100644 --- a/examples/undocumented/libshogun/converter_sobi_bss.cpp +++ b/examples/undocumented/libshogun/converter_sobi_bss.cpp @@ -91,7 +91,9 @@ void test() std::cout << "Separation error: " << sep_error << std::endl; SG_UNREF(sobi); + SG_UNREF(mixed_signals); SG_UNREF(signals); + return; } diff --git a/tests/unit/converter/ica/Jade_unittest.cc b/tests/unit/converter/ica/Jade_unittest.cc index 76269a8e1b6..05009ce46df 100644 --- a/tests/unit/converter/ica/Jade_unittest.cc +++ b/tests/unit/converter/ica/Jade_unittest.cc @@ -61,6 +61,7 @@ TEST(CJade, blind_source_separation) EXPECT_EQ(isperm,true); SG_UNREF(jade); + SG_UNREF(mixed_signals); SG_UNREF(signals); } diff --git a/tests/unit/converter/ica/SOBI_unittest.cc b/tests/unit/converter/ica/SOBI_unittest.cc index 09c75f3b301..8e0939d3879 100644 --- a/tests/unit/converter/ica/SOBI_unittest.cc +++ b/tests/unit/converter/ica/SOBI_unittest.cc @@ -56,6 +56,7 @@ TEST(CSOBI, blind_source_separation) EXPECT_LE(sep_error, 1e-6); SG_UNREF(sobi); + SG_UNREF(mixed_signals); SG_UNREF(signals); } From 27a4e9ec631dac742e66b831d452dcc143ad033e Mon Sep 17 00:00:00 2001 From: Kevin Date: Thu, 25 Jul 2013 11:32:07 -0400 Subject: [PATCH 3/3] changed AJD interface to pass by value --- .../mathematics/ajd/ApproxJointDiagonalizer.h | 4 ++-- src/shogun/mathematics/ajd/FFDiag.cpp | 2 +- src/shogun/mathematics/ajd/FFDiag.h | 8 ++++---- src/shogun/mathematics/ajd/JADiag.cpp | 2 +- src/shogun/mathematics/ajd/JADiag.h | 8 ++++---- src/shogun/mathematics/ajd/JADiagOrth.cpp | 13 +++---------- src/shogun/mathematics/ajd/JADiagOrth.h | 4 ++-- src/shogun/mathematics/ajd/UWedge.cpp | 2 +- src/shogun/mathematics/ajd/UWedge.h | 8 ++++---- 9 files changed, 22 insertions(+), 29 deletions(-) diff --git a/src/shogun/mathematics/ajd/ApproxJointDiagonalizer.h b/src/shogun/mathematics/ajd/ApproxJointDiagonalizer.h index 2f6fc115c09..3d84d71c3a9 100644 --- a/src/shogun/mathematics/ajd/ApproxJointDiagonalizer.h +++ b/src/shogun/mathematics/ajd/ApproxJointDiagonalizer.h @@ -54,8 +54,8 @@ class CApproxJointDiagonalizer : public CSGObject * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - virtual SGMatrix compute(const SGNDArray &C, - SGMatrix V0= SGMatrix(NULL,0,0), + virtual SGMatrix compute(SGNDArray C, + SGMatrix V0 = SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200) = 0; diff --git a/src/shogun/mathematics/ajd/FFDiag.cpp b/src/shogun/mathematics/ajd/FFDiag.cpp index 1adc82702a4..365ede1d0a4 100644 --- a/src/shogun/mathematics/ajd/FFDiag.cpp +++ b/src/shogun/mathematics/ajd/FFDiag.cpp @@ -16,7 +16,7 @@ using namespace shogun; void getW(double *C, int *ptN, int *ptK, double *W); -SGMatrix CFFDiag::diagonalize(const SGNDArray &C0, SGMatrix V0, +SGMatrix CFFDiag::diagonalize(SGNDArray C0, SGMatrix V0, double eps, int itermax) { int n = C0.dims[0]; diff --git a/src/shogun/mathematics/ajd/FFDiag.h b/src/shogun/mathematics/ajd/FFDiag.h index ddb1e392e94..85d27679b79 100644 --- a/src/shogun/mathematics/ajd/FFDiag.h +++ b/src/shogun/mathematics/ajd/FFDiag.h @@ -53,8 +53,8 @@ class CFFDiag : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - static SGMatrix diagonalize(const SGNDArray &C, - SGMatrix V0= SGMatrix(NULL,0,0), + static SGMatrix diagonalize(SGNDArray C, + SGMatrix V0 = SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200); @@ -65,8 +65,8 @@ class CFFDiag : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - virtual SGMatrix compute(const SGNDArray &C, - SGMatrix V0= SGMatrix(NULL,0,0), + virtual SGMatrix compute(SGNDArray C, + SGMatrix V0 = SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200) { diff --git a/src/shogun/mathematics/ajd/JADiag.cpp b/src/shogun/mathematics/ajd/JADiag.cpp index 47268b3af05..256cb975ddd 100644 --- a/src/shogun/mathematics/ajd/JADiag.cpp +++ b/src/shogun/mathematics/ajd/JADiag.cpp @@ -17,7 +17,7 @@ using namespace shogun; void jadiagw(double c[], double w[], int *ptn, int *ptm, double a[], double *logdet, double *decr, double *result); -SGMatrix CJADiag::diagonalize(const SGNDArray &C, SGMatrix V0, +SGMatrix CJADiag::diagonalize(SGNDArray C, SGMatrix V0, double eps, int itermax) { int d = C.dims[0]; diff --git a/src/shogun/mathematics/ajd/JADiag.h b/src/shogun/mathematics/ajd/JADiag.h index 73dc1159850..db36337cfda 100644 --- a/src/shogun/mathematics/ajd/JADiag.h +++ b/src/shogun/mathematics/ajd/JADiag.h @@ -53,8 +53,8 @@ class CJADiag : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - static SGMatrix diagonalize(const SGNDArray &C, - SGMatrix V0= SGMatrix(NULL,0,0), + static SGMatrix diagonalize(SGNDArray C, + SGMatrix V0 = SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200); @@ -65,8 +65,8 @@ class CJADiag : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - virtual SGMatrix compute(const SGNDArray &C, - SGMatrix V0= SGMatrix(NULL,0,0), + virtual SGMatrix compute(SGNDArray C, + SGMatrix V0 = SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200) { diff --git a/src/shogun/mathematics/ajd/JADiagOrth.cpp b/src/shogun/mathematics/ajd/JADiagOrth.cpp index 565d27e7556..fc60180d82c 100644 --- a/src/shogun/mathematics/ajd/JADiagOrth.cpp +++ b/src/shogun/mathematics/ajd/JADiagOrth.cpp @@ -19,18 +19,11 @@ void left_rot_stack(double *A, int M, int N, int K, int p, int q, double c, doub void right_rot_stack(double *A, int M, int N, int K, int p, int q, double c, double s); void left_rot_simple(double *A, int m, int n, int p, int q, double c, double s); -SGMatrix CJADiagOrth::diagonalize(const SGNDArray &C0, SGMatrix V0, +SGMatrix CJADiagOrth::diagonalize(SGNDArray C, SGMatrix V0, double eps, int itermax) { - int m = C0.dims[0]; - int L = C0.dims[2]; - - index_t * C_dims = SG_MALLOC(index_t, 3); - C_dims[0] = C0.dims[0]; - C_dims[1] = C0.dims[1]; - C_dims[2] = C0.dims[2]; - SGNDArray C(C_dims,3); - memcpy(C.array, C0.array, C0.dims[0]*C0.dims[1]*C0.dims[2]*sizeof(float64_t)); + int m = C.dims[0]; + int L = C.dims[2]; SGMatrix V; if (V0.num_rows != 0) diff --git a/src/shogun/mathematics/ajd/JADiagOrth.h b/src/shogun/mathematics/ajd/JADiagOrth.h index 1a2c6b4136e..e1a8ce78a3d 100644 --- a/src/shogun/mathematics/ajd/JADiagOrth.h +++ b/src/shogun/mathematics/ajd/JADiagOrth.h @@ -53,7 +53,7 @@ class CJADiagOrth : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - static SGMatrix diagonalize(const SGNDArray &C, + static SGMatrix diagonalize(SGNDArray C, SGMatrix V0 = SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200); @@ -65,7 +65,7 @@ class CJADiagOrth : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - virtual SGMatrix compute(const SGNDArray &C, + virtual SGMatrix compute(SGNDArray C, SGMatrix V0 = SGMatrix(NULL,0,0), double eps=CMath::MACHINE_EPSILON, int itermax=200) diff --git a/src/shogun/mathematics/ajd/UWedge.cpp b/src/shogun/mathematics/ajd/UWedge.cpp index ee2f69e6c70..629894fc846 100644 --- a/src/shogun/mathematics/ajd/UWedge.cpp +++ b/src/shogun/mathematics/ajd/UWedge.cpp @@ -14,7 +14,7 @@ typedef Matrix< float64_t, Dynamic, Dynamic, ColMajor > EMatrix; using namespace shogun; -SGMatrix CUWedge::diagonalize(const SGNDArray &C, SGMatrix V0, +SGMatrix CUWedge::diagonalize(SGNDArray C, SGMatrix V0, double eps, int itermax) { int d = C.dims[0]; diff --git a/src/shogun/mathematics/ajd/UWedge.h b/src/shogun/mathematics/ajd/UWedge.h index 5f4fc9799c2..c0443be8472 100644 --- a/src/shogun/mathematics/ajd/UWedge.h +++ b/src/shogun/mathematics/ajd/UWedge.h @@ -52,8 +52,8 @@ class CUWedge : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - static SGMatrix diagonalize(const SGNDArray &C, - SGMatrix V0= SGMatrix(NULL,0,0), + static SGMatrix diagonalize(SGNDArray C, + SGMatrix V0 = SGMatrix(NULL,0,0), double eps=1e-12, int itermax=200); @@ -64,8 +64,8 @@ class CUWedge : public CApproxJointDiagonalizer * @param itermax maximum number of iterations * @return V the matrix that best diagonalizes C */ - virtual SGMatrix compute(const SGNDArray &C, - SGMatrix V0= SGMatrix(NULL,0,0), + virtual SGMatrix compute(SGNDArray C, + SGMatrix V0 = SGMatrix(NULL,0,0), double eps=1e-12, int itermax=200) {