/* nag_trans_hessenberg_observer (g13ewc) Example Program. * * Copyright 1993 Numerical Algorithms Group * * Mark 3, 1993 * Mark 8 revised, 2004. */ #include #include #include #include #define A(I, J) a[(I) *tda + J] #define C(I, J) c[(I) *tdc + J] #define U(I, J) u[(I) *tdu + J] int main(void) { Integer exit_status = 0, i, j, n, p, tda, tdc, tdu; double *a = 0, *c = 0, one = 1.0, *u = 0, zero = 0.0; Nag_ObserverForm reduceto; NagError fail; INIT_FAIL(fail); printf( "nag_trans_hessenberg_observer (g13ewc) Example Program Results\n"); /* Skip the heading in the data file and read the data. */ scanf("%*[^\n]"); scanf("%ld%ld", &n, &p); if (n >= 1 || p >= 1) { if (!(a = NAG_ALLOC(n*n, double)) || !(c = NAG_ALLOC(p*n, double)) || !(u = NAG_ALLOC(n*n, double))) { printf("Allocation failure\n"); exit_status = -1; goto END; } tda = n; tdc = n; tdu = n; } else { printf("Invalid n or p.\n"); exit_status = 1; return exit_status; } reduceto = Nag_UH_Observer; for (j = 0; j < n; ++j) for (i = 0; i < n; ++i) scanf("%lf", &A(i, j)); for (i = 0; i < p; ++i) for (j = 0; j < n; ++j) scanf("%lf", &C(i, j)); if (u) /* Initialise U as the identity matrix. */ for (i = 0; i < n; ++i) { for (j = 0; j < n; ++j) U(i, j) = zero; U(i, i) = one; } /* Reduce the pair (A,C) to reduceto observer Hessenberg form. */ /* nag_trans_hessenberg_observer (g13ewc). * Unitary state-space transformation to reduce (AC) to * lower or upper observer Hessenberg form */ nag_trans_hessenberg_observer(n, p, reduceto, a, tda, c, tdc, u, tdu, &fail); if (fail.code != NE_NOERROR) { printf( "Error from nag_trans_hessenberg_observer (g13ewc).\n%s\n", fail.message); exit_status = 1; goto END; } printf("\nThe transformed state transition matrix is \n\n"); for (i = 0; i < n; ++i) { for (j = 0; j < n; ++j) printf("%8.4f ", A(i, j)); printf("\n"); } printf("\nThe transformed input matrix is \n\n"); for (i = 0; i < p; ++i) { for (j = 0; j < n; ++j) printf("%8.4f ", C(i, j)); printf("\n"); } if (u) { printf("\nThe transformation matrix that reduces (A,C) " "to observer Hessenberg form is \n\n"); for (i = 0; i < n; ++i) { for (j = 0; j < n; ++j) printf("%8.4f ", U(i, j)); printf("\n"); } } END: if (a) NAG_FREE(a); if (c) NAG_FREE(c); if (u) NAG_FREE(u); return exit_status; }