Paparazzi UAS  v7.0_unstable
Paparazzi is a free software Unmanned Aircraft System.
r8lib_min.c
Go to the documentation of this file.
1 /*
2  * This file is a modified subset of the R8lib from John Burkardt.
3  * http://people.sc.fsu.edu/~jburkardt/c_src/r8lib/r8lib.html
4  *
5  * It is the minimal set of functions from r8lib needed to use qr_solve.
6  *
7  * This code is distributed under the GNU LGPL license.
8  */
9 
10 # include "r8lib_min.h"
11 #include "std.h"
12 # include <stdlib.h>
13 # include <math.h>
14 
15 #define DEBUG_FPRINTF(...)
16 #define DEBUG_PRINT(...)
17 #define DEBUG_EXIT(...)
18 
19 void r8mat_copy_new ( int m, int n, float a1[], float a2[])
20 
21 /******************************************************************************/
22 /*
23  Purpose:
24 
25  R8MAT_COPY_NEW copies one R8MAT to a "new" R8MAT.
26 
27  Discussion:
28 
29  An R8MAT is a doubly dimensioned array of R8 values, stored as a vector
30  in column-major order.
31 
32  Licensing:
33 
34  This code is distributed under the GNU LGPL license.
35 
36  Modified:
37 
38  26 July 2008
39 
40  Author:
41 
42  John Burkardt
43 
44  Parameters:
45 
46  Input, int M, N, the number of rows and columns.
47 
48  Input, float A1[M*N], the matrix to be copied.
49 
50  Output, float R8MAT_COPY_NEW[M*N], the copy of A1.
51 */
52 {
53  int i;
54  int j;
55 
56  /*a2 = ( float * ) malloc ( m * n * sizeof ( float ) );*/
57 
58  for ( j = 0; j < n; j++ )
59  {
60  for ( i = 0; i < m; i++ )
61  {
62  a2[i+j*m] = a1[i+j*m];
63  }
64  }
65 }
66 /******************************************************************************/
67 
68 float r8_epsilon ( void )
69 
70 /******************************************************************************/
71 /*
72  Purpose:
73 
74  R8_EPSILON returns the R8 round off unit.
75 
76  Discussion:
77 
78  R8_EPSILON is a number R which is a power of 2 with the property that,
79  to the precision of the computer's arithmetic,
80  1 < 1 + R
81  but
82  1 = ( 1 + R / 2 )
83 
84  Licensing:
85 
86  This code is distributed under the GNU LGPL license.
87 
88  Modified:
89 
90  01 September 2012
91 
92  Author:
93 
94  John Burkardt
95 
96  Parameters:
97 
98  Output, float R8_EPSILON, the R8 round-off unit.
99 */
100 {
101  const float value = 1.192092896E-7;
102 
103  return value;
104 }
105 /******************************************************************************/
106 
107 float r8mat_amax ( int m, int n, float a[] )
108 
109 /******************************************************************************/
110 /*
111  Purpose:
112 
113  R8MAT_AMAX returns the maximum absolute value entry of an R8MAT.
114 
115  Discussion:
116 
117  An R8MAT is a doubly dimensioned array of R8 values, stored as a vector
118  in column-major order.
119 
120  Licensing:
121 
122  This code is distributed under the GNU LGPL license.
123 
124  Modified:
125 
126  07 September 2012
127 
128  Author:
129 
130  John Burkardt
131 
132  Parameters:
133 
134  Input, int M, the number of rows in A.
135 
136  Input, int N, the number of columns in A.
137 
138  Input, float A[M*N], the M by N matrix.
139 
140  Output, float R8MAT_AMAX, the maximum absolute value entry of A.
141 */
142 {
143  int i;
144  int j;
145  float value;
146 
147  value = fabs ( a[0+0*m] );
148 
149  for ( j = 0; j < n; j++ )
150  {
151  for ( i = 0; i < m; i++ )
152  {
153  if ( value < fabs ( a[i+j*m] ) )
154  {
155  value = fabs ( a[i+j*m] );
156  }
157  }
158  }
159  return value;
160 }
161 /******************************************************************************/
162 
163 float r8_sign ( float x )
164 
165 /******************************************************************************/
166 /*
167  Purpose:
168 
169  R8_SIGN returns the sign of an R8.
170 
171  Licensing:
172 
173  This code is distributed under the GNU LGPL license.
174 
175  Modified:
176 
177  08 May 2006
178 
179  Author:
180 
181  John Burkardt
182 
183  Parameters:
184 
185  Input, float X, the number whose sign is desired.
186 
187  Output, float R8_SIGN, the sign of X.
188 */
189 {
190  float value;
191 
192  if ( x < 0.0 )
193  {
194  value = - 1.0;
195  }
196  else
197  {
198  value = + 1.0;
199  }
200  return value;
201 }
202 /******************************************************************************/
203 
204 float r8_max ( float x, float y )
205 
206 /******************************************************************************/
207 /*
208  Purpose:
209 
210  R8_MAX returns the maximum of two R8's.
211 
212  Licensing:
213 
214  This code is distributed under the GNU LGPL license.
215 
216  Modified:
217 
218  07 May 2006
219 
220  Author:
221 
222  John Burkardt
223 
224  Parameters:
225 
226  Input, float X, Y, the quantities to compare.
227 
228  Output, float R8_MAX, the maximum of X and Y.
229 */
230 {
231  float value;
232 
233  if ( y < x )
234  {
235  value = x;
236  }
237  else
238  {
239  value = y;
240  }
241  return value;
242 }
243 /******************************************************************************/
244 
245 float *r8mat_l_solve ( int n, float a[], float b[] )
246 
247 /******************************************************************************/
248 /*
249  Purpose:
250 
251  R8MAT_L_SOLVE solves a lower triangular linear system.
252 
253  Discussion:
254 
255  An R8MAT is a doubly dimensioned array of R8 values, stored as a vector
256  in column-major order.
257 
258  Licensing:
259 
260  This code is distributed under the GNU LGPL license.
261 
262  Modified:
263 
264  07 June 2008
265 
266  Author:
267 
268  John Burkardt
269 
270  Parameters:
271 
272  Input, int N, the number of rows and columns of
273  the matrix A.
274 
275  Input, float A[N*N], the N by N lower triangular matrix.
276 
277  Input, float B[N], the right hand side of the linear system.
278 
279  Output, float R8MAT_L_SOLVE[N], the solution of the linear system.
280 */
281 {
282  float dot;
283  int i;
284  int j;
285  float *x;
286 
287  x = ( float * ) malloc ( n * sizeof ( float ) );
288 /*
289  Solve L * x = b.
290 */
291  for ( i = 0; i < n; i++ )
292  {
293  dot = 0.0;
294  for ( j = 0; j < i; j++ )
295  {
296  dot = dot + a[i+j*n] * x[j];
297  }
298  x[i] = ( b[i] - dot ) / a[i+i*n];
299  }
300 
301  return x;
302 }
303 /******************************************************************************/
304 
305 float *r8mat_lt_solve ( int n, float a[], float b[] )
306 
307 /******************************************************************************/
308 /*
309  Purpose:
310 
311  R8MAT_LT_SOLVE solves a transposed lower triangular linear system.
312 
313  Discussion:
314 
315  An R8MAT is a doubly dimensioned array of R8 values, stored as a vector
316  in column-major order.
317 
318  Given the lower triangular matrix A, the linear system to be solved is:
319 
320  A' * x = b
321 
322  Licensing:
323 
324  This code is distributed under the GNU LGPL license.
325 
326  Modified:
327 
328  08 April 2009
329 
330  Author:
331 
332  John Burkardt
333 
334  Parameters:
335 
336  Input, int N, the number of rows and columns of the matrix A.
337 
338  Input, float A[N*N], the N by N lower triangular matrix.
339 
340  Input, float B[N], the right hand side of the linear system.
341 
342  Output, float R8MAT_LT_SOLVE[N], the solution of the linear system.
343 */
344 {
345  int i;
346  int j;
347  float *x;
348 
349  x = ( float * ) malloc ( n * sizeof ( float ) );
350 
351  for ( j = n-1; 0 <= j; j-- )
352  {
353  x[j] = b[j];
354  for ( i = j+1; i < n; i++ )
355  {
356  x[j] = x[j] - x[i] * a[i+j*n];
357  }
358  x[j] = x[j] / a[j+j*n];
359  }
360 
361  return x;
362 }
363 /******************************************************************************/
364 
365 float *r8mat_mtv_new ( int m, int n, float a[], float x[] )
366 
367 /******************************************************************************/
368 /*
369  Purpose:
370 
371  R8MAT_MTV_NEW multiplies a transposed matrix times a vector.
372 
373  Discussion:
374 
375  An R8MAT is a doubly dimensioned array of R8 values, stored as a vector
376  in column-major order.
377 
378  For this routine, the result is returned as the function value.
379 
380  Licensing:
381 
382  This code is distributed under the GNU LGPL license.
383 
384  Modified:
385 
386  26 August 2011
387 
388  Author:
389 
390  John Burkardt
391 
392  Parameters:
393 
394  Input, int M, N, the number of rows and columns of the matrix.
395 
396  Input, float A[M,N], the M by N matrix.
397 
398  Input, float X[M], the vector to be multiplied by A.
399 
400  Output, float R8MAT_MTV_NEW[N], the product A'*X.
401 */
402 {
403  int i;
404  int j;
405  float *y;
406 
407  y = ( float * ) malloc ( n * sizeof ( float ) );
408 
409  for ( j = 0; j < n; j++ )
410  {
411  y[j] = 0.0;
412  for ( i = 0; i < m; i++ )
413  {
414  y[j] = y[j] + a[i+j*m] * x[i];
415  }
416  }
417 
418  return y;
419 }
420 /******************************************************************************/
421 
422 float r8vec_max ( int n, float r8vec[] )
423 
424 /******************************************************************************/
425 /*
426  Purpose:
427 
428  R8VEC_MAX returns the value of the maximum element in a R8VEC.
429 
430  Licensing:
431 
432  This code is distributed under the GNU LGPL license.
433 
434  Modified:
435 
436  05 May 2006
437 
438  Author:
439 
440  John Burkardt
441 
442  Parameters:
443 
444  Input, int N, the number of entries in the array.
445 
446  Input, float R8VEC[N], a pointer to the first entry of the array.
447 
448  Output, float R8VEC_MAX, the value of the maximum element. This
449  is set to 0.0 if N <= 0.
450 */
451 {
452  int i;
453  float value;
454 
455  if ( n <= 0 )
456  {
457  value = 0.0;
458  return value;
459  }
460 
461  value = r8vec[0];
462 
463  for ( i = 1; i < n; i++ )
464  {
465  if ( value < r8vec[i] )
466  {
467  value = r8vec[i];
468  }
469  }
470  return value;
471 }
472 /******************************************************************************/
473 
474 int i4_min ( int i1, int i2 )
475 
476 /******************************************************************************/
477 /*
478  Purpose:
479 
480  I4_MIN returns the smaller of two I4's.
481 
482  Licensing:
483 
484  This code is distributed under the GNU LGPL license.
485 
486  Modified:
487 
488  29 August 2006
489 
490  Author:
491 
492  John Burkardt
493 
494  Parameters:
495 
496  Input, int I1, I2, two integers to be compared.
497 
498  Output, int I4_MIN, the smaller of I1 and I2.
499 */
500 {
501  int value;
502 
503  if ( i1 < i2 )
504  {
505  value = i1;
506  }
507  else
508  {
509  value = i2;
510  }
511  return value;
512 }
513 /******************************************************************************/
514 
515 int i4_max ( int i1, int i2 )
516 
517 /******************************************************************************/
518 /*
519  Purpose:
520 
521  I4_MAX returns the maximum of two I4's.
522 
523  Licensing:
524 
525  This code is distributed under the GNU LGPL license.
526 
527  Modified:
528 
529  29 August 2006
530 
531  Author:
532 
533  John Burkardt
534 
535  Parameters:
536 
537  Input, int I1, I2, are two integers to be compared.
538 
539  Output, int I4_MAX, the larger of I1 and I2.
540 */
541 {
542  int value;
543 
544  if ( i2 < i1 )
545  {
546  value = i1;
547  }
548  else
549  {
550  value = i2;
551  }
552  return value;
553 }
554 /******************************************************************************/
float r8_max(float x, float y)
Definition: r8lib_min.c:204
float * r8mat_lt_solve(int n, float a[], float b[])
Definition: r8lib_min.c:305
int i4_max(int i1, int i2)
Definition: r8lib_min.c:515
float * r8mat_mtv_new(int m, int n, float a[], float x[])
Definition: r8lib_min.c:365
float r8_epsilon(void)
Definition: r8lib_min.c:68
float r8mat_amax(int m, int n, float a[])
Definition: r8lib_min.c:107
float r8vec_max(int n, float r8vec[])
Definition: r8lib_min.c:422
void r8mat_copy_new(int m, int n, float a1[], float a2[])
Definition: r8lib_min.c:19
float r8_sign(float x)
Definition: r8lib_min.c:163
float * r8mat_l_solve(int n, float a[], float b[])
Definition: r8lib_min.c:245
int i4_min(int i1, int i2)
Definition: r8lib_min.c:474
float b
Definition: wedgebug.c:202