← Back to branch summary

~ubuntu-branches/ubuntu/vivid/grass/vivid-proposed

« back to all changes in this revision

Viewing changes to lib/external/ccmath/qreval.c

  • Committer: Package Import Robot
  • Author(s): Bas Couwenberg
  • Date: 2015-02-20 23:12:08 UTC
  • mfrom: (8.2.6 experimental)
  • Revision ID: package-import@ubuntu.com-20150220231208-1u6qvqm84v430b10
Tags: 7.0.0-1~exp1
* New upstream release.
* Update python-ctypes-ternary.patch to use if/else instead of and/or.
* Drop check4dev patch, rely on upstream check.
* Add build dependency on libpq-dev to grass-dev for libpq-fe.h.
* Drop patches applied upstream, refresh remaining patches.
* Update symlinks for images switched from jpg to png.

Show diffs side-by-side

added added

removed removed

Lines of Context:
lib/external/ccmath/qreval.c
 
1
/*  qreval.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
#include "ccmath.h"
 
9
int qreval(double *ev, double *dp, int n)
 
10
{
 
11
    double cc, sc = 0.0, d, x, y, h, tzr = 1.e-15;
 
12
 
 
13
    int j, k, m, mqr = 8 * n;
 
14
 
 
15
    for (j = 0, m = n - 1;; ++j) {
 
16
        while (1) {
 
17
            if (m < 1)
 
18
                return 0;
 
19
            k = m - 1;
 
20
            if (fabs(dp[k]) <= fabs(ev[m]) * tzr)
 
21
                --m;
 
22
            else {
 
23
                x = (ev[k] - ev[m]) / 2.;
 
24
                h = sqrt(x * x + dp[k] * dp[k]);
 
25
                if (m > 1 && fabs(dp[m - 2]) > fabs(ev[k]) * tzr)
 
26
                    break;
 
27
                x += ev[m];
 
28
                ev[m--] = x - h;
 
29
                ev[m--] = x + h;
 
30
            }
 
31
        }
 
32
        if (j > mqr)
 
33
            return -1;
 
34
        if (x > 0.)
 
35
            d = ev[m] + x - h;
 
36
        else
 
37
            d = ev[m] + x + h;
 
38
        cc = 1.;
 
39
        y = 0.;
 
40
        ev[0] -= d;
 
41
        for (k = 0; k < m; ++k) {
 
42
            x = ev[k] * cc - y;
 
43
            y = dp[k] * cc;
 
44
            h = sqrt(x * x + dp[k] * dp[k]);
 
45
            if (k > 0)
 
46
                dp[k - 1] = sc * h;
 
47
            ev[k] = cc * h;
 
48
            cc = x / h;
 
49
            sc = dp[k] / h;
 
50
            ev[k + 1] -= d;
 
51
            y *= sc;
 
52
            ev[k] = cc * (ev[k] + y) + ev[k + 1] * sc * sc + d;
 
53
        }
 
54
        ev[k] = ev[k] * cc - y;
 
55
        dp[k - 1] = ev[k] * sc;
 
56
        ev[k] = ev[k] * cc + d;
 
57
    }
 
58
    return 0;
 
59
}