GRASS Programmer's Manual  6.5.svn(2014)-r66266
ccmath_grass_wrapper.c
Go to the documentation of this file.
1 #if defined(HAVE_CCMATH)
2 #include <ccmath.h>
3 #else
4 #include <grass/ccmath_grass.h>
5 #endif
6
163 int G_math_solv(double **a,double *b,int n)
164 {
165  return solv(a[0],b, n);
166 }
167
168
177  int G_math_solvps(double **a,double *b,int n)
178 {
179  return solvps(a[0], b,n);
180 }
181
182
193 void G_math_solvtd(double *a,double *b,double *c,double *x,int m)
194 {
195  solvtd(a, b, c, x, m);
196  return;
197 }
198
199
200 /*
201  \brief Solve an upper right triangular linear system T*x = b.
202
203  \param a = pointer to array of upper right triangular matrix T
204  \param b = pointer to array of system vector The computation overloads this with the solution vector x.
205  \param n = dimension (dim(a)=n*n,dim(b)=n)
206  \return value: f = status flag, with 0 -> normal exit, -1 -> system singular
207 */
208 int G_math_solvru(double **a,double *b,int n)
209 {
210  return solvru(a[0], b, n);
211 }
212
213
221 int G_math_minv(double **a,int n)
222 {
223  return minv(a[0], n);
224 }
225
226
235 int G_math_psinv(double **a,int n)
236 {
237  return psinv( a[0], n);
238 }
239
240
248 int G_math_ruinv(double **a,int n)
249 {
250  return ruinv(a[0], n);
251 }
252
253
254 /*
255 -----------------------------------------------------------------------------
256
257  Symmetric Eigensystem Analysis:
258 -----------------------------------------------------------------------------
259 */
268 void G_math_eigval(double **a,double *ev,int n)
269 {
270  eigval(a[0], ev, n);
271  return;
272 }
273
274
289 void G_math_eigen(double **a,double *ev,int n)
290 {
291  eigen(a[0], ev, n);
292  return;
293 }
294
295
296 /*
297  \brief Compute the maximum (absolute) eigenvalue and corresponding eigenvector of a real symmetric matrix A.
298
299
300  \param a = array containing symmetric input matrix A
301  \param u = array containing the n components of the eigenvector at exit (vector normalized to 1)
302  \param n = dimension of system
303  \return: ev = eigenvalue of A with maximum absolute value HUGE -> convergence failure
304 */
305 double G_math_evmax(double **a,double *u,int n)
306 {
307  return evmax(a[0], u, n);
308 }
309
310
311 /*
312 ------------------------------------------------------------------------------
313
314  Singular Value Decomposition:
315 ------------------------------------------------------------------------------
316
317  A number of versions of the Singular Value Decomposition (SVD)
318  are implemented in the library. They support the efficient
319  computation of this important factorization for a real m by n
320  matrix A. The general form of the SVD is
321
322  A = U*S*V~ with S = | D |
323  | 0 |
324
325  where U is an m by m orthogonal matrix, V is an n by n orthogonal matrix,
326  D is the n by n diagonal matrix of singular value, and S is the singular
327  m by n matrix produced by the transformation.
328
329  The singular values computed by these functions provide important
330  information on the rank of the matrix A, and on several matrix
331  norms of A. The number of non-zero singular values d[i] in D
332  equal to the rank of A. The two norm of A is
333
334  ||A|| = max(d[i]) , and the condition number is
335
336  k(A) = max(d[i])/min(d[i]) .
337
338  The Frobenius norm of the matrix A is
339
340  Fn(A) = Sum(i=0 to n-1) d[i]^2 .
341
342  Singular values consistent with zero are easily recognized, since
343  the decomposition algorithms have excellent numerical stability.
344  The value of a 'zero' d[i] is no larger than a few times the
345  computational rounding error e.
346
347  The matrix U1 is formed from the first n orthonormal column vectors
348  of U. U1[i,j] = U[i,j] for i = 1 to m and j = 1 to n. A singular
349  value decomposition of A can also be expressed in terms of the m by\
350  n matrix U1, with
351
352  A = U1*D*V~ .
353
354  SVD functions with three forms of output are provided. The first
355  form computes only the singular values, while the second computes
356  the singular values and the U and V orthogonal transformation
357  matrices. The third form of output computes singular values, the
358  V matrix, and saves space by overloading the input array with
359  the U1 matrix.
360
361  Two forms of decomposition algorithm are available for each of the
362  three output types. One is computationally efficient when m ~ n.
363  The second, distinguished by the prefix 'sv2' in the function name,
364  employs a two stage Householder reduction to accelerate computation
365  when m substantially exceeds n. Use of functions of the second form
366  is recommended for m > 2n.
367
368  Singular value output from each of the six SVD functions satisfies
369
370  d[i] >= 0 for i = 0 to n-1.
371 -------------------------------------------------------------------------------
372 */
373
374
386 int G_math_svdval(double *d,double **a,int m,int n)
387 {
388  return svdval(d, a[0], m, n);
389 }
390
391
402 int G_math_sv2val(double *d,double **a,int m,int n)
403 {
404  return sv2val(d, a[0], m, n);
405 }
406
407
408 /*
409  \brief Compute the singular value transformation S = U~*A*V.
410
411  \param d = pointer to double array of dimension n (output = singular values of A)
412  \param a = pointer to store of the m by n input matrix A (A is altered by the computation)
413  \param u = pointer to store for m by m orthogonal matrix U
414  \param v = pointer to store for n by n orthogonal matrix V
415  \param m = number of rows in A
416  \param n = number of columns in A (m>=n required)
417  \return value: status flag with: 0 -> success -1 -> input error m < n
418 */
419 int G_math_svduv(double *d,double **a,double **u,int m,double **v,int n)
420 {
421  return svduv(d, a[0], u[0], m, v[0], n);
422 }
423
424
436 int G_math_sv2uv(double *d,double **a,double **u,int m,double **v,int n)
437 {
438  return sv2uv(d, a[0], u[0], m, v[0], n);
439 }
440
441
455 int G_math_svdu1v(double *d,double **a,int m,double **v,int n)
456 {
457  return svdu1v(d, a[0], m, v[0], n);
458 }
int svduv(double *d, double *a, double *u, int m, double *v, int n)
Definition: svduv.c:10
int G_math_solvru(double **a, double *b, int n)
float b
Definition: named_colr.c:8
int solvps(double *s, double *x, int n)
Definition: solvps.c:9
int G_math_solvps(double **a, double *b, int n)
Solve a symmetric positive definite linear system S*x = b.
void G_math_eigen(double **a, double *ev, int n)
Compute the eigenvalues and eigenvectors of a real symmetric matrix A.
int svdval(double *d, double *a, int m, int n)
Definition: svdval.c:10
double G_math_evmax(double **a, double *u, int n)
void eigen(double *a, double *eval, int n)
Definition: eigen.c:10
void G_math_eigval(double **a, double *ev, int n)
Compute the eigenvalues of a real symmetric matrix A.
int G_math_sv2uv(double *d, double **a, double **u, int m, double **v, int n)
Compute the singular value transformation when m &gt;&gt; n.
void solvtd(double *a, double *b, double *c, double *x, int m)
Definition: solvtd.c:8
int G_math_sv2val(double *d, double **a, int m, int n)
Compute singular values when m &gt;&gt; n.
void G_math_solvtd(double *a, double *b, double *c, double *x, int m)
Solve a tridiagonal linear system M*x = y.
int sv2val(double *d, double *a, int m, int n)
Definition: sv2val.c:10
int minv(double *a, int n)
Definition: minv.c:10
int psinv(double *v, int n)
Definition: psinv.c:9
int G_math_svduv(double *d, double **a, double **u, int m, double **v, int n)
int G_math_psinv(double **a, int n)
Invert (in place) a symmetric real matrix, V -&gt; Inv(V).
int sv2uv(double *d, double *a, double *u, int m, double *v, int n)
Definition: sv2uv.c:10
int solvru(double *a, double *b, int n)
Definition: solvru.c:8
int svdu1v(double *d, double *a, int m, double *v, int n)
Definition: svdu1v.c:10
double evmax(double *a, double *u, int n)
Definition: evmax.c:10
int G_math_svdval(double *d, double **a, int m, int n)
Compute the singular values of a real m by n matrix A.
int G_math_svdu1v(double *d, double **a, int m, double **v, int n)
Compute the singular value transformation with A overloaded by the partial U-matrix.
int G_math_ruinv(double **a, int n)
Invert an upper right triangular matrix T -&gt; Inv(T).
int G_math_minv(double **a, int n)
Invert (in place) a general real matrix A -&gt; Inv(A).
int ruinv(double *a, int n)
Definition: ruinv.c:8
int solv(double *a, double *b, int n)
Definition: solv.c:10
int G_math_solv(double **a, double *b, int n)
Solve a general linear system A*x = b.
int n