/[ascend]/trunk/models/test/transamp.a4c
ViewVC logotype

Contents of /trunk/models/test/transamp.a4c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1255 - (show annotations) (download) (as text)
Sat Jan 27 10:56:58 2007 UTC (15 years, 5 months ago) by johnpye
File MIME type: text/x-ascend
File size: 3443 byte(s)
Small error with the chemakzo test case.
1 REQUIRE "ivpsystem.a4l";
2 REQUIRE "atoms.a4l";
3 (*
4 This is the 'Transistor amplifier' problem from the Test Set for IVP Solvers
5
6 The problem is a stiff Differential Algebraic Equation of index 1,
7 consisting of 8 equations. The problem originates from electrical circuit
8 analysis, and is a model for the transistor amplifier. A complete
9 description of the problem could be found in transamp.pdf
10
11 http://pitagora.dm.uniba.it/~testset/problems/transamp.php
12
13 We're not really concerned with the details of this problem. We're just
14 trying to determine that our solver gets the right answer. High-precision
15 solution results are provided online and we aim to reproduce those here.
16
17 Note that, according to the graphs given at the above link, the DASSL
18 solver is not so great at this particular problem. Not sure why this is,
19 but it was certainly observed that it wasn't possible to obtain results
20 at the same precision as obtained for the other examples from the Test Set.
21 *)
22 MODEL transamp;
23
24 U_b, U_F, alpha, beta, C[1..5], R[0..9] IS_A real_constant;
25 FOR k IN [1..5] CREATE
26 C[k] :== k * 1e-6;
27 END FOR;
28
29 R[0] :== 1000;
30 FOR k IN [1..9] CREATE
31 R[k] :== 9000;
32 END FOR;
33
34 U_b :== 6;
35 U_F :== 0.026;
36 alpha :== 0.99;
37 beta :== 1e-6;
38
39 y[1..8], dy_dt[1..8] IS_A solver_var;
40 t IS_A solver_var;
41
42 U_e IS_A solver_var;
43 U_e = 0.1 * sin(200{PI}*t);
44
45 g23, g56 IS_A solver_var;
46 (* note, comments about possible overflow with this fn *)
47 g23 = beta * (exp((y[2] - y[3])/U_F) - 1);
48 g56 = beta * (exp((y[5] - y[6])/U_F) - 1);
49
50 -C[1]*dy_dt[1] + C[1]*dy_dt[2] = (y[1] - U_e) / R[0];
51 C[1]*dy_dt[1] - C[1]*dy_dt[2] = y[2] / R[1] + (y[2] - U_b) / R[2] + (1 - alpha)*g23;
52 -C[2]*dy_dt[3] = y[3] / R[3] - g23;
53 -C[3]*dy_dt[4] + C[3]*dy_dt[5] = (y[4] - U_b) / R[4] + alpha*g23;
54 C[3]*dy_dt[4] - C[3]*dy_dt[5] = y[5] / R[5] + (y[5] - U_b) / R[6] + (1 - alpha)*g56;
55 -C[4]*dy_dt[6] = y[6] / R[7] - g56;
56 -C[5]*dy_dt[7] + C[5]*dy_dt[8] = (y[7] - U_b) / R[8] + alpha*g56;
57 C[5]*dy_dt[7] - C[5]*dy_dt[8] = y[8] / R[9];
58
59 METHODS
60 METHOD values;
61 y[1] := 0;
62 y[2] := U_b/(R[2]/R[1] + 1);
63 y[3] := U_b/(R[2]/R[1] + 1);
64 y[4] := U_b;
65 y[5] := U_b/(R[6]/R[5] + 1);
66 y[6] := U_b/(R[6]/R[5] + 1);
67 y[7] := U_b;
68 y[8] := 0;
69
70 dy_dt[3] := -y[2]/(C[2]*R[3]);
71 dy_dt[6] := -y[5]/(C[4]*R[7]);
72
73 (* these were already calculated, numerically *)
74 dy_dt[1] := 51.338775;
75 dy_dt[2] := 51.338775;
76 dy_dt[4] := -24.9757667;
77 dy_dt[5] := -24.9757667;
78 dy_dt[7] := -10.00564453;
79 dy_dt[8] := -10.00564453;
80
81 t := 0 {s};
82 END values;
83
84 METHOD specify;
85 FIX t; (* t is incident in the equations, so we need to explicitly fix it *)
86 END specify;
87
88 METHOD ode_init;
89 FOR i IN [1..8] DO
90 y[i].ode_id := i; y[i].ode_type := 1;
91 dy_dt[i].ode_id := i; dy_dt[i].ode_type := 2;
92 END FOR;
93 FOR i IN [1..8] DO
94 y[i].obs_id := i;
95 END FOR;
96 t.ode_type := -1;
97 END ode_init;
98
99 METHOD on_load;
100 RUN reset;
101 RUN values;
102 RUN ode_init;
103 END on_load;
104
105 METHOD self_test;
106 ASSERT abs( y[1] + 0.5562145012262709e-002 ) < 1e-7;
107 ASSERT abs( y[2] - 3.006522471903042 ) < 1e-7;
108 ASSERT abs( y[3] - 2.849958788608128 ) < 1e-7;
109 ASSERT abs( y[4] - 2.926422536206241 ) < 1e-7;
110 ASSERT abs( y[5] - 2.704617865010554 ) < 1e-7;
111 ASSERT abs( y[6] - 2.761837778393145 ) < 1e-6;
112 ASSERT abs( y[7] - 4.770927631616772 ) < 1e-7;
113 ASSERT abs( y[8] - 1.236995868091548 ) < 1e-7;
114 END self_test;
115 END transamp;

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