# Annotation of /trunk/models/johnpye/idakryx.a4c

Revision 956 - (hide annotations) (download) (as text)
Sat Dec 9 15:38:05 2006 UTC (13 years, 10 months ago) by johnpye
File MIME type: text/x-ascend
File size: 3011 byte(s)
```Removed references to 'ASC_BIG_BUGMAIL' and 'ASC_MILD_BUGMAIL'
Disabled stream redirection functions (they were making testing difficult)
Re-wrote the main test-runner to permit specific tests to be run from the commandline (ongoing)
Disabled readln tests that required stream redirection.
```
 1 johnpye 953 REQUIRE "ivpsystem.a4l"; 2 REQUIRE "atoms.a4l"; 3 4 (* 5 * ----------------------------------------------------------------- 6 * Example problem for IDA: 2D heat equation, serial, GMRES. 7 * 8 * This example solves a discretized 2D heat equation problem. 9 * This version uses the Krylov solver IDASpgmr. 10 * 11 * The DAE system solved is a spatial discretization of the PDE 12 * du/dt = d^2u/dx^2 + d^2u/dy^2 13 * on the unit square. The boundary condition is u = 0 on all edges. 14 * Initial conditions are given by u = 16 x (1 - x) y (1 - y). The 15 * PDE is treated with central differences on a uniform M x M grid. 16 * The values of u at the interior points satisfy ODEs, and 17 * equations u = 0 at the boundaries are appended, to form a DAE 18 * system of size N = M^2. Here M = 10. 19 * 20 * The system is solved with IDA/IDAS using the Krylov linear solver 21 * IDASPGMR. The preconditioner uses the diagonal elements of the 22 * Jacobian only. Routines for preconditioning, required by 23 * IDASPGMR, are supplied here. The constraints u >= 0 are posed 24 * for all components. Output is taken at t = 0, .01, .02, .04, 25 * ..., 10.24. Two cases are run -- with the Gram-Schmidt type 26 * being Modified in the first case, and Classical in the second. 27 * The second run uses IDAReInit and IDAReInitSpgmr. 28 * ----------------------------------------------------------------- 29 *) 30 MODEL idakryx; 31 32 M IS_A integer_constant; 33 johnpye 956 M :== 4; 34 johnpye 953 35 dx IS_A real_constant; 36 dx :== 1 / (M - 1); 37 38 coeff IS_A real_constant; 39 coeff :== 1/dx^2; 40 41 u[0..M-1][0..M-1] IS_A solver_var; 42 udot[0..M-1][0..M-1] IS_A solver_var; 43 44 t IS_A time; 45 46 (* central difference equations for internal points *) 47 FOR i IN [1..M-2] CREATE 48 FOR j IN [1..M-2] CREATE 49 udoteq[i][j]: udot[i][j] = coeff * ( 50 (* u_yy *) u[i-1][j] + u[i+1][j] - 2*u[i][j] 51 + (* u_xx *) u[i][j-1] + u[i][j+1] - 2*u[i][j] 52 ); 53 END FOR; 54 END FOR; 55 56 (* boundary conditions: constrained u = 0 *) 57 58 (* top and bottom *) 59 FOR i IN [0,M-1] CREATE 60 FOR j IN[0..M-1] CREATE 61 topbot[i][j]: u[i][j] = 0; 62 (* udot1[i][j]: udot[i][j] = 0; *) 63 END FOR; 64 END FOR; 65 66 (* sides *) 67 FOR i IN [1..M-2] CREATE 68 FOR j IN[0,M-1] CREATE 69 sides[i][j]: u[i][j] = 0; 70 (* udot2[i][j]: udot[i][j] = 0; *) 71 END FOR; 72 END FOR; 73 74 METHODS 75 METHOD on_load; 76 RUN specify; 77 RUN values; 78 RUN ode_init; 79 END on_load; 80 81 METHOD specify; 82 FIX udot[1..M-2][1..M-2]; 83 (* 84 FIX udot[0,M-1][0..M-1]; 85 FIX udot[1..M-2][0,M-1]; 86 *) 87 END specify; 88 89 90 METHOD values; 91 t := 0 {s}; 92 FOR i IN [0..M-1] DO 93 FOR j IN [0..M-1] DO 94 u[i][j] := 16 * (dx*j) * (1 - dx*j) * dx*i * (1 - dx*i); 95 END FOR; 96 END FOR; 97 johnpye 956 FOR i IN [1..M-2] DO 98 FOR j IN [1..M-2] DO 99 u[i][j] := 0; 100 END FOR; 101 END FOR; 102 johnpye 953 END values; 103 104 METHOD ode_init; 105 FREE udot[1..M-2][1..M-2]; 106 FREE udot[0,M-1][0..M-1]; 107 FREE udot[1..M-2][0,M-1]; 108 109 FOR i IN [0..M-1] DO 110 FOR j IN[0..M-1] DO 111 u[i][j].ode_id := i*M + j; 112 udot[i][j].ode_id := i*M + j; 113 114 u[i][j].ode_type := 1; 115 udot[i][j].ode_type := 2; 116 END FOR; 117 END FOR; 118 119 t.ode_type := -1; 120 121 END ode_init; 122 123 END idakryx;

 john.pye@anu.edu.au ViewVC Help Powered by ViewVC 1.1.22