http://github.com/b45ch1/hpsc_hanoi_2009_walter
C | 108 lines | 98 code | 1 blank | 9 comment | 8 complexity | dce60ff27967c900523946897f8f3956 MD5 | raw file
```  1/*        Generated by TAPENADE     (INRIA, Tropics team)
2    Tapenade 3.4 (r3375) - 10 Feb 2010 15:08
3*/
4#include "GlobalDeclarations_b.c"
5
6/*
7  Differentiation of inv in reverse (adjoint) mode:
8   gradient     of useful results: *qt
9   with respect to varying inputs: *qt *a
10   RW status of diff variables: *qt:in-out *a:out
11*/
12void inv_b(double *a, double *ab, double *qt, double *qtb, double *r, double *
13        rb, int na) {
14    int n, m, k;
15    double rnn, rnm;
16    double rnnb, rnmb;
17    int res;
18    int res0;
19    int res1;
20    double tmp;
21    int res2;
22    int res3;
23    double tmp0;
24    int res4;
25    int res5;
26    int res6;
27    int res7;
28    int res8;
29    double tmp1;
31    double tmpb;
32    double tmp0b;
33    double tmp1b;
34    pushreal8_(*r);
35    pushreal8_(*qt);
36    qr(a, qt, r, na);
37    for (n = na-1; n >= 0; --n) {
38        res = myindex(n, n, na);
39        pushreal8_(rnn);
40        rnn = r[res];
41        for (m = 0; m <= na-1; ++m) {
42            res0 = myindex(n, m, na);
43            res1 = myindex(n, m, na);
44            tmp = r[res1]/rnn;
45            pushreal8_(r[res0]);
46            r[res0] = tmp;
47            res2 = myindex(n, m, na);
48            res3 = myindex(n, m, na);
49            tmp0 = qt[res3]/rnn;
50            pushreal8_(qt[res2]);
51            qt[res2] = tmp0;
52        }
53        adFrom = n + 1;
54        for (m = adFrom; m <= na-1; ++m) {
55            res4 = myindex(n, m, na);
56            pushreal8_(rnm);
57            rnm = r[res4];
58            res5 = myindex(n, m, na);
59            pushreal8_(r[res5]);
60            r[res5] = 0;
61            for (k = 0; k <= na-1; ++k) {
62                res6 = myindex(n, k, na);
63                res7 = myindex(n, k, na);
64                res8 = myindex(m, k, na);
65                tmp1 = qt[res7] - qt[res8]*rnm;
66                pushreal8_(qt[res6]);
67                qt[res6] = tmp1;
68            }
69        }
71    }
72    *rb = 0.0;
73    for (n = 0; n <= na-1; ++n) {
75        for (m = na-1; m >= adFrom; --m) {
76            rnmb = 0.0;
77            for (k = na-1; k >= 0; --k) {
78                popreal8_(&qt[res6]);
79                tmp1b = qtb[res6];
80                qtb[res6] = 0.0;
81                qtb[res7] = qtb[res7] + tmp1b;
82                qtb[res8] = qtb[res8] - rnm*tmp1b;
83                rnmb = rnmb - qt[res8]*tmp1b;
84            }
85            popreal8_(&r[res5]);
86            rb[res5] = 0.0;
87            popreal8_(&rnm);
88            rb[res4] = rb[res4] + rnmb;
89        }
90        rnnb = 0.0;
91        for (m = na-1; m >= 0; --m) {
92            popreal8_(&qt[res2]);
93            tmp0b = qtb[res2];
94            qtb[res2] = 0.0;
95            qtb[res3] = qtb[res3] + tmp0b/rnn;
96            popreal8_(&r[res0]);
97            tmpb = rb[res0];
98            rnnb = rnnb - r[res1]*tmpb/(rnn*rnn) - qt[res3]*tmp0b/(rnn*rnn);
99            rb[res0] = 0.0;
100            rb[res1] = rb[res1] + tmpb/rnn;
101        }
102        popreal8_(&rnn);
103        rb[res] = rb[res] + rnnb;
104    }
105    popreal8_(qt);
106    popreal8_(r);
107    qr_b(a, ab, qt, qtb, r, rb, na);
108}
```