RcdMathLib_doc
Open Source Library for Linear and Non-linear Algebra
ieee_mobile_computing_non_lin_alg.c
Go to the documentation of this file.
1 /*
2  * Copyright (C) 2019 Zakaria Kasmi <zkasmi@inf.fu-berlin.de>
3  *
4  * This file is subject to the terms and conditions of the GNU Lesser General
5  * Public License v2.1. See the file LICENSE in the top level directory for more
6  * details.
7  */
8 
23 #include <stdio.h>
24 #include <math.h>
25 
26 #include "matrix.h"
27 #include "vector.h"
28 #include "modified_gauss_newton.h"
29 #include "levenberg_marquardt.h"
30 
46 void ieee_mobile_comp_get_exp_f(vector_t x_vec[], vector_t data_vec[],
47  vector_t f_vec[])
48 {
49  vector_t x1, x2;
50  uint8_t i;
51 
52  x1 = x_vec[0];
53  x2 = x_vec[1];
54  for (i = 1; i < 9; i++) {
55  f_vec[i - 1] = x1 * exp(i * x2) - data_vec[i - 1];
56  }
57 }
58 
83 {
84  vector_t x1, x2;
85  uint8_t i;
86 
87  x1 = x_vec[0];
88  x2 = x_vec[1];
89  for (i = 1; i < 9; i++) {
90  J[i - 1][0] = exp(i * x2);
91  J[i - 1][1] = i * x1 * exp(i * x2);
92  }
93 }
94 
95 /*
96  *
97  * SOLVING NONLINEAR LEAST-SQUARES PROBLEMS
98  * WITH THE GAUSS-NEWTON AND LEVENBERG-MARQUARDT METHODS
99  *
100  * ALFONSO CROEZE, LINDSEY PITTMAN, AND WINNIE REYNOLDS
101  *
102  */
103 void ieee_mobile_comp_exponential_data_test(void)
104 {
105  vector_t d_vec[] = { 8.3, 11.0, 14.7, 19.7, 26.7, 35.2, 44.4, 55.9 };
106  vector_t start_x_vec[] = { 6, .3 };
107  vector_t f_vec[8];
108  matrix_t J[8][2];
109 
110  // Parameters of GNM & LVM
111  matrix_t beta0 = 0.2;
112  matrix_t beta1 = 0.8;
113  matrix_t eps = 1e-3;
114  matrix_t tau = 1e-6;
115  matrix_t f_min = 1e-11;
116  uint8_t max_it_num = 3;
117  uint8_t iter_num = 0;
118  matrix_t est_x_vec[2];
119 
120  ieee_mobile_comp_get_exp_f(start_x_vec, d_vec, f_vec);
121  ieee_mobile_comp_get_exp_Jacobian(start_x_vec, J);
122 
123  puts("############### Test Gauss-Newton & LVM algorithms ###############");
124  puts("\nExponential Data:");
125  // Start timer
126  iter_num = modified_gauss_newton(8, 2, start_x_vec, d_vec, eps, f_min,
127  max_it_num, est_x_vec, &ieee_mobile_comp_get_exp_f,
129  // Stop timer
130  printf("start value = {%.6f, %.6f}\n", start_x_vec[0], start_x_vec[1]);
131  printf("Gauss-Newton = {%.6f, %.6f}\n", est_x_vec[0], est_x_vec[1]);
132  printf("iteration number = %d\n", iter_num);
133 
134  iter_num = 0;
135  vector_clear(2, est_x_vec);
136  // Start timer
137  iter_num = opt_levenberg_marquardt(8, 2, start_x_vec, d_vec, eps, tau,
138  beta0, beta1, max_it_num, est_x_vec,
141  // Stop timer
142  printf("Levenberg-Marquardt = {%.6f, %.6f}\n", est_x_vec[0], est_x_vec[1]);
143  printf("iteration number = %d\n", iter_num);
144 
145 }
146 
165  vector_t f_vec[])
166 {
167  vector_t x1, x2, x3, x4;
168  uint8_t i;
169 
170  x1 = x_vec[0];
171  x2 = x_vec[1];
172  x3 = x_vec[2];
173  x4 = x_vec[3];
174  for (i = 1; i < 13; i++) {
175  f_vec[i - 1] = x1 * sin(x2 * i + x3) + x4 - data_vec[i - 1];
176  }
177 }
178 
200 {
201  vector_t x1, x2, x3;
202  uint8_t i;
203 
204  x1 = x_vec[0];
205  x2 = x_vec[1];
206  x3 = x_vec[2];
207  for (i = 1; i < 13; i++) {
208  J[i - 1][0] = sin(i * x2 + x3);
209  J[i - 1][1] = i * x1 * cos(i * x2 + x3);
210  J[i - 1][2] = x1 * cos(i * x2 + x3);
211  J[i - 1][3] = 1;
212  }
213 }
214 
215 void ieee_mobile_comp_sinusoidal_data_test(void)
216 {
217  vector_t d_vec[] = { 61, 65, 72, 78, 85, 90, 92, 92, 88, 81, 72, 63 };
218  vector_t start_x_vec[] = { 17, 0.5, 10.5, 77 };
219  vector_t f_vec[12];
220  matrix_t J[12][4];
221 
222  // Parameters of GNM & LVM
223  matrix_t beta0 = 0.2;
224  matrix_t beta1 = 0.8;
225  matrix_t eps = 1e-3;
226  matrix_t tau = 1e-6;
227  matrix_t f_min = 1e-11;
228  uint8_t max_it_num = 2;
229  uint8_t iter_num = 0;
230  matrix_t est_x_vec[4];
231 
232  ieee_mobile_comp_get_sin_f(start_x_vec, d_vec, f_vec);
233  ieee_mobile_comp_get_sin_Jacobian(start_x_vec, J);
234 
235 // puts("############### Test Gauss-Newton & LVM algorithms ###############");
236  puts("\nSinusoidal Data:");
237  // Start timer
238  iter_num = modified_gauss_newton(12, 4, start_x_vec, d_vec, eps, f_min,
239  max_it_num, est_x_vec, &ieee_mobile_comp_get_sin_f,
241  // Stop timer
242  printf("start value = {%.6f, %.6f, %.6f, %.6f}\n", start_x_vec[0],
243  start_x_vec[1], start_x_vec[2], start_x_vec[3]);
244  printf("Gauss-Newton = {%.6f, %.6f, %.6f, %.6f}\n", est_x_vec[0],
245  est_x_vec[1], est_x_vec[2], est_x_vec[3]);
246  printf("iteration number = %d\n", iter_num);
247 
248  iter_num = 0;
249  vector_clear(4, est_x_vec);
250  // Start timer
251 iter_num = opt_levenberg_marquardt(12, 4, start_x_vec, d_vec, eps, tau,
252  beta0, beta1, max_it_num, est_x_vec, &ieee_mobile_comp_get_sin_f,
254 
255  // Stop timer
256  printf("Levenberg-Marquardt = {%.6f, %.6f, %.6f, %.6f}\n", est_x_vec[0],
257  est_x_vec[1], est_x_vec[2], est_x_vec[3]);
258  printf("iteration number = %d\n", iter_num);
259 
260 }
modified_gauss_newton
uint8_t modified_gauss_newton(uint8_t f_length, uint8_t n, vector_t x0_vec[n], vector_t data_vec[f_length], matrix_t eps, matrix_t fmin, uint8_t max_iter_num, vector_t est_x_vec[n], void(*get_f_error)(vector_t x0_vec[], vector_t data_vec[], vector_t f_vec[]), void(*get_jacobian)(vector_t x0_vec[], matrix_t J[][n]))
Implements the modified Gauss–Newton algorithm.
Definition: modified_gauss_newton.c:32
vector_t
#define vector_t
Define the data type of the vector elements.
Definition: vector.h:33
ieee_mobile_comp_get_exp_Jacobian
void ieee_mobile_comp_get_exp_Jacobian(vector_t x_vec[], matrix_t J[][2])
Calculate the Jacobian matrix using exponential data.
Definition: ieee_mobile_computing_non_lin_alg.c:82
matrix.h
Matrix computations.
ieee_mobile_comp_get_sin_Jacobian
void ieee_mobile_comp_get_sin_Jacobian(vector_t x_vec[], matrix_t J[][4])
Calculate the Jacobian matrix using sinusoidal data.
Definition: ieee_mobile_computing_non_lin_alg.c:199
levenberg_marquardt.h
Implement the Levenberg–Marquardt (LVM) algorithm.
vector_clear
void vector_clear(uint8_t size, vector_t arr[])
Clear all the elements of the vector.
Definition: vector.c:32
opt_levenberg_marquardt
uint8_t opt_levenberg_marquardt(uint8_t f_length, uint8_t n, vector_t x0_vec[n], vector_t data_vec[f_length], matrix_t eps, matrix_t tau, matrix_t beta0, matrix_t beta1, uint8_t max_iter_num, vector_t est_x_vec[n], void(*get_f_error)(vector_t x0_vec[], vector_t data_vec[], vector_t f_vec[]), void(*get_jacobian)(vector_t x0_vec[], matrix_t J[][n]))
Implements the Levenberg–Marquardt (LVM) algorithm.
Definition: levenberg_marquardt.c:143
matrix_t
#define matrix_t
Define the data type of the matrix elements.
Definition: matrix.h:38
ieee_mobile_comp_get_exp_f
void ieee_mobile_comp_get_exp_f(vector_t x_vec[], vector_t data_vec[], vector_t f_vec[])
Calculate the error vector using exponential data.
Definition: ieee_mobile_computing_non_lin_alg.c:46
ieee_mobile_comp_get_sin_f
void ieee_mobile_comp_get_sin_f(vector_t x_vec[], vector_t data_vec[], vector_t f_vec[])
Calculate the error vector using sinusoidal data.
Definition: ieee_mobile_computing_non_lin_alg.c:164
modified_gauss_newton.h
Implement the Gauss–Newton algorithm.
vector.h
Vector computations.