atovm.c 1.2 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344
  1. /* atovm.c CCMATH mathematics library source code.
  2. *
  3. * Copyright (C) 2000 Daniel A. Atkinson All rights reserved.
  4. * This code may be redistributed under the terms of the GNU library
  5. * public license (LGPL). ( See the lgpl.license file for details.)
  6. * ------------------------------------------------------------------------
  7. */
  8. void atovm(double *v, int n)
  9. {
  10. double *p0, *q0, *p, *q, *qq;
  11. double h, s;
  12. int i, j, k, mm;
  13. q0 = v + n * n - 1;
  14. *q0 = 1.;
  15. q0 -= n + 1;
  16. p0 = v + n * n - n - n - 1;
  17. for (i = n - 2, mm = 1; i >= 0; --i, p0 -= n + 1, q0 -= n + 1, ++mm) {
  18. if (i && *(p0 - 1) != 0.) {
  19. for (j = 0, p = p0, h = 1.; j < mm; ++j, ++p)
  20. h += *p * *p;
  21. h = *(p0 - 1);
  22. *q0 = 1. - h;
  23. for (j = 0, q = q0 + n, p = p0; j < mm; ++j, q += n)
  24. *q = -h * *p++;
  25. for (k = i + 1, q = q0 + 1; k < n; ++k) {
  26. for (j = 0, qq = q + n, p = p0, s = 0.; j < mm; ++j, qq += n)
  27. s += *qq * *p++;
  28. s *= h;
  29. for (j = 0, qq = q + n, p = p0; j < mm; ++j, qq += n)
  30. *qq -= s * *p++;
  31. *q++ = -s;
  32. }
  33. }
  34. else {
  35. *q0 = 1.;
  36. for (j = 0, p = q0 + 1, q = q0 + n; j < mm; ++j, q += n)
  37. *q = *p++ = 0.;
  38. }
  39. }
  40. }