1 /*
2 * Copyright (c) 2019 Kevin Townsend
3 *
4 * SPDX-License-Identifier: Apache-2.0
5 */
6
7 #include <stdio.h>
8 #include <stdlib.h>
9 #include <string.h>
10 #include <zephyr/kernel.h>
11 #include <zsl/zsl.h>
12 #include <zsl/matrices.h>
13 #include "zsl/vectors.h"
14
15 #ifndef EPSILON
16 /* Define EPSILON for floating point comparisons. */
17 #define EPSILON 1e-6
18 #endif
19
svd_demo(void)20 void svd_demo(void)
21 {
22 size_t q = 18;
23 size_t p = 3;
24
25 zsl_real_t vi[18 * 3] = {
26 3.2, 5.0, -4.2,
27 -3.1, 5.4, -2.2,
28 7.7, 8.7, 9.9,
29 8.9, 0.6, 5.4,
30 3.7, -3.3, 7.2,
31 -5.4, -0.6, 8.4,
32 2.4, 5.1, -5.7,
33 6.9, -2.1, 0.4,
34 6.4, -9.1, -4.4,
35 0.1, 7.4, 0.0,
36 6.1, -2.3, 5.5,
37 6.1, -8.4, -6.6,
38 7.1, 2.3, 4.1,
39 -0.8, -4.7, 9.9,
40 -3.1, 1.2, 5.2,
41 6.4, -6.3, 5.2,
42 8.8, 7.3, 4.2,
43 -4.7, 0.0, -0.4
44 };
45
46 /* Large stack space warning. */
47 printf("WARNING: Operations involving large matrices and eigenvectors\n");
48 printf(" require a large amount of stack memory. This sample\n");
49 printf(" may not work on small embedded devices, and may\n");
50 printf(" cause a stack overflow.\n\n");
51
52 /* Define the 18x3 matrix, assigning the values from 'vi'. */
53 ZSL_MATRIX_DEF(mep, q, p);
54 zsl_mtx_from_arr(&mep, vi);
55
56 /* Define the three SVD output matrices. */
57 ZSL_MATRIX_DEF(u, mep.sz_rows, mep.sz_rows);
58 ZSL_MATRIX_DEF(e, mep.sz_rows, mep.sz_cols);
59 ZSL_MATRIX_DEF(v, mep.sz_cols, mep.sz_cols);
60
61 /* Define placeholder matrices for SVD and related operations. */
62 ZSL_MATRIX_DEF(ue, mep.sz_rows, mep.sz_cols);
63 ZSL_MATRIX_DEF(vt, mep.sz_cols, mep.sz_cols);
64 ZSL_MATRIX_DEF(uevt, mep.sz_rows, mep.sz_cols);
65
66 /* Display the input matrix for reference. */
67 printf("INPUT (18x3 MATRIX)\n");
68 printf("-------------------\n");
69 zsl_mtx_print(&mep);
70 printf("\n");
71
72 /* Calculate and display the SVD using 150 iterations. */
73 zsl_mtx_svd(&mep, &u, &e, &v, 150);
74 printf("SVD OUTPUT\n");
75 printf("----------\n");
76 printf("U matrix (18x18):\n");
77 zsl_mtx_print(&u);
78 printf("\n");
79 printf("Sigma matrix (18x3):\n");
80 zsl_mtx_print(&e);
81 printf("\n");
82 printf("V matrix: (3x3)\n");
83 zsl_mtx_print(&v);
84 printf("\n");
85
86 /* Multiply u, sigma and v(transposed) from the SVD calculations,
87 * which shoud return the original matrix. */
88 printf("VERIFICATION\n");
89 printf("------------\n");
90 zsl_mtx_mult(&u, &e, &ue);
91 zsl_mtx_trans(&v, &vt);
92 zsl_mtx_mult(&ue, &vt, &uevt);
93 printf("U * Sigma * V^t:\n");
94 zsl_mtx_print(&uevt);
95 printf("\n");
96 }
97
main(void)98 int main(void)
99 {
100 printf("\n\nzscilib SVD demo\n\n");
101
102 svd_demo();
103
104 return 0;
105 }
106
