| 1 | /* Generated by TAPENADE (INRIA, Ecuador team)
|
|---|
| 2 | Tapenade 3.12 (r6213) - 13 Oct 2016 10:54
|
|---|
| 3 | */
|
|---|
| 4 | #include "GlobalDeclarations_bv.c"
|
|---|
| 5 | #include "DIFFSIZES.inc"
|
|---|
| 6 | /* Hint: NBDirsMax should be the maximum number of differentiation directions
|
|---|
| 7 | */
|
|---|
| 8 |
|
|---|
| 9 | /*
|
|---|
| 10 | Differentiation of o_fcn in reverse (adjoint) mode:
|
|---|
| 11 | gradient of useful results: *obj
|
|---|
| 12 | with respect to varying inputs: *obj x[0:6-1]
|
|---|
| 13 | RW status of diff variables: *obj:in-out x[0:6-1]:out
|
|---|
| 14 | Plus diff mem management of: obj:in x:in
|
|---|
| 15 | */
|
|---|
| 16 | // = 2.0/sqrt3;
|
|---|
| 17 | //double tsqrt3 = 2.0/sqrt3;
|
|---|
| 18 | void o_fcn_bv(double *obj, double (*objb)[NBDirsMax], double x[6], double xb[6
|
|---|
| 19 | ][NBDirsMax], int nbdirs) {
|
|---|
| 20 | /* const
|
|---|
| 21 | static */
|
|---|
| 22 | double matr[4], f;
|
|---|
| 23 | double matrb[4][NBDirsMax], fb[NBDirsMax];
|
|---|
| 24 | /* static */
|
|---|
| 25 | double g;
|
|---|
| 26 | double gb[NBDirsMax];
|
|---|
| 27 | double retval;
|
|---|
| 28 | int nd;
|
|---|
| 29 | int ii1;
|
|---|
| 30 | double tempb[NBDirsMax];
|
|---|
| 31 | double tempb0[NBDirsMax];
|
|---|
| 32 | double tempb1[NBDirsMax];
|
|---|
| 33 | double o_fcn;
|
|---|
| 34 | /* Calculate M = A*inv(W). */
|
|---|
| 35 | matr[0] = x[1] - x[0];
|
|---|
| 36 | matr[1] = (2.0*x[2]-x[1]-x[0])*.577350269189625797959429519858;
|
|---|
| 37 | matr[2] = x[4] - x[3];
|
|---|
| 38 | matr[3] = (2.0*x[5]-x[4]-x[3])*.577350269189625797959429519858;
|
|---|
| 39 | /* Calculate det(M). */
|
|---|
| 40 | g = matr[0]*matr[3] - matr[1]*matr[2];
|
|---|
| 41 | if (g <= 1.00000e-14) {
|
|---|
| 42 | for (nd = 0; nd < nbdirs; ++nd) {
|
|---|
| 43 | gb[nd] = *objb[nd];
|
|---|
| 44 | *objb[nd] = 0.0;
|
|---|
| 45 | }
|
|---|
| 46 | for (nd = 0; nd < nbdirs; ++nd)
|
|---|
| 47 | for (ii1 = 0; ii1 < 4; ++ii1)
|
|---|
| 48 | matrb[ii1][nd] = 0.0;
|
|---|
| 49 | } else {
|
|---|
| 50 | /* Calculate norm(M). */
|
|---|
| 51 | f = matr[0]*matr[0] + matr[1]*matr[1] + matr[2]*matr[2] + matr[3]*matr
|
|---|
| 52 | [3];
|
|---|
| 53 | /* Calculate objective function. */
|
|---|
| 54 | for (nd = 0; nd < nbdirs; ++nd) {
|
|---|
| 55 | tempb1[nd] = .500000000000000000000000000000*(*objb[nd])/g;
|
|---|
| 56 | fb[nd] = tempb1[nd];
|
|---|
| 57 | gb[nd] = -(f*tempb1[nd]/g);
|
|---|
| 58 | *objb[nd] = 0.0;
|
|---|
| 59 | for (ii1 = 0; ii1 < 4; ++ii1)
|
|---|
| 60 | matrb[ii1][nd] = 0.0;
|
|---|
| 61 | matrb[0][nd] = matrb[0][nd] + 2*matr[0]*fb[nd];
|
|---|
| 62 | matrb[1][nd] = matrb[1][nd] + 2*matr[1]*fb[nd];
|
|---|
| 63 | matrb[2][nd] = matrb[2][nd] + 2*matr[2]*fb[nd];
|
|---|
| 64 | matrb[3][nd] = matrb[3][nd] + 2*matr[3]*fb[nd];
|
|---|
| 65 | }
|
|---|
| 66 | }
|
|---|
| 67 | for (nd = 0; nd < nbdirs; ++nd) {
|
|---|
| 68 | matrb[0][nd] = matrb[0][nd] + matr[3]*gb[nd];
|
|---|
| 69 | matrb[3][nd] = matrb[3][nd] + matr[0]*gb[nd];
|
|---|
| 70 | matrb[1][nd] = matrb[1][nd] - matr[2]*gb[nd];
|
|---|
| 71 | matrb[2][nd] = matrb[2][nd] - matr[1]*gb[nd];
|
|---|
| 72 | for (ii1 = 0; ii1 < 6; ++ii1)
|
|---|
| 73 | xb[ii1][nd] = 0.0;
|
|---|
| 74 | tempb[nd] = .577350269189625797959429519858*matrb[3][nd];
|
|---|
| 75 | xb[5][nd] = xb[5][nd] + 2.0*tempb[nd];
|
|---|
| 76 | xb[4][nd] = xb[4][nd] - tempb[nd];
|
|---|
| 77 | xb[3][nd] = xb[3][nd] - tempb[nd];
|
|---|
| 78 | matrb[3][nd] = 0.0;
|
|---|
| 79 | xb[4][nd] = xb[4][nd] + matrb[2][nd];
|
|---|
| 80 | xb[3][nd] = xb[3][nd] - matrb[2][nd];
|
|---|
| 81 | matrb[2][nd] = 0.0;
|
|---|
| 82 | tempb0[nd] = .577350269189625797959429519858*matrb[1][nd];
|
|---|
| 83 | xb[2][nd] = xb[2][nd] + 2.0*tempb0[nd];
|
|---|
| 84 | xb[1][nd] = xb[1][nd] - tempb0[nd];
|
|---|
| 85 | xb[0][nd] = xb[0][nd] - tempb0[nd];
|
|---|
| 86 | matrb[1][nd] = 0.0;
|
|---|
| 87 | xb[1][nd] = xb[1][nd] + matrb[0][nd];
|
|---|
| 88 | xb[0][nd] = xb[0][nd] - matrb[0][nd];
|
|---|
| 89 | }
|
|---|
| 90 | }
|
|---|
| 91 | /****************************************************************************
|
|---|
| 92 | This set of functions reference triangular elements to an equilateral
|
|---|
| 93 | triangle. The input are the coordinates in the following order:
|
|---|
| 94 | [x1 x2 x3 y1 y2 y3]
|
|---|
| 95 | A zero return value indicates success, while a nonzero value indicates
|
|---|
| 96 | failure.
|
|---|
| 97 | ***************************************************************************
|
|---|
| 98 | Not all compilers substitute out constants (especially the square root).
|
|---|
| 99 | Therefore, they are substituted out manually. The values below were
|
|---|
| 100 | calculated on a solaris machine using long doubles. I believe they are
|
|---|
| 101 | accurate.
|
|---|
| 102 | ***************************************************************************
|
|---|
| 103 | #define sqrt3 .577350269189625797959429519858 1.0/sqrt(3.0)
|
|---|
| 104 | #define tsqrt3 1.15470053837925159591885903972 2.0/sqrt(3.0) */
|
|---|
| 105 | double tsqrt3;
|
|---|
