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

Contents of /trunk/test.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1270 - (show annotations) (download) (as text)
Sun Feb 4 01:01:18 2007 UTC (15 years, 10 months ago) by johnpye
File MIME type: text/x-python
File size: 34870 byte(s)
Reduced some debug output.
Added ATOM for specific_energy_rate.
Fixed default_self routines for DSG models.
Added testdsgsatrepeat test to check above.
1 #!/usr/bin/env python
2 # ASCEND modelling environment
3 # Copyright (C) 2006, 2007 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,parameters={}):
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.setSolver(ascpy.Solver(solvername))
53 for k,v in parameters.iteritems():
54 M.setParameter(k,v)
55 M.solve(ascpy.Solver(solvername),ascpy.SolverReporter())
56 M.run(T.getMethod('self_test'))
57 return M
58
59 class TestCompiler(Ascend):
60
61 def _run(self,filen,modeln=""):
62 self.L.load('test/compiler/%s.a4c' % filen)
63 T = self.L.findType('%s%s' % (filen,modeln))
64 M = T.getSimulation('sim')
65 M.build()
66
67 def _runfail(self,filen,n,msg="failed"):
68 try:
69 self._run(filen,'fail%d' % n)
70 except Exception,e:
71 print "(EXPECTED) ERROR: %s" % e
72 return
73 self.fail(msg)
74
75
76 def testloading(self):
77 """library startup"""
78 pass
79
80 def testsystema4l(self):
81 """loading system.a4l?"""
82 self.L.load('system.a4l')
83
84 def testatomsa4l(self):
85 """loading atoms.a4l?"""
86 self.L.load('atoms.a4l')
87
88 def testmissingreq(self):
89 """flagging a missing REQUIRE"""
90 self._runfail('missingreq',1)
91
92 def testmissingsubreq(self):
93 """flagging a subsidiary missing REQUIRE"""
94 self._runfail('missingreq',1)
95
96 class TestSolver(AscendSelfTester):
97
98 def testlog10(self):
99 self._run('testlog10')
100
101 def testrootsofpoly(self):
102 self._run('roots_of_poly',filename="roots_of_poly.a4c")
103
104 def testcollapsingcan(self):
105 self._run('collapsingcan',filename="collapsingcan.a4c")
106
107 def testdistancecalc(self):
108 self._run('distance_calc',filename="distance_calc.a4c")
109
110 def testconopt(self):
111 self._run('testconopt',"CONOPT")
112
113 def testcmslv2(self):
114 self._run('testcmslv2',"CMSlv")
115
116 def testsunpos1(self):
117 self._run('example_1_6_1',"QRSlv","johnpye/sunpos.a4c")
118
119 def testsunpos2(self):
120 self._run('example_1_6_2',"QRSlv","johnpye/sunpos.a4c")
121
122 def testsunpos3(self):
123 self._run('example_1_7_1',"QRSlv","johnpye/sunpos.a4c")
124
125 def testsunpos4(self):
126 self._run('example_1_7_2',"QRSlv","johnpye/sunpos.a4c")
127
128 def testsunpos5(self):
129 self._run('example_1_7_3',"QRSlv","johnpye/sunpos.a4c")
130
131 def testsunpos6(self):
132 self._run('example_1_8_1',"QRSlv","johnpye/sunpos.a4c")
133
134 def testinstanceas(self):
135 M = self._run('example_1_6_1',"QRSlv","johnpye/sunpos.a4c")
136 self.assertAlmostEqual( float(M.t_solar), M.t_solar.as("s"))
137 self.assertAlmostEqual( float(M.t_solar)/3600, M.t_solar.as("h"))
138
139
140 class TestCMSlv(AscendSelfTester):
141 def testheatex(self):
142 self._run('heatex',"CMSlv","heatex.a4c")
143 def testphaseeq(self):
144 self._run('phaseq',"CMSlv","phaseq.a4c")
145 def testpipeline(self):
146 self._run('pipeline',"CMSlv","pipeline.a4c"
147 ,{'infinity':3.2e9}
148 )
149 def testrachford(self):
150 self._run('rachford',"CMSlv","rachford.a4c")
151 def testlinmassbal(self):
152 self._run('linmassbal',"CMSlv","linmassbal.a4c")
153
154
155 class TestMatrix(AscendSelfTester):
156 def testlog10(self):
157 M = self._run('testlog10')
158 print M.getMatrix().write(sys.stderr,"mmio")
159
160
161 class TestIntegrator(Ascend):
162
163 def testListIntegrators(self):
164 I = ascpy.Integrator.getEngines()
165 s1 = sorted([str(i) for i in I.values()])
166 s2 = sorted(['IDA','LSODE','AWW'])
167 assert s1==s2
168
169 # this routine is reused by both testIDA and testLSODE
170 def _testIntegrator(self,integratorname):
171 self.L.load('johnpye/shm.a4c')
172 M = self.L.findType('shm').getSimulation('sim')
173 M.setSolver(ascpy.Solver('QRSlv'))
174 P = M.getParameters()
175 M.setParameter('feastol',1e-12)
176 print M.getChildren()
177 assert float(M.x) == 10.0
178 assert float(M.v) == 0.0
179 t_end = math.pi
180
181 I = ascpy.Integrator(M)
182 I.setReporter(ascpy.IntegratorReporterNull(I))
183 I.setEngine(integratorname);
184 I.setLinearTimesteps(ascpy.Units("s"), 0.0, t_end, 100);
185 I.setMinSubStep(0.0001); # these limits are required by IDA at present (numeric diff)
186 I.setMaxSubStep(0.1);
187 I.setInitialSubStep(0.001);
188 I.setMaxSubSteps(200);
189 if(integratorname=='IDA'):
190 I.setParameter('autodiff',False)
191 for p in M.getParameters():
192 print p.getName(),"=",p.getValue()
193 I.analyse();
194 I.solve();
195 print "At end of simulation,"
196 print "x = %f" % M.x
197 print "v = %f" % M.v
198 assert abs(float(M.x) + 10) < 1e-2
199 assert abs(float(M.v)) < 1e-2
200 assert I.getNumObservedVars() == 3
201
202 def testInvalidIntegrator(self):
203 self.L.load('johnpye/shm.a4c')
204 M = self.L.findType('shm').getSimulation('sim')
205 M.setSolver(ascpy.Solver('QRSlv'))
206 I = ascpy.Integrator(M)
207 try:
208 I.setEngine('___NONEXISTENT____')
209 except RuntimeError:
210 return
211 self.fail("setEngine did not raise error!")
212
213 def testLSODE(self):
214 self._testIntegrator('LSODE')
215
216 def testIDA(self):
217 self._testIntegrator('IDA')
218
219 def testparameters(self):
220 self.L.load('johnpye/shm.a4c')
221 M = self.L.findType('shm').getSimulation('sim')
222 M.build()
223 I = ascpy.Integrator(M)
224 I.setEngine('IDA')
225 P = I.getParameters()
226 for p in P:
227 print p.getName(),"=",p.getValue()
228 assert len(P)==11
229 assert P[0].isStr()
230 assert P[0].getName()=="linsolver"
231 assert P[0].getValue()=='SPGMR'
232 assert P[2].getName()=="autodiff"
233 assert P[2].getValue()==True
234 assert P[7].getName()=="atolvect"
235 assert P[7].getBoolValue() == True
236 P[2].setBoolValue(False)
237 assert P[2].getBoolValue()==False
238 I.setParameters(P)
239 assert I.getParameterValue('autodiff')==False
240 I.setParameter('autodiff',True)
241 try:
242 v = I.getParameterValue('nonexist')
243 except KeyError:
244 pass
245 else:
246 self.fail('Failed to trip invalid Integrator parameter')
247
248 class TestLSODE(Ascend):
249
250 def testzill(self):
251 self.L.load('johnpye/zill.a4c')
252 T = self.L.findType('zill')
253 M = T.getSimulation('sim')
254 M.setSolver(ascpy.Solver('QRSlv'))
255 I = ascpy.Integrator(M)
256 I.setEngine('LSODE')
257 I.setMinSubStep(1e-7)
258 I.setMaxSubStep(0.001)
259 I.setMaxSubSteps(10000)
260 I.setReporter(ascpy.IntegratorReporterConsole(I))
261 I.setLinearTimesteps(ascpy.Units(), 1.0, 1.5, 5)
262 I.analyse()
263 I.solve()
264 M.run(T.getMethod('self_test'))
265
266 def testnewton(self):
267 sys.stderr.write("STARTING TESTNEWTON\n")
268 self.L.load('johnpye/newton.a4c')
269 T = self.L.findType('newton')
270 M = T.getSimulation('sim')
271 M.solve(ascpy.Solver("QRSlv"),ascpy.SolverReporter())
272 I = ascpy.Integrator(M)
273 I.setEngine('LSODE')
274 I.setParameter('rtolvect',False)
275 I.setParameter('rtol',1e-7)
276 I.setParameter('atolvect',False)
277 I.setParameter('atol',1e-7)
278 I.setMinSubStep(1e-7)
279 I.setMaxSubStep(0.001)
280 I.setMaxSubSteps(10000)
281
282 I.setReporter(ascpy.IntegratorReporterConsole(I))
283 I.setLinearTimesteps(ascpy.Units("s"), 0, 2*float(M.v)/float(M.g), 2)
284 I.analyse()
285 I.solve()
286 print "At end of simulation,"
287 print "x = %f" % M.x
288 print "v = %f" % M.v
289 M.run(T.getMethod('self_test'))
290
291 def testlotka(self):
292 self.L.load('johnpye/lotka.a4c')
293 M = self.L.findType('lotka').getSimulation('sim')
294 M.setSolver(ascpy.Solver("QRSlv"))
295 I = ascpy.Integrator(M)
296 I.setEngine('LSODE')
297 I.setReporter(ascpy.IntegratorReporterConsole(I))
298 I.setLinearTimesteps(ascpy.Units("s"), 0, 200, 5)
299 I.analyse()
300 print "Number of vars = %d" % I.getNumVars()
301 assert I.getNumVars()==2
302 I.solve()
303 assert I.getNumObservedVars() == 3;
304 assert abs(M.R - 832) < 1.0
305 assert abs(M.F - 21.36) < 0.1
306
307 #-------------------------------------------------------------------------------
308 # Testing of a external blackbox functions
309
310 class TestBlackBox(AscendSelfTester):
311 def testparsefail0(self):
312 try:
313 self.L.load('test/blackbox/parsefail0.a4c')
314 self.fail("parsefail0 should not have loaded without errors")
315 except:
316 pass
317
318 def testparsefail1(self):
319 try:
320 self.L.load('test/blackbox/parsefail1.a4c')
321 self.fail("parsefail1 should not have loaded without errors")
322 except:
323 pass
324
325 def testparsefail2(self):
326 try:
327 self.L.load('test/blackbox/parsefail2.a4c')
328 self.fail("parsefail2 should not have loaded without errors")
329 except:
330 pass
331
332 def testparsefail3(self):
333 try:
334 self.L.load('test/blackbox/parsefail3.a4c')
335 self.fail("parsefail3 should not have loaded without errors")
336 except:
337 pass
338
339 def testparsefail4(self):
340 try:
341 self.L.load('test/blackbox/parsefail4.a4c')
342 self.fail("parsefail4 should not have loaded")
343 except:
344 pass
345
346 def testfail1(self):
347 """Mismatched arg counts check-- tests bbox, not ascend."""
348 self.L.load('test/blackbox/fail1.a4c')
349 try:
350 M = self.L.findType('fail1').getSimulation('sim')
351 self.fail("expected exception was not raised")
352 except RuntimeError,e:
353 print "Caught exception '%s', assumed ok" % e
354
355 def testfail2(self):
356 """Incorrect data arg check -- tests bbox, not ascend"""
357 self.L.load('test/blackbox/fail2.a4c')
358 try:
359 M = self.L.findType('fail2').getSimulation('sim')
360 self.fail("expected exception was not raised")
361 except RuntimeError,e:
362 print "Caught exception '%s', assumed ok (should mention errors during instantiation)" % e
363
364 def testpass1(self):
365 """simple single bbox forward solve"""
366 M = self._run('pass1',filename='test/blackbox/pass.a4c')
367
368 def testpass2(self):
369 """simple single bbox reverse solve"""
370 M = self._run('pass2',filename='test/blackbox/pass.a4c')
371
372 def testpass3(self):
373 """simple double bbox solve"""
374 M = self._run('pass3',filename='test/blackbox/pass3.a4c')
375
376 def testpass4(self):
377 """simple double bbox reverse solve"""
378 M = self._run('pass4',filename='test/blackbox/pass3.a4c')
379
380 def testpass5(self):
381 M = self._run('pass5',filename='test/blackbox/pass5.a4c')
382
383 def testpass6(self):
384 M = self._run('pass6',filename='test/blackbox/pass5.a4c')
385
386 def testpass7(self):
387 M = self._run('pass7',filename='test/blackbox/passmerge.a4c')
388
389 def testpass8(self):
390 M = self._run('pass8',filename='test/blackbox/passmerge.a4c')
391
392 def testpass9(self):
393 M = self._run('pass9',filename='test/blackbox/passmerge.a4c')
394
395 def testpass10(self):
396 M = self._run('pass10',filename='test/blackbox/passmerge.a4c')
397
398 def testpass11(self):
399 M = self._run('pass11',filename='test/blackbox/passmerge.a4c')
400
401 def testpass12(self):
402 M = self._run('pass12',filename='test/blackbox/passmerge.a4c')
403
404 # this test doesn't work: 'system is inconsistent' -- and structurally singular
405 # def testpass13(self):
406 # """cross-merged input/output solve"""
407 # M = self._run('pass13',filename='test/blackbox/passmerge.a4c')
408
409 def testpass14(self):
410 """cross-merged input/output reverse solve"""
411 M = self._run('pass14',filename='test/blackbox/passmerge.a4c')
412
413 def testpass20(self):
414 M = self._run('pass20',filename='test/blackbox/passarray.a4c')
415
416 def testparsefail21(self):
417 """dense array of black boxes wrong syntax"""
418 try:
419 self.L.load('test/blackbox/parsefail21.a4c')
420 self.fail("parsefail21 should not have loaded without errors")
421 except:
422 pass
423
424 def testpass22(self):
425 M = self._run('pass22',filename='test/blackbox/passarray.a4c')
426
427 def testpass23(self):
428 M = self._run('pass23',filename='test/blackbox/passarray.a4c')
429
430 def testpass61(self):
431 M = self._run('pass61',filename='test/blackbox/reinstantiate.a4c')
432
433 def testpass62(self):
434 M = self._run('pass62',filename='test/blackbox/reinstantiate.a4c')
435
436 def testpass64(self):
437 M = self._run('pass64',filename='test/blackbox/reinstantiate.a4c')
438
439 def testpass65(self):
440 M = self._run('pass65',filename='test/blackbox/reinstantiate.a4c')
441
442 def testpass66(self):
443 M = self._run('pass66',filename='test/blackbox/reinstantiate.a4c')
444
445 def testpass67(self):
446 M = self._run('pass67',filename='test/blackbox/reinstantiate.a4c')
447
448 class TestExtFn(AscendSelfTester):
449 def testextfntest(self):
450 M = self._run('extfntest',filename='johnpye/extfn/extfntest.a4c')
451 self.assertAlmostEqual(M.y, 2);
452 self.assertAlmostEqual(M.x, 1);
453 self.assertAlmostEqual(M.y, M.x + 1);
454
455 def testextrelfor(self):
456 M = self._run('extrelfor',filename='johnpye/extfn/extrelfor.a4c')
457
458 ## @TODO fix bug with badly-named bbox rel in a loop (Ben, maybe)
459 # def testextrelforbadnaming(self):
460 # self.L.load('johnpye/extfn/extrelforbadnaming.a4c')
461 # T = self.L.findType('extrelfor')
462 # M = T.getSimulation('sim')
463 # M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
464 # print "x[1] = %f" % M.x[1]
465 # print "x[2] = %f" % M.x[2]
466 # print "x[3] = %f" % M.x[3]
467 # print "x[4] = %f" % M.x[4]
468 # print "x[5] = %f" % M.x[5]
469 # M.run(T.getMethod('self_test'))
470
471 def testextrelrepeat(self):
472 M = self._run('extrelrepeat',filename='johnpye/extfn/extrelrepeat.a4c')
473
474 #-------------------------------------------------------------------------------
475 # Testing of Sensitivity module
476
477 class TestSensitivity(AscendSelfTester):
478 def test1(self):
479 self.L.load('sensitivity_test.a4c')
480 T = self.L.findType('sensitivity_test')
481 M = T.getSimulation('sim',False)
482 M.run(T.getMethod('on_load'))
483 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
484 M.run(T.getMethod('analyse'))
485 M.run(T.getMethod('self_test'))
486
487 # def testall(self):
488 # self.L.load('sensitivity_test.a4c')
489 # T = self.L.findType('sensitivity_test_all')
490 # M = T.getSimulation('sim',False)
491 # M.run(T.getMethod('on_load'))
492 # M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
493 # M.run(T.getMethod('analyse'))
494 # M.run(T.getMethod('self_test'))
495 # CAUSES CRASH
496
497 #-------------------------------------------------------------------------------
498 # Testing of a ExtPy - external python methods
499
500 class TestExtPy(AscendSelfTester):
501 def test1(self):
502 self.L.load('johnpye/extpy/extpytest.a4c')
503 T = self.L.findType('extpytest')
504 M = T.getSimulation('sim')
505 M.run(T.getMethod('self_test'))
506
507 def test2(self):
508 self.L.load('johnpye/extpy/extpytest.a4c')
509 T = self.L.findType('extpytest')
510 M = T.getSimulation('sim')
511 M.run(T.getMethod('pythonthing'))
512 M.run(T.getMethod('pythonthing'))
513 M.run(T.getMethod('pythonthing'))
514 M.run(T.getMethod('pythonthing'))
515 M.run(T.getMethod('pythonthing'))
516 M.run(T.getMethod('pythonthing'))
517 M.run(T.getMethod('pythonthing'))
518 # causes crash!
519
520 #-------------------------------------------------------------------------------
521 # Testing of saturated steam properties library (iapwssatprops.a4c)
522
523 class TestSteam(AscendSelfTester):
524
525 def testiapwssatprops1(self):
526 M = self._run('testiapwssatprops1',filename='steam/iapwssatprops.a4c')
527 def testiapwssatprops2(self):
528 M = self._run('testiapwssatprops2',filename='steam/iapwssatprops.a4c')
529 def testiapwssatprops3(self):
530 M = self._run('testiapwssatprops3',filename='steam/iapwssatprops.a4c')
531
532 # test the stream model basically works
533 def testsatsteamstream(self):
534 M = self._run('satsteamstream',filename='steam/satsteamstream.a4c')
535
536 # test that we can solve in terms of various (rho,u)
537 def testsatuv(self):
538 self.L.load('steam/iapwssat.a4c')
539 T = self.L.findType('testiapwssatuv')
540 M = T.getSimulation('sim',False)
541 M.run(T.getMethod('on_load'))
542 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
543 print "p = %f bar" % M.p.as('bar');
544 print "T = %f C" % (M.T.as('K') - 273.15);
545 print "x = %f" % M.x;
546 M.run(T.getMethod('self_test'))
547 M.run(T.getMethod('values2'))
548 # M.v.setRealValueWithUnits(1.0/450,"m^3/kg");
549 # M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
550 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
551 print "p = %f bar" % M.p.as('bar');
552 print "T = %f C" % (M.T.as('K') - 273.15);
553 print "x = %f" % M.x;
554 M.run(T.getMethod('self_test2'))
555
556
557 ## @TODO fix error capture from bounds checking during initialisation
558 # def testiapwssat1(self):
559 # M = self._run('testiapwssat1',filename='steam/iapwssat.a4c')
560
561 def testdsgsat(self):
562 self.L.load('steam/dsgsat3.a4c')
563 T = self.L.findType('dsgsat3')
564 M = T.getSimulation('sim',False)
565 M.run(T.getMethod('on_load'))
566 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
567 self.assertAlmostEqual(M.dTw_dt[2],0.0);
568 M.run(T.getMethod('configure_dynamic'))
569 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
570 return M
571
572 def testdsgsatrepeat(self):
573 self.L.load('steam/dsgsat3.a4c')
574 T = self.L.findType('dsgsat3')
575 M = T.getSimulation('sim',False)
576 M.run(T.getMethod('on_load'))
577 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
578 M.run(T.getMethod('on_load'))
579 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
580 M.run(T.getMethod('on_load'))
581 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
582
583 def teststeadylsode(self):
584 "test that steady conditions are stable with LSODE"
585 M = self.testdsgsat()
586 #M.qdot_s.setRealValueWithUnits(1000,"W/m")
587 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
588 #M.setParameter('
589 I = ascpy.Integrator(M)
590 I.setEngine('LSODE')
591 I.setReporter(ascpy.IntegratorReporterConsole(I))
592 I.setLinearTimesteps(ascpy.Units("s"), 0, 5, 1)
593 I.analyse()
594 I.solve()
595
596 # DOESN'T CONVERGE
597 # def testpeturblsode(self):
598 # "test that steady conditions are stable with LSODE"
599 # M = self.testdsgsat()
600 # # here is the peturbation...
601 # M.qdot_s.setRealValueWithUnits(1000,"W/m")
602 # M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
603 # I = ascpy.Integrator(M)
604 # I.setEngine('LSODE')
605 # I.setReporter(ascpy.IntegratorReporterConsole(I))
606 # I.setLinearTimesteps(ascpy.Units("s"), 0, 5, 1)
607 # I.analyse()
608 # I.solve()
609
610 def teststeadyida(self):
611 M = self.testdsgsat()
612 self.assertAlmostEqual(M.dTw_dt[2],0.0)
613 Tw1 = float(M.T_w[2])
614 T = self.L.findType('dsgsat3')
615 M.run(T.getMethod('free_states'))
616 I = ascpy.Integrator(M)
617 I.setEngine('IDA')
618 I.setParameter('linsolver','DENSE')
619 I.setParameter('safeeval',True)
620 I.setParameter('rtol',1e-5)
621 I.setInitialSubStep(0.01)
622 I.setMaxSubSteps(100)
623 I.setReporter(ascpy.IntegratorReporterConsole(I))
624 I.setLinearTimesteps(ascpy.Units("s"), 0, 3600, 5)
625 try:
626 I.analyse()
627 except Exception,e:
628 print "ERROR: %s" % e
629 I.writeDebug(sys.stdout)
630
631 I.solve()
632 self.assertAlmostEqual(float(M.T_w[2]),Tw1)
633 M.qdot_s.setRealValueWithUnits(1000,"W/m")
634 self.assertAlmostEqual(M.qdot_s.as("W/m"),1000)
635 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
636 print "dTw/dt = %f" % M.dTw_dt[2]
637 self.assertNotAlmostEqual(M.dTw_dt[2],0.0)
638
639 def testpeturbida(self):
640 M = self.testdsgsat()
641 self.assertAlmostEqual(M.dTw_dt[2],0.0)
642 T = self.L.findType('dsgsat3')
643 M.run(T.getMethod('free_states'))
644 # here is the peturbation...
645 M.qdot_s.setRealValueWithUnits(2000,"W/m")
646 # IDA has its own initial conditions solver, so no need to call QRSlv here
647 I = ascpy.Integrator(M)
648 I.setEngine('IDA')
649 I.setParameter('linsolver','DENSE')
650 I.setParameter('safeeval',True)
651 I.setParameter('rtol',1e-5)
652 I.setParameter('atolvect',False)
653 I.setParameter('atol',1e-5)
654 I.setInitialSubStep(0.0001)
655 I.setReporter(ascpy.IntegratorReporterConsole(I))
656 I.setLogTimesteps(ascpy.Units("s"), 1, 5, 20)
657 I.analyse()
658 I.solve()
659
660 #-------------------------------------------------------------------------------
661 # Testing of freesteam external steam properties functions
662
663 with_freesteam = True
664 try:
665 # we assume that if the freesteam python module is installed, the ASCEND
666 # external library will also be.
667 import freesteam
668 have_freesteam = True
669 except ImportError,e:
670 have_freesteam = False
671
672 if with_freesteam and have_freesteam:
673 class TestFreesteam(AscendSelfTester):
674 # def testfreesteamtest(self):
675 # """run the self-test cases bundled with freesteam"""
676 # self._run('testfreesteam',filename='testfreesteam.a4c')
677
678 def testload(self):
679 """check that we can load 'thermalequilibrium2' (IMPORT "freesteam", etc)"""
680 self.L.load('johnpye/thermalequilibrium2.a4c')
681
682 def testinstantiate(self):
683 """load an instantiate 'thermalequilibrium2'"""
684 self.testload()
685 M = self.L.findType('thermalequilibrium2').getSimulation('sim')
686 return M
687
688 def testintegrate(self):
689 """integrate transfer of heat from one mass of water/steam to another
690 according to Newton's law of cooling"""
691 M = self.testinstantiate()
692 M.setSolver(ascpy.Solver("QRSlv"))
693 I = ascpy.Integrator(M)
694 I.setEngine('LSODE')
695 I.setReporter(ascpy.IntegratorReporterConsole(I))
696 I.setLinearTimesteps(ascpy.Units("s"), 0, 3000, 30)
697 I.setMinSubStep(0.01)
698 I.setInitialSubStep(1)
699 I.analyse()
700 print "Number of vars = %d" % I.getNumVars()
701 assert I.getNumVars()==2
702 I.solve()
703 assert I.getNumObservedVars() == 3;
704 print "S[1].T = %f K" % M.S[1].T
705 print "S[2].T = %f K" % M.S[2].T
706 print "Q = %f W" % M.Q
707 self.assertAlmostEqual(float(M.S[1].T),506.77225109,4);
708 self.assertAlmostEqual(float(M.S[2].T),511.605173967,5);
709 self.assertAlmostEqual(float(M.Q),-48.32922877329,3);
710 self.assertAlmostEqual(float(M.t),3000);
711 print "Note that the above values have not been verified analytically"
712
713 def testcollapsingcan2(self):
714 """ solve the collapsing can model using IAPWS-IF97 steam props """
715 M = self._run("collapsingcan2",filename="collapsingcan2.a4c");
716
717 #-------------------------------------------------------------------------------
718 # Testing of IDA's analysis module
719
720 class TestIDA(Ascend):
721 def _run(self,filen,modeln=""):
722 self.L.load('test/ida/%s.a4c' % filen)
723 T = self.L.findType('%s%s' % (filen,modeln))
724 M = T.getSimulation('sim')
725 M.build()
726 I = ascpy.Integrator(M)
727 I.setEngine('IDA')
728 I.analyse()
729 return M;
730
731 def _runfail(self,filen,n,msg="failed"):
732 try:
733 self._run(filen,'fail%d' % n)
734 except Exception,e:
735 print "(EXPECTED) ERROR: %s" % e
736 return
737 self.fail(msg)
738
739 def testsinglederiv(self):
740 self._run('singlederiv')
741
742 def testsinglederivfail1(self):
743 self._runfail('singlederiv',1
744 ,"t.ode_id=-1 did not trigger error")
745
746 def testsinglederivfail2(self):
747 self._runfail('singlederiv',2
748 ,"dy_dt.ode_id=2 did not trigger error")
749
750 def testsinglederivfail3(self):
751 self._runfail('singlederiv',3
752 ,"dy_dt.ode_type=3 did not trigger error")
753
754 def testsinglederivfail4(self):
755 self._runfail('singlederiv',4
756 ,"duplicate ode_type=1 did not trigger error")
757
758 def testsinglederivfail5(self):
759 self._runfail('singlederiv',5
760 ,"duplicate ode_type=1 did not trigger error")
761
762 def testsinglederivfail6(self):
763 self._runfail('singlederiv',6
764 ,"duplicate ode_type=1 did not trigger error")
765
766 def testtwoderiv(self):
767 self._run('twoderiv')
768
769 def testtwoderivfail1(self):
770 self._runfail('twoderiv',1)
771
772 def testtwoderivfail2(self):
773 self._runfail('twoderiv',2)
774
775 def testtwoderivfail3(self):
776 self._runfail('twoderiv',3)
777 def testtwoderivfail4(self):
778 self._runfail('twoderiv',4)
779 def testtwoderivfail5(self):
780 self._runfail('twoderiv',5)
781
782 def testnoderivs(self):
783 self._runfail('noderivs',1)
784
785 def testnoindeps(self):
786 self._runfail('indeps',1)
787
788 def testtwoindeps(self):
789 self._runfail('indeps',2)
790
791 def testfixedvars(self):
792 self._run('fixedvars')
793
794 def testfixedvars1(self):
795 self._run('fixedvars',1)
796
797 def testfixedvars2(self):
798 self._run('fixedvars',2)
799
800 def testfixedvars3(self):
801 self._run('fixedvars',3)
802
803 def testincidence(self):
804 self._run('incidence')
805
806 def testincidence1(self):
807 self._run('incidence',1)
808 def testincidence2(self):
809 self._run('incidence',2)
810 def testincidence3(self):
811 M = self._run('incidence',3)
812
813 def testincidence4(self):
814 self._run('incidence',4)
815 def testincidencefail5(self):
816 self._runfail('incidence',5)
817
818 # doesn't work yet:
819 # def testincidence5(self):
820 # self._run('incidence',5)
821
822
823 #-------------------------------------------------------------------------------
824 # Testing of IDA models using DENSE linear solver
825
826 class TestIDADENSE(Ascend):
827 """IDA DAE integrator, DENSE linear solver"""
828
829 def testlotka(self):
830 self.L.load('johnpye/lotka.a4c')
831 M = self.L.findType('lotka').getSimulation('sim')
832 M.setSolver(ascpy.Solver("QRSlv"))
833 I = ascpy.Integrator(M)
834 I.setEngine('IDA')
835 I.setReporter(ascpy.IntegratorReporterConsole(I))
836 I.setLinearTimesteps(ascpy.Units("s"), 0, 200, 5);
837 I.setParameter('linsolver','DENSE')
838 I.setParameter('rtol',1e-8);
839 I.analyse()
840 assert I.getNumVars()==2
841 assert abs(M.R - 1000) < 1e-300
842 I.solve()
843 assert I.getNumObservedVars() == 3
844 assert abs(M.R - 832) < 1.0
845 assert abs(M.F - 21.36) < 0.1
846
847 def testdenx(self):
848 print "-----------------------------====="
849 self.L.load('johnpye/idadenx.a4c')
850 M = self.L.findType('idadenx').getSimulation('sim')
851 M.setSolver(ascpy.Solver("QRSlv"))
852 I = ascpy.Integrator(M)
853 I.setEngine('IDA')
854 I.setParameter('calcic','YA_YDP')
855 I.setParameter('linsolver','DENSE')
856 I.setParameter('safeeval',True)
857 I.setReporter(ascpy.IntegratorReporterConsole(I))
858 I.setLogTimesteps(ascpy.Units("s"), 0.4, 4e10, 11)
859 I.setMaxSubStep(0);
860 I.setInitialSubStep(0)
861 I.setMaxSubSteps(0)
862 I.setParameter('autodiff',True)
863 I.analyse()
864 I.solve()
865 assert abs(float(M.y1) - 5.1091e-08) < 2e-9
866 assert abs(float(M.y2) - 2.0437e-13) < 2e-14
867 assert abs(float(M.y3) - 1.0) < 1e-5
868
869 def testhires(self):
870 self.L.load('test/hires.a4c')
871 T = self.L.findType('hires')
872 M = T.getSimulation('sim')
873 M.setSolver(ascpy.Solver('QRSlv'))
874 I = ascpy.Integrator(M)
875 I.setEngine('IDA')
876 I.setParameter('linsolver','DENSE')
877 I.setParameter('rtol',1.1e-15)
878 I.setParameter('atolvect',0)
879 I.setParameter('atol',1.1e-15)
880 I.setReporter(ascpy.IntegratorReporterConsole(I))
881 I.setLogTimesteps(ascpy.Units(""), 1, 321.8122, 5)
882 I.setInitialSubStep(1e-5)
883 I.setMaxSubSteps(10000)
884 I.analyse()
885 I.solve()
886 for i in range(8):
887 print "y[%d] = %.20g" % (i+1, M.y[i+1])
888 M.run(T.getMethod('self_test'))
889
890 def testchemakzo(self):
891 self.L.load('test/chemakzo.a4c')
892 T = self.L.findType('chemakzo')
893 M = T.getSimulation('sim')
894 M.setSolver(ascpy.Solver('QRSlv'))
895 I = ascpy.Integrator(M)
896 I.setEngine('IDA')
897 I.setParameter('linsolver','DENSE')
898 I.setParameter('rtol',1e-15)
899 I.setParameter('atolvect',0)
900 I.setParameter('atol',1e-15)
901 I.setReporter(ascpy.IntegratorReporterConsole(I))
902 I.setLinearTimesteps(ascpy.Units("s"), 1, 180, 5)
903 I.setInitialSubStep(1e-13)
904 I.setMaxSubSteps(10000)
905 I.analyse()
906 I.solve()
907 for i in range(6):
908 print "y[%d] = %.20g" % (i+1, M.y[i+1])
909 M.run(T.getMethod('self_test'))
910
911 def testtransamp(self):
912 self.L.load('test/transamp.a4c')
913 T = self.L.findType('transamp')
914 M = T.getSimulation('sim')
915 M.setSolver(ascpy.Solver('QRSlv'))
916 I = ascpy.Integrator(M)
917 I.setEngine('IDA')
918 I.setParameter('linsolver','DENSE')
919 I.setParameter('rtol',1e-7)
920 I.setParameter('atolvect',0)
921 I.setParameter('atol',1e-7)
922 I.setReporter(ascpy.IntegratorReporterConsole(I))
923 I.setLinearTimesteps(ascpy.Units("s"), 0.05, 0.2, 20)
924 I.setInitialSubStep(0.00001)
925 I.setMaxSubSteps(10000)
926 I.analyse()
927 I.solve()
928 for i in range(6):
929 print "y[%d] = %.20g" % (i+1, M.y[i+1])
930 M.run(T.getMethod('self_test'))
931
932 # MODEL FAILS ANALYSIS: we need to add support for non-incident differential vars
933 # def testpollution(self):
934 # self.L.load('test/pollution.a4c')
935 # T = self.L.findType('pollution')
936 # M = T.getSimulation('sim')
937 # M.setSolver(ascpy.Solver('QRSlv'))
938 # I = ascpy.Integrator(M)
939 # I.setEngine('IDA')
940 # I.setParameter('linsolver','DENSE')
941 # I.setParameter('rtol',1.1e-15)
942 # I.setParameter('atolvect',0)
943 # I.setParameter('atol',1.1e-15)
944 # I.setReporter(ascpy.IntegratorReporterConsole(I))
945 # I.setLogTimesteps(ascpy.Units("s"), 1, 60, 5)
946 # I.setInitialSubStep(1e-5)
947 # I.setMaxSubSteps(10000)
948 # I.analyse()
949 # I.solve()
950 # for i in range(20):
951 # print "y[%d] = %.20g" % (i+1, M.y[i+1])
952 # M.run(T.getMethod('self_test'))
953
954 ## @TODO fails during IDACalcIC (model too big?)
955 # def testkryx(self):
956 # self.L.load('johnpye/idakryx.a4c')
957 # ascpy.getCompiler().setUseRelationSharing(False)
958 # M = self.L.findType('idakryx').getSimulation('sim')
959 # M.setSolver(ascpy.Solver('QRSlv'))
960 # M.build()
961 # I = ascpy.Integrator(M)
962 # I.setEngine('IDA')
963 # I.setReporter(ascpy.IntegratorReporterConsole(I))
964 # I.setParameter('linsolver','DENSE')
965 # I.setParameter('maxl',8)
966 # I.setParameter('gsmodified',False)
967 # I.setParameter('autodiff',True)
968 # I.setParameter('rtol',0)
969 # I.setParameter('atol',1e-3);
970 # I.setParameter('atolvect',False)
971 # I.setParameter('calcic','YA_YDP')
972 # I.analyse()
973 # I.setLogTimesteps(ascpy.Units("s"), 0.01, 10.24, 11)
974 # I.solve()
975 # assert abs(M.u[2][2].getValue()) < 1e-5
976
977 #-------------------------------------------------------------------------------
978 # Testing of IDA models using SPGMR linear solver (Krylov)
979
980 # these tests are disabled until SPGMR preconditioning has been implemented
981 class TestIDASPGMR:#(Ascend):
982 def testlotka(self):
983 self.L.load('johnpye/lotka.a4c')
984 M = self.L.findType('lotka').getSimulation('sim')
985 M.setSolver(ascpy.Solver("QRSlv"))
986 I = ascpy.Integrator(M)
987 I.setEngine('IDA')
988 I.setReporter(ascpy.IntegratorReporterConsole(I))
989 I.setLinearTimesteps(ascpy.Units("s"), 0, 200, 5)
990 I.setParameter('rtol',1e-8)
991 I.analyse()
992 assert I.getNumVars()==2
993 assert abs(M.R - 1000) < 1e-300
994 I.solve()
995 assert I.getNumObservedVars() == 3
996 assert abs(M.R - 832) < 1.0
997 assert abs(M.F - 21.36) < 0.1
998
999
1000 def testkryx(self):
1001 self.L.load('johnpye/idakryx.a4c')
1002 M = self.L.findType('idakryx').getSimulation('sim')
1003 M.build()
1004 I = ascpy.Integrator(M)
1005 I.setEngine('IDA')
1006 I.setReporter(ascpy.IntegratorReporterConsole(I))
1007 I.setParameter('linsolver','SPGMR')
1008 I.setParameter('prec','JACOBI')
1009 I.setParameter('maxl',8)
1010 I.setParameter('gsmodified',False)
1011 I.setParameter('autodiff',True)
1012 I.setParameter('gsmodified',True)
1013 I.setParameter('rtol',0)
1014 I.setParameter('atol',1e-3);
1015 I.setParameter('atolvect',False)
1016 I.setParameter('calcic','Y')
1017 I.analyse()
1018 I.setLogTimesteps(ascpy.Units("s"), 0.01, 10.24, 10);
1019 print M.udot[1][3]
1020 I.solve()
1021 assert 0
1022
1023 def testzill(self):
1024 self.L.load('johnpye/zill.a4c')
1025 T = self.L.findType('zill')
1026 M = T.getSimulation('sim')
1027 M.setSolver(ascpy.Solver('QRSlv'))
1028 I = ascpy.Integrator(M)
1029 I.setEngine('IDA')
1030 I.setParameter('safeeval',False)
1031 I.setMinSubStep(1e-7)
1032 I.setMaxSubStep(0.001)
1033 I.setMaxSubSteps(10000)
1034 I.setReporter(ascpy.IntegratorReporterConsole(I))
1035 I.setLinearTimesteps(ascpy.Units(), 1.0, 1.5, 5)
1036 I.analyse()
1037 I.solve()
1038 M.run(T.getMethod('self_test'))
1039
1040 def testdenxSPGMR(self):
1041 self.L.load('johnpye/idadenx.a4c')
1042 M = self.L.findType('idadenx').getSimulation('sim')
1043 M.setSolver(ascpy.Solver('QRSlv'))
1044 I = ascpy.Integrator(M)
1045 I.setEngine('IDA')
1046 I.setReporter(ascpy.IntegratorReporterConsole(I))
1047 I.setLogTimesteps(ascpy.Units("s"), 0.4, 4e10, 11)
1048 I.setMaxSubStep(0);
1049 I.setInitialSubStep(0);
1050 I.setMaxSubSteps(0);
1051 I.setParameter('autodiff',True)
1052 I.setParameter('linsolver','SPGMR')
1053 I.setParameter('gsmodified',False)
1054 I.setParameter('maxncf',10)
1055 I.analyse()
1056 I.solve()
1057 assert abs(float(M.y1) - 5.1091e-08) < 1e-10
1058 assert abs(float(M.y2) - 2.0437e-13) < 1e-15
1059 assert abs(float(M.y3) - 1.0) < 1e-5
1060
1061 # move code above down here if you want to temporarily avoid testing it
1062 class NotToBeTested:
1063 def nothing(self):
1064 pass
1065
1066 def testnewton(self):
1067 sys.stderr.write("STARTING TESTNEWTON\n")
1068 self.L.load('johnpye/newton.a4c')
1069 T = self.L.findType('newton')
1070 M = T.getSimulation('sim')
1071 M.solve(ascpy.Solver("QRSlv"),ascpy.SolverReporter())
1072 I = ascpy.Integrator(M)
1073 I.setEngine('IDA')
1074 I.setParameter('linsolver','DENSE')
1075 I.setParameter('safeeval',True)
1076 I.setParameter('rtol',1e-8)
1077 I.setMaxSubStep(0.001)
1078 I.setMaxSubSteps(10000)
1079
1080 I.setReporter(ascpy.IntegratorReporterConsole(I))
1081 I.setLinearTimesteps(ascpy.Units("s"), 0, 2*float(M.v)/float(M.g), 2)
1082 I.analyse()
1083 I.solve()
1084 print "At end of simulation,"
1085 print "x = %f" % M.x
1086 print "v = %f" % M.v
1087 M.run(T.getMethod('self_test'))
1088
1089 if __name__=='__main__':
1090 # a whole bag of tricks to make sure we get the necessary dirs in our ascend, python and ld path vars
1091 restart = 0
1092
1093 if platform.system()=="Windows":
1094 LD_LIBRARY_PATH="PATH"
1095 SEP = ";"
1096 else:
1097 LD_LIBRARY_PATH="LD_LIBRARY_PATH"
1098 SEP = ":"
1099
1100 freesteamdir = os.path.expanduser("~/freesteam/ascend")
1101 modeldirs = [os.path.abspath(os.path.join(sys.path[0],"models")),os.path.abspath(freesteamdir)]
1102 if not os.environ.get('ASCENDLIBRARY'):
1103 os.environ['ASCENDLIBRARY'] = SEP.join(modeldirs)
1104 restart = 1
1105 else:
1106 envmodelsdir = [os.path.abspath(i) for i in os.environ['ASCENDLIBRARY'].split(SEP)]
1107 for l in modeldirs:
1108 if l in envmodelsdir[len(modeldirs):]:
1109 envmodelsdir.remove(l)
1110 restart = 1
1111 for l in modeldirs:
1112 if l not in envmodelsdir:
1113 envmodelsdir.insert(0,l)
1114 restart = 1
1115 os.environ['ASCENDLIBRARY'] = SEP.join(envmodelsdir)
1116
1117 libdirs = ["pygtk","."]
1118 libdirs = [os.path.normpath(os.path.join(sys.path[0],l)) for l in libdirs]
1119 if not os.environ.get(LD_LIBRARY_PATH):
1120 os.environ[LD_LIBRARY_PATH]=SEP.join(libdirs)
1121 restart = 1
1122 else:
1123 envlibdirs = [os.path.normpath(i) for i in os.environ[LD_LIBRARY_PATH].split(SEP)]
1124 for l in libdirs:
1125 if l in envlibdirs[len(libdirs):]:
1126 envlibdirs.remove(l)
1127 restart = 1
1128 for l in libdirs:
1129 if l not in envlibdirs:
1130 envlibdirs.insert(0,l)
1131 restart = 1
1132 os.environ[LD_LIBRARY_PATH] = SEP.join(envlibdirs)
1133
1134 pypath = os.path.normpath(os.path.join(sys.path[0],"pygtk"))
1135 if not os.environ.get('PYTHONPATH'):
1136 os.environ['PYTHONPATH']=pypath
1137 else:
1138 envpypath = os.environ['PYTHONPATH'].split(SEP)
1139 if pypath not in envpypath:
1140 envpypath.insert(0,pypath)
1141 os.environ['PYTHONPATH']=SEP.join(envpypath)
1142 restart = 1
1143
1144 if restart:
1145 script = os.path.join(sys.path[0],"test.py")
1146 print "Restarting with..."
1147 print " export LD_LIBRARY_PATH=%s" % os.environ.get(LD_LIBRARY_PATH)
1148 print " export PYTHONPATH=%s" % os.environ.get('PYTHONPATH')
1149 print " export ASCENDLIBRARY=%s" % os.environ.get('ASCENDLIBRARY')
1150
1151 os.execvp("python",[script] + sys.argv)
1152
1153 import ascpy
1154
1155 try:
1156 import cunit
1157 except:
1158 pass
1159
1160 atexit.register(ascpy.shutdown)
1161 #suite = unittest.TestSuite()
1162 #suite = unittest.defaultTestLoader.loadTestsFromName('__main__')
1163 #unittest.TextTestRunner(verbosity=2).run(suite)
1164 unittest.main()

Properties

Name Value
svn:executable *

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