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

Contents of /trunk/test.py

Parent Directory Parent Directory | Revision Log Revision Log


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

Properties

Name Value
svn:executable *

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