/[ascend]/trunk/test.py
ViewVC logotype

Contents of /trunk/test.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1201 - (show annotations) (download) (as text)
Mon Jan 22 13:26:06 2007 UTC (12 years, 9 months ago) by johnpye
File MIME type: text/x-python
File size: 29345 byte(s)
Added some more test cases for DAE problem analysis.
Little bit of tidying up in the slv*.c files.
1 #!/usr/bin/env python
2 # ASCEND modelling environment
3 # Copyright (C) 2006 Carnegie Mellon University
4 #
5 # This program is free software; you can redistribute it and/or modify
6 # it under the terms of the GNU General Public License as published by
7 # the Free Software Foundation; either version 2, or (at your option)
8 # any later version.
9 #
10 # This program is distributed in the hope that it will be useful,
11 # but WITHOUT ANY WARRANTY; without even the implied warranty of
12 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 # GNU General Public License for more details.
14 #
15 # You should have received a copy of the GNU General Public License
16 # along with this program; if not, write to the Free Software
17 # Foundation, Inc., 59 Temple Place - Suite 330,
18 # Boston, MA 02111-1307, USA.
19
20 # This script gives a test suite for the high-level interface of ASCEND via
21 # Python. It is also planned to be a wrapper for the CUnit test suite, although
22 # this is still experimental.
23
24 import unittest
25 import os, sys
26 import math
27 import atexit
28
29 import platform
30 if platform.system() != "Windows":
31 import dl
32 sys.setdlopenflags(dl.RTLD_GLOBAL|dl.RTLD_NOW)
33
34 class Ascend(unittest.TestCase):
35
36 def setUp(self):
37 import ascpy
38 self.L = ascpy.Library()
39
40 def tearDown(self):
41 self.L.clear()
42 del self.L
43
44 class AscendSelfTester(Ascend):
45
46 def _run(self,modelname,solvername="QRSlv",filename=None):
47 if filename==None:
48 filename = 'johnpye/%s.a4c' % modelname
49 self.L.load(filename)
50 T = self.L.findType(modelname)
51 M = T.getSimulation('sim')
52 M.build()
53 M.solve(ascpy.Solver(solvername),ascpy.SolverReporter())
54 M.run(T.getMethod('self_test'))
55 return M
56
57 class TestCompiler(Ascend):
58
59 def testloading(self):
60 pass
61
62 def testsystema4l(self):
63 self.L.load('system.a4l')
64
65 def testatomsa4l(self):
66 self.L.load('atoms.a4l')
67
68 class TestSolver(AscendSelfTester):
69
70 def testlog10(self):
71 self._run('testlog10')
72
73 def testrootsofpoly(self):
74 self._run('roots_of_poly',filename="roots_of_poly.a4c")
75
76 def testcollapsingcan(self):
77 self._run('collapsingcan',filename="collapsingcan.a4c")
78
79 def testdistancecalc(self):
80 self._run('distance_calc',filename="distance_calc.a4c")
81
82 def testconopt(self):
83 self._run('testconopt',"CONOPT")
84
85 def testcmslv2(self):
86 self._run('testcmslv2',"CMSlv")
87
88 def testsunpos1(self):
89 self._run('example_1_6_1',"QRSlv","johnpye/sunpos.a4c")
90
91 def testsunpos2(self):
92 self._run('example_1_6_2',"QRSlv","johnpye/sunpos.a4c")
93
94 def testsunpos3(self):
95 self._run('example_1_7_1',"QRSlv","johnpye/sunpos.a4c")
96
97 def testsunpos4(self):
98 self._run('example_1_7_2',"QRSlv","johnpye/sunpos.a4c")
99
100 def testsunpos5(self):
101 self._run('example_1_7_3',"QRSlv","johnpye/sunpos.a4c")
102
103 def testsunpos6(self):
104 self._run('example_1_8_1',"QRSlv","johnpye/sunpos.a4c")
105
106 def testinstanceas(self):
107 M = self._run('example_1_6_1',"QRSlv","johnpye/sunpos.a4c")
108 self.assertAlmostEqual( float(M.t_solar), M.t_solar.as("s"))
109 self.assertAlmostEqual( float(M.t_solar)/3600, M.t_solar.as("h"))
110
111 class TestMatrix(AscendSelfTester):
112 def testlog10(self):
113 M = self._run('testlog10')
114 print M.getMatrix().write(sys.stderr,"mmio")
115
116
117 class TestIntegrator(Ascend):
118
119 def testListIntegrators(self):
120 I = ascpy.Integrator.getEngines()
121 s1 = sorted([str(i) for i in I.values()])
122 s2 = sorted(['IDA','LSODE','AWW'])
123 assert s1==s2
124
125 # this routine is reused by both testIDA and testLSODE
126 def _testIntegrator(self,integratorname):
127 self.L.load('johnpye/shm.a4c')
128 M = self.L.findType('shm').getSimulation('sim')
129 M.setSolver(ascpy.Solver('QRSlv'))
130 P = M.getParameters()
131 M.setParameter('feastol',1e-12)
132 print M.getChildren()
133 assert float(M.x) == 10.0
134 assert float(M.v) == 0.0
135 t_end = math.pi
136
137 I = ascpy.Integrator(M)
138 I.setReporter(ascpy.IntegratorReporterNull(I))
139 I.setEngine(integratorname);
140 I.setLinearTimesteps(ascpy.Units("s"), 0.0, t_end, 100);
141 I.setMinSubStep(0.0001); # these limits are required by IDA at present (numeric diff)
142 I.setMaxSubStep(0.1);
143 I.setInitialSubStep(0.001);
144 I.setMaxSubSteps(200);
145 if(integratorname=='IDA'):
146 I.setParameter('autodiff',False)
147 for p in M.getParameters():
148 print p.getName(),"=",p.getValue()
149 I.analyse();
150 I.solve();
151 print "At end of simulation,"
152 print "x = %f" % M.x
153 print "v = %f" % M.v
154 assert abs(float(M.x) + 10) < 1e-2
155 assert abs(float(M.v)) < 1e-2
156 assert I.getNumObservedVars() == 3
157
158 def testInvalidIntegrator(self):
159 self.L.load('johnpye/shm.a4c')
160 M = self.L.findType('shm').getSimulation('sim')
161 M.setSolver(ascpy.Solver('QRSlv'))
162 I = ascpy.Integrator(M)
163 try:
164 I.setEngine('___NONEXISTENT____')
165 except RuntimeError:
166 return
167 self.fail("setEngine did not raise error!")
168
169 def testLSODE(self):
170 self._testIntegrator('LSODE')
171
172 def testIDA(self):
173 self._testIntegrator('IDA')
174
175 def testparameters(self):
176 self.L.load('johnpye/shm.a4c')
177 M = self.L.findType('shm').getSimulation('sim')
178 M.build()
179 I = ascpy.Integrator(M)
180 I.setEngine('IDA')
181 P = I.getParameters()
182 for p in P:
183 print p.getName(),"=",p.getValue()
184 assert len(P)==11
185 assert P[0].isStr()
186 assert P[0].getName()=="linsolver"
187 assert P[0].getValue()=='SPGMR'
188 assert P[2].getName()=="autodiff"
189 assert P[2].getValue()==True
190 assert P[7].getName()=="atolvect"
191 assert P[7].getBoolValue() == True
192 P[2].setBoolValue(False)
193 assert P[2].getBoolValue()==False
194 I.setParameters(P)
195 assert I.getParameterValue('autodiff')==False
196 I.setParameter('autodiff',True)
197 try:
198 v = I.getParameterValue('nonexist')
199 except KeyError:
200 pass
201 else:
202 self.fail('Failed to trip invalid Integrator parameter')
203
204 class TestLSODE(Ascend):
205
206 def testzill(self):
207 self.L.load('johnpye/zill.a4c')
208 T = self.L.findType('zill')
209 M = T.getSimulation('sim')
210 M.setSolver(ascpy.Solver('QRSlv'))
211 I = ascpy.Integrator(M)
212 I.setEngine('LSODE')
213 I.setMinSubStep(1e-7)
214 I.setMaxSubStep(0.001)
215 I.setMaxSubSteps(10000)
216 I.setReporter(ascpy.IntegratorReporterConsole(I))
217 I.setLinearTimesteps(ascpy.Units(), 1.0, 1.5, 5)
218 I.analyse()
219 I.solve()
220 M.run(T.getMethod('self_test'))
221
222 def testnewton(self):
223 sys.stderr.write("STARTING TESTNEWTON\n")
224 self.L.load('johnpye/newton.a4c')
225 T = self.L.findType('newton')
226 M = T.getSimulation('sim')
227 M.solve(ascpy.Solver("QRSlv"),ascpy.SolverReporter())
228 I = ascpy.Integrator(M)
229 I.setEngine('LSODE')
230 I.setParameter('rtolvect',False)
231 I.setParameter('rtol',1e-7)
232 I.setParameter('atolvect',False)
233 I.setParameter('atol',1e-7)
234 I.setMinSubStep(1e-7)
235 I.setMaxSubStep(0.001)
236 I.setMaxSubSteps(10000)
237
238 I.setReporter(ascpy.IntegratorReporterConsole(I))
239 I.setLinearTimesteps(ascpy.Units("s"), 0, 2*float(M.v)/float(M.g), 2)
240 I.analyse()
241 I.solve()
242 print "At end of simulation,"
243 print "x = %f" % M.x
244 print "v = %f" % M.v
245 M.run(T.getMethod('self_test'))
246
247 def testlotka(self):
248 self.L.load('johnpye/lotka.a4c')
249 M = self.L.findType('lotka').getSimulation('sim')
250 M.setSolver(ascpy.Solver("QRSlv"))
251 I = ascpy.Integrator(M)
252 I.setEngine('LSODE')
253 I.setReporter(ascpy.IntegratorReporterConsole(I))
254 I.setLinearTimesteps(ascpy.Units("s"), 0, 200, 5)
255 I.analyse()
256 print "Number of vars = %d" % I.getNumVars()
257 assert I.getNumVars()==2
258 I.solve()
259 assert I.getNumObservedVars() == 3;
260 assert abs(M.R - 832) < 1.0
261 assert abs(M.F - 21.36) < 0.1
262
263 #-------------------------------------------------------------------------------
264 # Testing of a external blackbox functions
265
266 class TestBlackBox(AscendSelfTester):
267 def testparsefail0(self):
268 try:
269 self.L.load('test/blackbox/parsefail0.a4c')
270 self.fail("parsefail0 should not have loaded without errors")
271 except:
272 pass
273
274 def testparsefail1(self):
275 try:
276 self.L.load('test/blackbox/parsefail1.a4c')
277 self.fail("parsefail1 should not have loaded without errors")
278 except:
279 pass
280
281 def testparsefail2(self):
282 try:
283 self.L.load('test/blackbox/parsefail2.a4c')
284 self.fail("parsefail2 should not have loaded without errors")
285 except:
286 pass
287
288 def testparsefail3(self):
289 try:
290 self.L.load('test/blackbox/parsefail3.a4c')
291 self.fail("parsefail3 should not have loaded without errors")
292 except:
293 pass
294
295 def testparsefail4(self):
296 try:
297 self.L.load('test/blackbox/parsefail4.a4c')
298 self.fail("parsefail4 should not have loaded")
299 except:
300 pass
301
302 def testfail1(self):
303 """Mismatched arg counts check-- tests bbox, not ascend."""
304 self.L.load('test/blackbox/fail1.a4c')
305 try:
306 M = self.L.findType('fail1').getSimulation('sim')
307 self.fail("expected exception was not raised")
308 except RuntimeError,e:
309 print "Caught exception '%s', assumed ok" % e
310
311 def testfail2(self):
312 """Incorrect data arg check -- tests bbox, not ascend"""
313 self.L.load('test/blackbox/fail2.a4c')
314 try:
315 M = self.L.findType('fail2').getSimulation('sim')
316 self.fail("expected exception was not raised")
317 except RuntimeError,e:
318 print "Caught exception '%s', assumed ok (should mention errors during instantiation)" % e
319
320 def testpass1(self):
321 """simple single bbox forward solve"""
322 M = self._run('pass1',filename='test/blackbox/pass.a4c')
323
324 def testpass2(self):
325 """simple single bbox reverse solve"""
326 M = self._run('pass2',filename='test/blackbox/pass.a4c')
327
328 def testpass3(self):
329 """simple double bbox solve"""
330 M = self._run('pass3',filename='test/blackbox/pass3.a4c')
331
332 def testpass4(self):
333 """simple double bbox reverse solve"""
334 M = self._run('pass4',filename='test/blackbox/pass3.a4c')
335
336 def testpass5(self):
337 M = self._run('pass5',filename='test/blackbox/pass5.a4c')
338
339 def testpass6(self):
340 M = self._run('pass6',filename='test/blackbox/pass5.a4c')
341
342 def testpass7(self):
343 M = self._run('pass7',filename='test/blackbox/passmerge.a4c')
344
345 def testpass8(self):
346 M = self._run('pass8',filename='test/blackbox/passmerge.a4c')
347
348 def testpass9(self):
349 M = self._run('pass9',filename='test/blackbox/passmerge.a4c')
350
351 def testpass10(self):
352 M = self._run('pass10',filename='test/blackbox/passmerge.a4c')
353
354 def testpass11(self):
355 M = self._run('pass11',filename='test/blackbox/passmerge.a4c')
356
357 def testpass12(self):
358 M = self._run('pass12',filename='test/blackbox/passmerge.a4c')
359
360 # this test doesn't work: 'system is inconsistent' -- and structurally singular
361 # def testpass13(self):
362 # """cross-merged input/output solve"""
363 # M = self._run('pass13',filename='test/blackbox/passmerge.a4c')
364
365 def testpass14(self):
366 """cross-merged input/output reverse solve"""
367 M = self._run('pass14',filename='test/blackbox/passmerge.a4c')
368
369 def testpass20(self):
370 M = self._run('pass20',filename='test/blackbox/passarray.a4c')
371
372 def testparsefail21(self):
373 """dense array of black boxes wrong syntax"""
374 try:
375 self.L.load('test/blackbox/parsefail21.a4c')
376 self.fail("parsefail21 should not have loaded without errors")
377 except:
378 pass
379
380 def testpass22(self):
381 M = self._run('pass22',filename='test/blackbox/passarray.a4c')
382
383 def testpass23(self):
384 M = self._run('pass23',filename='test/blackbox/passarray.a4c')
385
386 def testpass61(self):
387 M = self._run('pass61',filename='test/blackbox/reinstantiate.a4c')
388
389 def testpass62(self):
390 M = self._run('pass62',filename='test/blackbox/reinstantiate.a4c')
391
392 def testpass64(self):
393 M = self._run('pass64',filename='test/blackbox/reinstantiate.a4c')
394
395 def testpass65(self):
396 M = self._run('pass65',filename='test/blackbox/reinstantiate.a4c')
397
398 def testpass66(self):
399 M = self._run('pass66',filename='test/blackbox/reinstantiate.a4c')
400
401 def testpass67(self):
402 M = self._run('pass67',filename='test/blackbox/reinstantiate.a4c')
403
404 class TestExtFn(AscendSelfTester):
405 def testextfntest(self):
406 M = self._run('extfntest',filename='johnpye/extfn/extfntest.a4c')
407 self.assertAlmostEqual(M.y, 2);
408 self.assertAlmostEqual(M.x, 1);
409 self.assertAlmostEqual(M.y, M.x + 1);
410
411 def testextrelfor(self):
412 M = self._run('extrelfor',filename='johnpye/extfn/extrelfor.a4c')
413
414 ## @TODO fix bug with badly-named bbox rel in a loop (Ben, maybe)
415 # def testextrelforbadnaming(self):
416 # self.L.load('johnpye/extfn/extrelforbadnaming.a4c')
417 # T = self.L.findType('extrelfor')
418 # M = T.getSimulation('sim')
419 # M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
420 # print "x[1] = %f" % M.x[1]
421 # print "x[2] = %f" % M.x[2]
422 # print "x[3] = %f" % M.x[3]
423 # print "x[4] = %f" % M.x[4]
424 # print "x[5] = %f" % M.x[5]
425 # M.run(T.getMethod('self_test'))
426
427 def testextrelrepeat(self):
428 M = self._run('extrelrepeat',filename='johnpye/extfn/extrelrepeat.a4c')
429
430 #-------------------------------------------------------------------------------
431 # Testing of Sensitivity module
432
433 class TestSensitivity(AscendSelfTester):
434 def test1(self):
435 self.L.load('sensitivity_test.a4c')
436 T = self.L.findType('sensitivity_test')
437 M = T.getSimulation('sim',False)
438 M.run(T.getMethod('on_load'))
439 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
440 M.run(T.getMethod('analyse'))
441 M.run(T.getMethod('self_test'))
442
443 # def testall(self):
444 # self.L.load('sensitivity_test.a4c')
445 # T = self.L.findType('sensitivity_test_all')
446 # M = T.getSimulation('sim',False)
447 # M.run(T.getMethod('on_load'))
448 # M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
449 # M.run(T.getMethod('analyse'))
450 # M.run(T.getMethod('self_test'))
451 # CAUSES CRASH
452
453 #-------------------------------------------------------------------------------
454 # Testing of a ExtPy - external python methods
455
456 class TestExtPy(AscendSelfTester):
457 def test1(self):
458 self.L.load('johnpye/extpy/extpytest.a4c')
459 T = self.L.findType('extpytest')
460 M = T.getSimulation('sim')
461 M.run(T.getMethod('self_test'))
462
463 def test2(self):
464 self.L.load('johnpye/extpy/extpytest.a4c')
465 T = self.L.findType('extpytest')
466 M = T.getSimulation('sim')
467 M.run(T.getMethod('pythonthing'))
468 M.run(T.getMethod('pythonthing'))
469 M.run(T.getMethod('pythonthing'))
470 M.run(T.getMethod('pythonthing'))
471 # causes crash!
472
473 #-------------------------------------------------------------------------------
474 # Testing of saturated steam properties library (iapwssatprops.a4c)
475
476 class TestSteam(AscendSelfTester):
477
478 def testiapwssatprops1(self):
479 M = self._run('testiapwssatprops1',filename='steam/iapwssatprops.a4c')
480 def testiapwssatprops2(self):
481 M = self._run('testiapwssatprops2',filename='steam/iapwssatprops.a4c')
482 def testiapwssatprops3(self):
483 M = self._run('testiapwssatprops3',filename='steam/iapwssatprops.a4c')
484
485 # test the stream model basically works
486 def testsatsteamstream(self):
487 M = self._run('satsteamstream',filename='steam/satsteamstream.a4c')
488
489 # test that we can solve in terms of various (rho,u)
490 def testsatuv(self):
491 self.L.load('steam/iapwssat.a4c')
492 T = self.L.findType('testiapwssatuv')
493 M = T.getSimulation('sim',False)
494 M.run(T.getMethod('on_load'))
495 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
496 print "p = %f bar" % M.p.as('bar');
497 print "T = %f C" % (M.T.as('K') - 273.15);
498 print "x = %f" % M.x;
499 M.run(T.getMethod('self_test'))
500 M.run(T.getMethod('values2'))
501 # M.v.setRealValueWithUnits(1.0/450,"m^3/kg");
502 # M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
503 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
504 print "p = %f bar" % M.p.as('bar');
505 print "T = %f C" % (M.T.as('K') - 273.15);
506 print "x = %f" % M.x;
507 M.run(T.getMethod('self_test2'))
508
509
510 ## @TODO fix error capture from bounds checking during initialisation
511 # def testiapwssat1(self):
512 # M = self._run('testiapwssat1',filename='steam/iapwssat.a4c')
513
514 def testdsgsat(self):
515 self.L.load('steam/dsgsat3.a4c')
516 T = self.L.findType('dsgsat3')
517 M = T.getSimulation('sim',False)
518 M.run(T.getMethod('on_load'))
519 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
520 self.assertAlmostEqual(M.dTw_dt[2],0.0);
521 M.run(T.getMethod('configure_dynamic'))
522 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
523 return M
524
525 def testintegLSODE(self):
526 M = self.testdsgsat()
527 M.qdot_s.setRealValueWithUnits(1000,"W/m")
528 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
529 #M.setParameter('
530 I = ascpy.Integrator(M)
531 I.setEngine('LSODE')
532 I.setParameter('meth','AM')
533 I.setParameter('maxord',12)
534 I.setReporter(ascpy.IntegratorReporterConsole(I))
535 I.setLinearTimesteps(ascpy.Units("s"), 0, 5, 1)
536 I.analyse()
537 I.solve()
538
539 def testintegIDA(self):
540 M = self.testdsgsat()
541 self.assertAlmostEqual(M.dTw_dt[2],0.0)
542 Tw1 = float(M.T_w[2])
543 T = self.L.findType('dsgsat3')
544 M.run(T.getMethod('free_states'))
545 I = ascpy.Integrator(M)
546 I.setEngine('IDA')
547 I.setParameter('linsolver','DENSE')
548 I.setParameter('safeeval',True)
549 I.setParameter('rtol',1e-8)
550 I.setInitialSubStep(0.01)
551 I.setMinSubStep(0.001)
552 I.setMaxSubSteps(100)
553 I.setReporter(ascpy.IntegratorReporterConsole(I))
554 I.setLinearTimesteps(ascpy.Units("s"), 0, 3600, 100)
555 I.analyse()
556 I.solve()
557 self.assertAlmostEqual(float(M.T_w[2]),Tw1)
558 M.qdot_s.setRealValueWithUnits(1000,"W/m")
559 self.assertAlmostEqual(M.qdot_s.as("W/m"),1000)
560 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
561 self.assertNotAlmostEqual(M.dTw_dt[2],0.0)
562 # I = ascpy.Integrator(M)
563 # I.setEngine('LSODE')
564 # I.setReporter(ascpy.IntegratorReporterConsole(I))
565 # I.setReporter(ascpy.IntegratorReporterConsole(I))
566 # I.setLinearTimesteps(ascpy.Units("s"), 0, 5, 100)
567 # I.setMinSubStep(0.0001)
568 # I.setMaxSubStep(100)
569 # I.setInitialSubStep(0.1)
570 # I.analyse()
571 # I.solve()
572
573 #-------------------------------------------------------------------------------
574 # Testing of freesteam external steam properties functions
575
576 with_freesteam = True
577 try:
578 # we assume that if the freesteam python module is installed, the ASCEND
579 # external library will also be.
580 import freesteam
581 have_freesteam = True
582 except ImportError,e:
583 have_freesteam = False
584
585 if with_freesteam and have_freesteam:
586 class TestFreesteam(AscendSelfTester):
587 # def testfreesteamtest(self):
588 # """run the self-test cases bundled with freesteam"""
589 # self._run('testfreesteam',filename='testfreesteam.a4c')
590
591 def testload(self):
592 """check that we can load 'thermalequilibrium2' (IMPORT "freesteam", etc)"""
593 self.L.load('johnpye/thermalequilibrium2.a4c')
594
595 def testinstantiate(self):
596 """load an instantiate 'thermalequilibrium2'"""
597 self.testload()
598 M = self.L.findType('thermalequilibrium2').getSimulation('sim')
599 return M
600
601 def testintegrate(self):
602 """integrate transfer of heat from one mass of water/steam to another
603 according to Newton's law of cooling"""
604 M = self.testinstantiate()
605 M.setSolver(ascpy.Solver("QRSlv"))
606 I = ascpy.Integrator(M)
607 I.setEngine('LSODE')
608 I.setReporter(ascpy.IntegratorReporterConsole(I))
609 I.setLinearTimesteps(ascpy.Units("s"), 0, 3000, 30)
610 I.setMinSubStep(0.01)
611 I.setInitialSubStep(1)
612 I.analyse()
613 print "Number of vars = %d" % I.getNumVars()
614 assert I.getNumVars()==2
615 I.solve()
616 assert I.getNumObservedVars() == 3;
617 print "S[1].T = %f K" % M.S[1].T
618 print "S[2].T = %f K" % M.S[2].T
619 print "Q = %f W" % M.Q
620 self.assertAlmostEqual(float(M.S[1].T),506.77225109,4);
621 self.assertAlmostEqual(float(M.S[2].T),511.605173967,5);
622 self.assertAlmostEqual(float(M.Q),-48.32922877329,3);
623 self.assertAlmostEqual(float(M.t),3000);
624 print "Note that the above values have not been verified analytically"
625
626 def testcollapsingcan2(self):
627 """ solve the collapsing can model using IAPWS-IF97 steam props """
628 M = self._run("collapsingcan2",filename="collapsingcan2.a4c");
629
630 #-------------------------------------------------------------------------------
631 # Testing of IDA's analysis module
632
633 class TestIDA(Ascend):
634 def _run(self,filen,modeln=""):
635 self.L.load('test/ida/%s.a4c' % filen)
636 T = self.L.findType('%s%s' % (filen,modeln))
637 M = T.getSimulation('sim')
638 M.build()
639 I = ascpy.Integrator(M)
640 I.setEngine('IDA')
641 I.analyse()
642
643 def _runfail(self,filen,n,msg="failed"):
644 try:
645 self._run(filen,'fail%d' % n)
646 except Exception,e:
647 print "(EXPECTED) ERROR: %s" % e
648 return
649 self.fail(msg)
650
651 def testsinglederiv(self):
652 self._run('singlederiv')
653
654 def testsinglederivfail1(self):
655 self._runfail('singlederiv',1
656 ,"t.ode_id=-1 did not trigger error")
657
658 def testsinglederivfail2(self):
659 self._runfail('singlederiv',2
660 ,"dy_dt.ode_id=2 did not trigger error")
661
662 def testsinglederivfail3(self):
663 self._runfail('singlederiv',3
664 ,"dy_dt.ode_type=3 did not trigger error")
665
666 def testsinglederivfail4(self):
667 self._runfail('singlederiv',4
668 ,"duplicate ode_type=1 did not trigger error")
669
670 def testsinglederivfail5(self):
671 self._runfail('singlederiv',5
672 ,"duplicate ode_type=1 did not trigger error")
673
674 def testsinglederivfail6(self):
675 self._runfail('singlederiv',6
676 ,"duplicate ode_type=1 did not trigger error")
677
678 def testtwoderiv(self):
679 self._run('twoderiv')
680
681 def testtwoderivfail1(self):
682 self._runfail('twoderiv',1)
683
684 def testtwoderivfail2(self):
685 self._runfail('twoderiv',2)
686
687 def testtwoderivfail3(self):
688 self._runfail('twoderiv',3)
689 def testtwoderivfail4(self):
690 self._runfail('twoderiv',4)
691 def testtwoderivfail5(self):
692 self._runfail('twoderiv',5)
693
694 def testnoderivs(self):
695 self._runfail('noderivs',1)
696
697 def testnoindeps(self):
698 self._runfail('indeps',1)
699
700 def testtwoindeps(self):
701 self._runfail('indeps',2)
702
703 def testfixedvars(self):
704 self._run('fixedvars')
705
706 def testfixedvars1(self):
707 self._run('fixedvars',1)
708
709 def testfixedvars2(self):
710 self._run('fixedvars',2)
711
712 def testfixedvars3(self):
713 self._run('fixedvars',3)
714
715
716 #-------------------------------------------------------------------------------
717 # Testing of IDA models using DENSE linear solver
718
719 class TestIDADENSE(Ascend):
720 """IDA DAE integrator, DENSE linear solver"""
721
722 def testnewton(self):
723 sys.stderr.write("STARTING TESTNEWTON\n")
724 self.L.load('johnpye/newton.a4c')
725 T = self.L.findType('newton')
726 M = T.getSimulation('sim')
727 M.solve(ascpy.Solver("QRSlv"),ascpy.SolverReporter())
728 I = ascpy.Integrator(M)
729 I.setEngine('IDA')
730 I.setParameter('linsolver','DENSE')
731 I.setParameter('safeeval',True)
732 I.setParameter('rtol',1e-8)
733 I.setMaxSubStep(0.001)
734 I.setMaxSubSteps(10000)
735
736 I.setReporter(ascpy.IntegratorReporterConsole(I))
737 I.setLinearTimesteps(ascpy.Units("s"), 0, 2*float(M.v)/float(M.g), 2)
738 I.analyse()
739 I.solve()
740 print "At end of simulation,"
741 print "x = %f" % M.x
742 print "v = %f" % M.v
743 M.run(T.getMethod('self_test'))
744
745 def testlotka(self):
746 self.L.load('johnpye/lotka.a4c')
747 M = self.L.findType('lotka').getSimulation('sim')
748 M.setSolver(ascpy.Solver("QRSlv"))
749 I = ascpy.Integrator(M)
750 I.setEngine('IDA')
751 I.setReporter(ascpy.IntegratorReporterConsole(I))
752 I.setLinearTimesteps(ascpy.Units("s"), 0, 200, 5);
753 I.setParameter('linsolver','DENSE')
754 I.setParameter('rtol',1e-8);
755 I.analyse()
756 assert I.getNumVars()==2
757 assert abs(M.R - 1000) < 1e-300
758 I.solve()
759 assert I.getNumObservedVars() == 3
760 assert abs(M.R - 832) < 1.0
761 assert abs(M.F - 21.36) < 0.1
762
763 def testdenx(self):
764 print "-----------------------------====="
765 self.L.load('johnpye/idadenx.a4c')
766 M = self.L.findType('idadenx').getSimulation('sim')
767 M.setSolver(ascpy.Solver("QRSlv"))
768 I = ascpy.Integrator(M)
769 I.setEngine('IDA')
770 I.setParameter('calcic','YA_YPD')
771 I.setParameter('linsolver','DENSE')
772 I.setReporter(ascpy.IntegratorReporterConsole(I))
773 I.setLogTimesteps(ascpy.Units("s"), 0.4, 4e10, 11)
774 I.setMaxSubStep(0);
775 I.setInitialSubStep(0)
776 I.setMaxSubSteps(0);
777 I.setParameter('autodiff',True)
778 I.analyse()
779 I.solve()
780 assert abs(float(M.y1) - 5.1091e-08) < 2e-9
781 assert abs(float(M.y2) - 2.0437e-13) < 2e-14
782 assert abs(float(M.y3) - 1.0) < 1e-5
783
784 ## @TODO fails during IDACalcIC (model too big?)
785 # def testkryx(self):
786 # self.L.load('johnpye/idakryx.a4c')
787 # ascpy.getCompiler().setUseRelationSharing(False)
788 # M = self.L.findType('idakryx').getSimulation('sim')
789 # M.setSolver(ascpy.Solver('QRSlv'))
790 # M.build()
791 # I = ascpy.Integrator(M)
792 # I.setEngine('IDA')
793 # I.setReporter(ascpy.IntegratorReporterConsole(I))
794 # I.setParameter('linsolver','DENSE')
795 # I.setParameter('maxl',8)
796 # I.setParameter('gsmodified',False)
797 # I.setParameter('autodiff',True)
798 # I.setParameter('rtol',0)
799 # I.setParameter('atol',1e-3);
800 # I.setParameter('atolvect',False)
801 # I.setParameter('calcic','YA_YDP')
802 # I.analyse()
803 # I.setLogTimesteps(ascpy.Units("s"), 0.01, 10.24, 11)
804 # I.solve()
805 # assert abs(M.u[2][2].getValue()) < 1e-5
806
807 #-------------------------------------------------------------------------------
808 # Testing of IDA models using SPGMR linear solver (Krylov)
809
810 # these tests are disabled until SPGMR preconditioning has been implemented
811 class TestIDASPGMR:#(Ascend):
812 def testlotka(self):
813 self.L.load('johnpye/lotka.a4c')
814 M = self.L.findType('lotka').getSimulation('sim')
815 M.setSolver(ascpy.Solver("QRSlv"))
816 I = ascpy.Integrator(M)
817 I.setEngine('IDA')
818 I.setReporter(ascpy.IntegratorReporterConsole(I))
819 I.setLinearTimesteps(ascpy.Units("s"), 0, 200, 5)
820 I.setParameter('rtol',1e-8)
821 I.analyse()
822 assert I.getNumVars()==2
823 assert abs(M.R - 1000) < 1e-300
824 I.solve()
825 assert I.getNumObservedVars() == 3
826 assert abs(M.R - 832) < 1.0
827 assert abs(M.F - 21.36) < 0.1
828
829
830 def testkryx(self):
831 self.L.load('johnpye/idakryx.a4c')
832 M = self.L.findType('idakryx').getSimulation('sim')
833 M.build()
834 I = ascpy.Integrator(M)
835 I.setEngine('IDA')
836 I.setReporter(ascpy.IntegratorReporterConsole(I))
837 I.setParameter('linsolver','SPGMR')
838 I.setParameter('prec','JACOBI')
839 I.setParameter('maxl',8)
840 I.setParameter('gsmodified',False)
841 I.setParameter('autodiff',True)
842 I.setParameter('gsmodified',True)
843 I.setParameter('rtol',0)
844 I.setParameter('atol',1e-3);
845 I.setParameter('atolvect',False)
846 I.setParameter('calcic','Y')
847 I.analyse()
848 I.setLogTimesteps(ascpy.Units("s"), 0.01, 10.24, 10);
849 print M.udot[1][3]
850 I.solve()
851 assert 0
852
853 def testzill(self):
854 self.L.load('johnpye/zill.a4c')
855 T = self.L.findType('zill')
856 M = T.getSimulation('sim')
857 M.setSolver(ascpy.Solver('QRSlv'))
858 I = ascpy.Integrator(M)
859 I.setEngine('IDA')
860 I.setParameter('safeeval',False)
861 I.setMinSubStep(1e-7)
862 I.setMaxSubStep(0.001)
863 I.setMaxSubSteps(10000)
864 I.setReporter(ascpy.IntegratorReporterConsole(I))
865 I.setLinearTimesteps(ascpy.Units(), 1.0, 1.5, 5)
866 I.analyse()
867 I.solve()
868 M.run(T.getMethod('self_test'))
869
870 def testdenxSPGMR(self):
871 self.L.load('johnpye/idadenx.a4c')
872 M = self.L.findType('idadenx').getSimulation('sim')
873 M.setSolver(ascpy.Solver('QRSlv'))
874 I = ascpy.Integrator(M)
875 I.setEngine('IDA')
876 I.setReporter(ascpy.IntegratorReporterConsole(I))
877 I.setLogTimesteps(ascpy.Units("s"), 0.4, 4e10, 11)
878 I.setMaxSubStep(0);
879 I.setInitialSubStep(0);
880 I.setMaxSubSteps(0);
881 I.setParameter('autodiff',True)
882 I.setParameter('linsolver','SPGMR')
883 I.setParameter('gsmodified',False)
884 I.setParameter('maxncf',10)
885 I.analyse()
886 I.solve()
887 assert abs(float(M.y1) - 5.1091e-08) < 1e-10
888 assert abs(float(M.y2) - 2.0437e-13) < 1e-15
889 assert abs(float(M.y3) - 1.0) < 1e-5
890
891 # move code above down here if you want to temporarily avoid testing it
892 class NotToBeTested:
893 def nothing(self):
894 pass
895
896 if __name__=='__main__':
897 # a whole bag of tricks to make sure we get the necessary dirs in our ascend, python and ld path vars
898 restart = 0
899
900 if platform.system()=="Windows":
901 LD_LIBRARY_PATH="PATH"
902 SEP = ";"
903 else:
904 LD_LIBRARY_PATH="LD_LIBRARY_PATH"
905 SEP = ":"
906
907 freesteamdir = os.path.expanduser("~/freesteam/ascend")
908 modeldirs = [os.path.abspath(os.path.join(sys.path[0],"models")),os.path.abspath(freesteamdir)]
909 if not os.environ.get('ASCENDLIBRARY'):
910 os.environ['ASCENDLIBRARY'] = SEP.join(modeldirs)
911 restart = 1
912 else:
913 envmodelsdir = [os.path.abspath(i) for i in os.environ['ASCENDLIBRARY'].split(SEP)]
914 for l in modeldirs:
915 if l in envmodelsdir[len(modeldirs):]:
916 envmodelsdir.remove(l)
917 restart = 1
918 for l in modeldirs:
919 if l not in envmodelsdir:
920 envmodelsdir.insert(0,l)
921 restart = 1
922 os.environ['ASCENDLIBRARY'] = SEP.join(envmodelsdir)
923
924 libdirs = ["pygtk","."]
925 libdirs = [os.path.normpath(os.path.join(sys.path[0],l)) for l in libdirs]
926 if not os.environ.get(LD_LIBRARY_PATH):
927 os.environ[LD_LIBRARY_PATH]=SEP.join(libdirs)
928 restart = 1
929 else:
930 envlibdirs = [os.path.normpath(i) for i in os.environ[LD_LIBRARY_PATH].split(SEP)]
931 for l in libdirs:
932 if l in envlibdirs[len(libdirs):]:
933 envlibdirs.remove(l)
934 restart = 1
935 for l in libdirs:
936 if l not in envlibdirs:
937 envlibdirs.insert(0,l)
938 restart = 1
939 os.environ[LD_LIBRARY_PATH] = SEP.join(envlibdirs)
940
941 pypath = os.path.normpath(os.path.join(sys.path[0],"pygtk"))
942 if not os.environ.get('PYTHONPATH'):
943 os.environ['PYTHONPATH']=pypath
944 else:
945 envpypath = os.environ['PYTHONPATH'].split(SEP)
946 if pypath not in envpypath:
947 envpypath.insert(0,pypath)
948 os.environ['PYTHONPATH']=SEP.join(envpypath)
949 restart = 1
950
951 if restart:
952 script = os.path.join(sys.path[0],"test.py")
953 print "Restarting with..."
954 print " export LD_LIBRARY_PATH=%s" % os.environ.get(LD_LIBRARY_PATH)
955 print " export PYTHONPATH=%s" % os.environ.get('PYTHONPATH')
956 print " export ASCENDLIBRARY=%s" % os.environ.get('ASCENDLIBRARY')
957
958 if sys.argv[0] !="/usr/bin/gdb":
959 os.execvp("python",[script] + sys.argv)
960 else:
961 print "Can't restart from inside GDB"
962
963 import ascpy
964
965 try:
966 import cunit
967 except:
968 pass
969
970 atexit.register(ascpy.shutdown)
971 #suite = unittest.TestSuite()
972 #suite = unittest.defaultTestLoader.loadTestsFromName('__main__')
973 #unittest.TextTestRunner(verbosity=2).run(suite)
974 unittest.main()

Properties

Name Value
svn:executable *

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