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

Contents of /trunk/test.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2329 - (show annotations) (download) (as text)
Wed Dec 22 12:52:47 2010 UTC (6 years, 10 months ago) by jpye
File MIME type: text/x-python
File size: 50038 byte(s)
Suppressing some console output.
Added test case that current crashes ASCEND/IPOPT with a memory error.
Issue with library.cpp not containing any way to free libascend memory.
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 TestSystem(AscendSelfTester):
178
179 def testwritegraph(self):
180 M = self._run('testlog10')
181 M.run(self.L.findType('testlog10').getMethod('on_load'))
182 if platform.system!="Windows":
183 f = file('temp.png','wb')
184 M.write(f,"dot")
185 f.close()
186 else:
187 self.fail("not implemented on windows")
188
189 class TestSolver(AscendSelfTester):
190
191 def testlog10(self):
192 M = self._run('testlog10')
193
194 def testrootsofpoly(self):
195 self._run('roots_of_poly',filename="roots_of_poly.a4c")
196
197 def testcollapsingcan(self):
198 self._run('collapsingcan',filename="collapsingcan.a4c")
199
200 def testdistancecalc(self):
201 self._run('distance_calc',filename="distance_calc.a4c")
202
203 def testconopt(self):
204 self._run('conopttest',"CONOPT",filename="conopttest.a4c")
205
206 def testcmslv2(self):
207 self._run('testcmslv2',"CMSlv")
208
209 def testsunpos1(self):
210 self._run('example_1_6_1',"QRSlv","johnpye/sunpos.a4c")
211
212 def testsunpos2(self):
213 self._run('example_1_6_2',"QRSlv","johnpye/sunpos.a4c")
214
215 def testsunpos3(self):
216 self._run('example_1_7_1',"QRSlv","johnpye/sunpos.a4c")
217
218 def testsunpos4(self):
219 self._run('example_1_7_2',"QRSlv","johnpye/sunpos.a4c")
220
221 def testsunpos5(self):
222 self._run('example_1_7_3',"QRSlv","johnpye/sunpos.a4c")
223
224 def testsunpos6(self):
225 self._run('example_1_8_1',"QRSlv","johnpye/sunpos.a4c")
226
227 def testinstanceas(self):
228 M = self._run('example_1_6_1',"QRSlv","johnpye/sunpos.a4c")
229 self.assertAlmostEqual( float(M.t_solar), M.t_solar.to("s"))
230 self.assertAlmostEqual( float(M.t_solar)/3600, M.t_solar.to("h"))
231
232 def testrelinclude(self):
233 self.L.load('test/relinclude.a4c')
234 T = self.L.findType('relinclude')
235 M = T.getSimulation('sim')
236 M.eq1.setIncluded(True)
237 M.eq2.setIncluded(False)
238 M.eq3.setIncluded(False)
239 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
240 self.assertAlmostEqual( float(M.z), 2.0)
241 M.eq1.setIncluded(False)
242 M.eq2.setIncluded(True)
243 M.eq3.setIncluded(False)
244 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
245 self.assertAlmostEqual( float(M.z), 4.0)
246 M.eq1.setIncluded(False)
247 M.eq2.setIncluded(False)
248 M.eq3.setIncluded(True)
249 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
250 self.assertAlmostEqual( float(M.z), 4.61043629206)
251
252
253 class TestBinTokens(AscendSelfTester):
254
255 def test1(self):
256 ascpy.getCompiler().setBinaryCompilation(True)
257 self.L.load('johnpye/testlog10.a4c')
258 T = self.L.findType('testlog10')
259 M = T.getSimulation('sim')
260 M.build()
261 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
262
263 class TestLRSlv(AscendSelfTester):
264 def testonerel(self):
265 self._run('onerel',"LRSlv","test/lrslv/onerel.a4c")
266
267 # need to migrate to 'FIX boolvar', currently not supported...
268 # def testonerel(self):
269 # self._run('onerel',"LRSlv","test/lrslv/onerel.a4c")
270
271 def testsequencecrash(self):
272 try:
273 self._run('onerel',"LRSlv","test/lrslv/sequencecrash.a4c")
274 except:
275 pass
276 # it just has to not crash, that's all
277
278 def testsequence(self):
279 self._run('onerel',"LRSlv","test/lrslv/sequence.a4c")
280
281
282 class TestCMSlv(AscendSelfTester):
283 def testsonic(self):
284 M = self._run('sonic',"CMSlv","sonic.a4c")
285 assert(M.sonic_flow.getBoolValue())
286
287 # other side of boundary...
288 M.D.setRealValueWithUnits(4.,"cm")
289 T = self.L.findType('sonic')
290 M.solve(ascpy.Solver('CMSlv'),ascpy.SolverReporter())
291 M.run(T.getMethod('self_test'))
292 assert(not M.sonic_flow.getBoolValue())
293
294 def testheatex(self):
295 self._run('heatex',"CMSlv","heatex.a4c")
296 def testphaseeq(self):
297 self._run('phaseq',"CMSlv","phaseq.a4c")
298 def testpipeline(self):
299 self._run('pipeline',"CMSlv","pipeline.a4c"
300 ,{'infinity':3.2e9}
301 )
302 def testrachford(self):
303 self._run('rachford',"CMSlv","rachford.a4c")
304 def testlinmassbal(self):
305 self._run('linmassbal',"CMSlv","linmassbal.a4c")
306
307
308 class TestMatrix(AscendSelfTester):
309 def testlog10(self):
310 M = self._run('testlog10')
311 print "FETCHING MATRIX................."
312 X = M.getMatrix()
313 # this stuff crashes Windows because the FILE* structure used by Python is not the same
314 # as used by MinGW...
315 #print "GOT MATRIX"
316 #sys.stderr.flush()
317 #sys.stdout.flush()
318 #F = os.tmpfile()
319 #X.write(F.fileno,"mmio")
320 #F.seek(0)
321 #print F.read()
322
323 class TestIntegrator(Ascend):
324
325 def testListIntegrators(self):
326 I = ascpy.Integrator.getEngines()
327 s1 = sorted([str(i) for i in I])
328 assert 'IDA' in s1
329 assert 'LSODE' in s1
330
331 # this routine is reused by both testIDA and testLSODE
332 def _testIntegrator(self,integratorname):
333 self.L.load('johnpye/shm.a4c')
334 M = self.L.findType('shm').getSimulation('sim')
335 M.setSolver(ascpy.Solver('QRSlv'))
336 P = M.getParameters()
337 M.setParameter('feastol',1e-12)
338 print M.getChildren()
339 assert float(M.x) == 10.0
340 assert float(M.v) == 0.0
341 t_end = math.pi
342
343 I = ascpy.Integrator(M)
344 I.setReporter(ascpy.IntegratorReporterNull(I))
345 I.setEngine(integratorname);
346 I.setLinearTimesteps(ascpy.Units("s"), 0.0, t_end, 100);
347 I.setMinSubStep(0.0001); # these limits are required by IDA at present (numeric diff)
348 I.setMaxSubStep(0.1);
349 I.setInitialSubStep(0.001);
350 I.setMaxSubSteps(200);
351 if(integratorname=='IDA'):
352 I.setParameter('autodiff',False)
353 for p in M.getParameters():
354 print p.getName(),"=",p.getValue()
355 I.analyse();
356 I.solve();
357 print "At end of simulation,"
358 print "x = %f" % M.x
359 print "v = %f" % M.v
360 assert abs(float(M.x) + 10) < 1e-2
361 assert abs(float(M.v)) < 1e-2
362 assert I.getNumObservedVars() == 3
363
364 def testInvalidIntegrator(self):
365 self.L.load('johnpye/shm.a4c')
366 M = self.L.findType('shm').getSimulation('sim')
367 M.setSolver(ascpy.Solver('QRSlv'))
368 I = ascpy.Integrator(M)
369 try:
370 I.setEngine('___NONEXISTENT____')
371 except IndexError:
372 return
373 self.fail("setEngine did not raise error!")
374
375 def testLSODE(self):
376 self._testIntegrator('LSODE')
377
378 def testIDA(self):
379 self._testIntegrator('IDA')
380
381 def testparameters(self):
382 self.L.load('johnpye/shm.a4c')
383 M = self.L.findType('shm').getSimulation('sim')
384 M.build()
385 I = ascpy.Integrator(M)
386 I.setEngine('IDA')
387 P = I.getParameters()
388 for p in P:
389 print p.getName(),"=",p.getValue()
390 assert len(P)==12
391 assert P[0].isStr()
392 assert P[0].getName()=="linsolver"
393 assert P[0].getValue()=='DENSE'
394 assert P[2].getName()=="maxord"
395 assert P[3].getName()=="autodiff"
396 assert P[3].getValue()==True
397 assert P[8].getName()=="atolvect"
398 assert P[8].getBoolValue() == True
399 P[3].setBoolValue(False)
400 assert P[3].getBoolValue()==False
401 I.setParameters(P)
402 assert I.getParameterValue('autodiff')==False
403 I.setParameter('autodiff',True)
404 try:
405 v = I.getParameterValue('nonexist')
406 except KeyError:
407 pass
408 else:
409 self.fail('Failed to trip invalid Integrator parameter')
410
411 class TestLSODE(Ascend):
412
413 def testzill(self):
414 self.L.load('johnpye/zill.a4c')
415 T = self.L.findType('zill')
416 M = T.getSimulation('sim')
417 M.setSolver(ascpy.Solver('QRSlv'))
418 I = ascpy.Integrator(M)
419 I.setEngine('LSODE')
420 I.setMinSubStep(1e-7)
421 I.setMaxSubStep(0.001)
422 I.setMaxSubSteps(10000)
423 I.setReporter(ascpy.IntegratorReporterConsole(I))
424 I.setLinearTimesteps(ascpy.Units(), 1.0, 1.5, 5)
425 I.analyse()
426 I.solve()
427 M.run(T.getMethod('self_test'))
428
429 def testnewton(self):
430 sys.stderr.write("STARTING TESTNEWTON\n")
431 self.L.load('johnpye/newton.a4c')
432 T = self.L.findType('newton')
433 M = T.getSimulation('sim')
434 M.solve(ascpy.Solver("QRSlv"),ascpy.SolverReporter())
435 I = ascpy.Integrator(M)
436 I.setEngine('LSODE')
437 I.setParameter('rtolvect',False)
438 I.setParameter('rtol',1e-7)
439 I.setParameter('atolvect',False)
440 I.setParameter('atol',1e-7)
441 I.setMinSubStep(1e-7)
442 I.setMaxSubStep(0.001)
443 I.setMaxSubSteps(10000)
444
445 I.setReporter(ascpy.IntegratorReporterConsole(I))
446 I.setLinearTimesteps(ascpy.Units("s"), 0, 2*float(M.v)/float(M.g), 2)
447 I.analyse()
448 I.solve()
449 print "At end of simulation,"
450 print "x = %f" % M.x
451 print "v = %f" % M.v
452 M.run(T.getMethod('self_test'))
453
454 def testlotka(self):
455 self.L.load('johnpye/lotka.a4c')
456 M = self.L.findType('lotka').getSimulation('sim')
457 M.setSolver(ascpy.Solver("QRSlv"))
458 I = ascpy.Integrator(M)
459 I.setEngine('LSODE')
460 I.setReporter(ascpy.IntegratorReporterConsole(I))
461 I.setLinearTimesteps(ascpy.Units("s"), 0, 200, 5)
462 I.analyse()
463 print "Number of vars = %d" % I.getNumVars()
464 assert I.getNumVars()==2
465 I.solve()
466 assert I.getNumObservedVars() == 3;
467 assert abs(M.R - 832) < 1.0
468 assert abs(M.F - 21.36) < 0.1
469
470 def testwritegraph(self):
471 self.L.load('johnpye/lotka.a4c')
472 M = self.L.findType('lotka').getSimulation('sim')
473 F = file('lotka.png','w')
474 M.build()
475 M.write(F,"dot")
476
477
478 #-------------------------------------------------------------------------------
479 # Testing of a external blackbox functions
480
481 class TestBlackBox(AscendSelfTester):
482 def testparsefail0(self):
483 try:
484 self.L.load('test/blackbox/parsefail0.a4c')
485 self.fail("parsefail0 should not have loaded without errors")
486 except:
487 pass
488
489 def testparsefail1(self):
490 try:
491 self.L.load('test/blackbox/parsefail1.a4c')
492 self.fail("parsefail1 should not have loaded without errors")
493 except:
494 pass
495
496 def testparsefail2(self):
497 try:
498 self.L.load('test/blackbox/parsefail2.a4c')
499 self.fail("parsefail2 should not have loaded without errors")
500 except:
501 pass
502
503 def testparsefail3(self):
504 try:
505 self.L.load('test/blackbox/parsefail3.a4c')
506 self.fail("parsefail3 should not have loaded without errors")
507 except:
508 pass
509
510 def testparsefail4(self):
511 try:
512 self.L.load('test/blackbox/parsefail4.a4c')
513 self.fail("parsefail4 should not have loaded")
514 except:
515 pass
516
517 def testfail1(self):
518 """Mismatched arg counts check-- tests bbox, not ascend."""
519 self.L.load('test/blackbox/fail1.a4c')
520 try:
521 M = self.L.findType('fail1').getSimulation('sim')
522 self.fail("expected exception was not raised")
523 except RuntimeError,e:
524 print "Caught exception '%s', assumed ok" % e
525
526 def testfail2(self):
527 """Incorrect data arg check -- tests bbox, not ascend"""
528 self.L.load('test/blackbox/fail2.a4c')
529 try:
530 M = self.L.findType('fail2').getSimulation('sim')
531 self.fail("expected exception was not raised")
532 except RuntimeError,e:
533 print "Caught exception '%s', assumed ok (should mention errors during instantiation)" % e
534
535 def testpass1(self):
536 """simple single bbox forward solve"""
537 M = self._run('pass1',filename='test/blackbox/pass.a4c')
538
539 def testpass2(self):
540 """simple single bbox reverse solve"""
541 M = self._run('pass2',filename='test/blackbox/pass.a4c')
542
543 def testpass3(self):
544 """simple double bbox solve"""
545 M = self._run('pass3',filename='test/blackbox/pass3.a4c')
546
547 def testpass4(self):
548 """simple double bbox reverse solve"""
549 M = self._run('pass4',filename='test/blackbox/pass3.a4c')
550
551 def testpass5(self):
552 M = self._run('pass5',filename='test/blackbox/pass5.a4c')
553
554 def testpass6(self):
555 M = self._run('pass6',filename='test/blackbox/pass5.a4c')
556
557 def testpass7(self):
558 M = self._run('pass7',filename='test/blackbox/passmerge.a4c')
559
560 def testpass8(self):
561 M = self._run('pass8',filename='test/blackbox/passmerge.a4c')
562
563 def testpass9(self):
564 M = self._run('pass9',filename='test/blackbox/passmerge.a4c')
565
566 def testpass10(self):
567 M = self._run('pass10',filename='test/blackbox/passmerge.a4c')
568
569 def testpass11(self):
570 M = self._run('pass11',filename='test/blackbox/passmerge.a4c')
571
572 def testpass12(self):
573 M = self._run('pass12',filename='test/blackbox/passmerge.a4c')
574
575 # this test doesn't work: 'system is inconsistent' -- and structurally singular
576 # def testpass13(self):
577 # """cross-merged input/output solve"""
578 # M = self._run('pass13',filename='test/blackbox/passmerge.a4c')
579
580 def testpass14(self):
581 """cross-merged input/output reverse solve"""
582 M = self._run('pass14',filename='test/blackbox/passmerge.a4c')
583
584 def testpass20(self):
585 M = self._run('pass20',filename='test/blackbox/passarray.a4c')
586
587 def testparsefail21(self):
588 """dense array of black boxes wrong syntax"""
589 try:
590 self.L.load('test/blackbox/parsefail21.a4c')
591 self.fail("parsefail21 should not have loaded without errors")
592 except:
593 pass
594
595 def testpass22(self):
596 M = self._run('pass22',filename='test/blackbox/passarray.a4c')
597
598 def testpass23(self):
599 M = self._run('pass23',filename='test/blackbox/passarray.a4c')
600
601 def testpass61(self):
602 M = self._run('pass61',filename='test/blackbox/reinstantiate.a4c')
603
604 def testpass62(self):
605 M = self._run('pass62',filename='test/blackbox/reinstantiate.a4c')
606
607 def testpass64(self):
608 M = self._run('pass64',filename='test/blackbox/reinstantiate.a4c')
609
610 def testpass65(self):
611 M = self._run('pass65',filename='test/blackbox/reinstantiate.a4c')
612
613 def testpass66(self):
614 M = self._run('pass66',filename='test/blackbox/reinstantiate.a4c')
615
616 def testpass67(self):
617 M = self._run('pass67',filename='test/blackbox/reinstantiate.a4c')
618
619 class TestExtFn(AscendSelfTester):
620 def testextfntest(self):
621 M = self._run('extfntest',filename='johnpye/extfn/extfntest.a4c')
622 self.assertAlmostEqual(M.y, 2);
623 self.assertAlmostEqual(M.x, 1);
624 self.assertAlmostEqual(M.y, M.x + 1);
625
626 def testextrelfor(self):
627 M = self._run('extrelfor',filename='johnpye/extfn/extrelfor.a4c')
628
629 ## @TODO fix bug with badly-named bbox rel in a loop (Ben, maybe)
630 # def testextrelforbadnaming(self):
631 # self.L.load('johnpye/extfn/extrelforbadnaming.a4c')
632 # T = self.L.findType('extrelfor')
633 # M = T.getSimulation('sim')
634 # M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
635 # print "x[1] = %f" % M.x[1]
636 # print "x[2] = %f" % M.x[2]
637 # print "x[3] = %f" % M.x[3]
638 # print "x[4] = %f" % M.x[4]
639 # print "x[5] = %f" % M.x[5]
640 # M.run(T.getMethod('self_test'))
641
642 def testextrelrepeat(self):
643 M = self._run('extrelrepeat',filename='johnpye/extfn/extrelrepeat.a4c')
644
645 #-------------------------------------------------------------------------------
646 # Testing of Sensitivity module
647
648 class TestSensitivity(AscendSelfTester):
649 def test1(self):
650 self.L.load('sensitivity_test.a4c')
651 T = self.L.findType('sensitivity_test')
652 M = T.getSimulation('sim',False)
653 M.run(T.getMethod('on_load'))
654 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
655 M.run(T.getMethod('analyse'))
656 M.run(T.getMethod('self_test'))
657
658 # def testall(self):
659 # self.L.load('sensitivity_test.a4c')
660 # T = self.L.findType('sensitivity_test_all')
661 # M = T.getSimulation('sim',False)
662 # M.run(T.getMethod('on_load'))
663 # M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
664 # M.run(T.getMethod('analyse'))
665 # M.run(T.getMethod('self_test'))
666 # CAUSES CRASH
667
668 #-------------------------------------------------------------------------------
669 # Testing of a ExtPy - external python methods
670
671 class TestExtPy(AscendSelfTester):
672 def test1(self):
673 self.L.load('johnpye/extpy/extpytest.a4c')
674 T = self.L.findType('extpytest')
675 M = T.getSimulation('sim')
676 M.run(T.getMethod('self_test'))
677
678 def test2(self):
679 self.L.load('johnpye/extpy/extpytest.a4c')
680 T = self.L.findType('extpytest')
681 M = T.getSimulation('sim')
682 M.run(T.getMethod('pythonthing'))
683 M.run(T.getMethod('pythonthing'))
684 M.run(T.getMethod('pythonthing'))
685 M.run(T.getMethod('pythonthing'))
686 M.run(T.getMethod('pythonthing'))
687 M.run(T.getMethod('pythonthing'))
688 M.run(T.getMethod('pythonthing'))
689 # causes crash!
690
691 #-------------------------------------------------------------------------------
692 # Testing of saturated steam properties library (iapwssatprops.a4c)
693
694 class TestSteam(AscendSelfTester):
695
696 def testiapwssatprops1(self):
697 M = self._run('testiapwssatprops1',filename='steam/iapwssatprops.a4c')
698 def testiapwssatprops2(self):
699 M = self._run('testiapwssatprops2',filename='steam/iapwssatprops.a4c')
700 def testiapwssatprops3(self):
701 M = self._run('testiapwssatprops3',filename='steam/iapwssatprops.a4c')
702
703 # test the stream model basically works
704 def testsatsteamstream(self):
705 M = self._run('satsteamstream',filename='steam/satsteamstream.a4c')
706
707 # test that we can solve in terms of various (rho,u)
708 def testsatuv(self):
709 self.L.load('steam/iapwssat.a4c')
710 T = self.L.findType('testiapwssatuv')
711 M = T.getSimulation('sim',False)
712 M.run(T.getMethod('on_load'))
713 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
714 print "p = %f bar" % M.p.to('bar');
715 print "T = %f C" % (M.T.to('K') - 273.15);
716 print "x = %f" % M.x;
717 M.run(T.getMethod('self_test'))
718 M.run(T.getMethod('values2'))
719 # M.v.setRealValueWithUnits(1.0/450,"m^3/kg");
720 # M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
721 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
722 print "p = %f bar" % M.p.to('bar');
723 print "T = %f C" % (M.T.to('K') - 273.15);
724 print "x = %f" % M.x;
725 M.run(T.getMethod('self_test2'))
726
727
728 ## @TODO fix error capture from bounds checking during initialisation
729 # def testiapwssat1(self):
730 # M = self._run('testiapwssat1',filename='steam/iapwssat.a4c')
731
732 def testdsgsat(self):
733 self.L.load('steam/dsgsat3.a4c')
734 T = self.L.findType('dsgsat3')
735 M = T.getSimulation('sim',False)
736 M.run(T.getMethod('on_load'))
737 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
738 self.assertAlmostEqual(M.dTw_dt[2],0.0);
739 M.run(T.getMethod('configure_dynamic'))
740 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
741 return M
742
743 def testdsgsatrepeat(self):
744 self.L.load('steam/dsgsat3.a4c')
745 T = self.L.findType('dsgsat3')
746 M = T.getSimulation('sim',False)
747 M.run(T.getMethod('on_load'))
748 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
749 M.run(T.getMethod('on_load'))
750 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
751 M.run(T.getMethod('on_load'))
752 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
753
754 def testvary(self):
755 self.L.load('steam/dsgsat3.a4c')
756 T = self.L.findType('dsgsat3')
757 M = T.getSimulation('sim',False)
758 M.run(T.getMethod('on_load'))
759 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
760 print "----- setting qdot_s -----"
761 M.qdot_s.setRealValueWithUnits(1000,"W/m")
762 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
763 print "----- setting qdot_s -----"
764 M.qdot_s.setRealValueWithUnits(2000,"W/m")
765 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
766
767 def teststeadylsode(self):
768 "test that steady conditions are stable with LSODE"
769 M = self.testdsgsat()
770 #M.qdot_s.setRealValueWithUnits(1000,"W/m")
771 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
772 #M.setParameter('
773 I = ascpy.Integrator(M)
774 I.setEngine('LSODE')
775 I.setReporter(ascpy.IntegratorReporterConsole(I))
776 I.setLinearTimesteps(ascpy.Units("s"), 0, 3600, 10)
777 I.analyse()
778 I.solve()
779
780 # def testpeturblsode(self):
781 # "test that steady conditions are stable with LSODE"
782 # M = self.testdsgsat()
783 # # here is the peturbation...
784 # M.qdot_s.setRealValueWithUnits(1000,"W/m")
785 # M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
786 # I = ascpy.Integrator(M)
787 # I.setEngine('LSODE')
788 # I.setReporter(ascpy.IntegratorReporterConsole(I))
789 # I.setLinearTimesteps(ascpy.Units("s"), 0, 5, 1)
790 # I.analyse()
791 # I.solve()
792
793 def teststeadyida(self):
794 M = self.testdsgsat()
795 self.assertAlmostEqual(M.dTw_dt[2],0.0)
796 Tw1 = float(M.T_w[2])
797 T = self.L.findType('dsgsat3')
798 M.run(T.getMethod('free_states'))
799 I = ascpy.Integrator(M)
800 I.setEngine('IDA')
801 I.setParameter('linsolver','DENSE')
802 I.setParameter('safeeval',True)
803 I.setParameter('rtol',1e-4)
804 I.setParameter('atolvect',False)
805 I.setParameter('atol',1e-4)
806 I.setParameter('maxord',3)
807 I.setInitialSubStep(0.001)
808 I.setReporter(ascpy.IntegratorReporterConsole(I))
809 I.setLinearTimesteps(ascpy.Units("s"), 0, 3600, 10)
810 I.analyse()
811 I.solve()
812 self.assertAlmostEqual(float(M.T_w[2]),Tw1)
813 M.qdot_s.setRealValueWithUnits(1000,"W/m")
814 self.assertAlmostEqual(M.qdot_s.to("W/m"),1000)
815 M.solve(ascpy.Solver('QRSlv'),ascpy.SolverReporter())
816 print "dTw/dt = %f" % M.dTw_dt[2]
817 self.assertNotAlmostEqual(M.dTw_dt[2],0.0)
818 F=file('dsgsat.dot','w')
819 M.write(F,'dot')
820
821 def testpeturbida(self):
822 M = self.testdsgsat()
823 self.assertAlmostEqual(M.dTw_dt[2],0.0)
824 T = self.L.findType('dsgsat3')
825 M.run(T.getMethod('free_states'))
826 # here is the peturbation...
827 qdot_s = float(M.qdot_s)
828 print "OLD QDOT_S = %f" % qdot_s
829 M.qdot_s.setRealValueWithUnits(6000,"W/m")
830 # IDA has its own initial conditions solver, so no need to call QRSlv here
831 I = ascpy.Integrator(M)
832 I.setEngine('IDA')
833 I.setParameter('linsolver','DENSE')
834 I.setReporter(ascpy.IntegratorReporterConsole(I))
835 #I.setLinearTimesteps(ascpy.Units("s"), 0,300,300)
836 I.setLogTimesteps(ascpy.Units("s"), 0.009, 1200, 150)
837 I.analyse()
838 F = file('ga.mm','w')
839 I.writeMatrix(F,'dg/dz')
840 F = file('gd.mm','w')
841 I.writeMatrix(F,'dg/dx')
842 F = file('fa.mm','w')
843 I.writeMatrix(F,'df/dz')
844 F = file('fd.mm','w')
845 I.writeMatrix(F,'df/dx')
846 F = file('fdp.mm','w')
847 I.writeMatrix(F,"df/dx'")
848 I.solve()
849
850 #-------------------------------------------------------------------------------
851 # Testing of freesteam external steam properties functions
852
853 with_freesteam = True
854 try:
855 # we assume that if the freesteam python module is installed, the ASCEND
856 # external library will also be.
857 import freesteam
858 have_freesteam = True
859 except ImportError,e:
860 have_freesteam = False
861
862 if with_freesteam and have_freesteam:
863 class TestFreesteam(AscendSelfTester):
864 # def testfreesteamtest(self):
865 # """run the self-test cases bundled with freesteam"""
866 # self._run('testfreesteam',filename='testfreesteam.a4c')
867
868 def testload(self):
869 """check that we can load 'thermalequilibrium2' (IMPORT "freesteam", etc)"""
870 self.L.load('johnpye/thermalequilibrium2.a4c')
871
872 def testinstantiate(self):
873 """load an instantiate 'thermalequilibrium2'"""
874 self.testload()
875 M = self.L.findType('thermalequilibrium2').getSimulation('sim')
876 return M
877
878 def testintegrate(self):
879 """integrate transfer of heat from one mass of water/steam to another
880 according to Newton's law of cooling"""
881 M = self.testinstantiate()
882 M.setSolver(ascpy.Solver("QRSlv"))
883 I = ascpy.Integrator(M)
884 I.setEngine('LSODE')
885 I.setReporter(ascpy.IntegratorReporterConsole(I))
886 I.setLinearTimesteps(ascpy.Units("s"), 0, 3000, 30)
887 I.setMinSubStep(0.01)
888 I.setInitialSubStep(1)
889 I.analyse()
890 print "Number of vars = %d" % I.getNumVars()
891 assert I.getNumVars()==2
892 I.solve()
893 assert I.getNumObservedVars() == 3;
894 print "S[1].T = %f K" % M.S[1].T
895 print "S[2].T = %f K" % M.S[2].T
896 print "Q = %f W" % M.Q
897 self.assertAlmostEqual(float(M.S[1].T),506.77225109,4);
898 self.assertAlmostEqual(float(M.S[2].T),511.605173967,5);
899 self.assertAlmostEqual(float(M.Q),-48.32922877329,3);
900 self.assertAlmostEqual(float(M.t),3000);
901 print "Note that the above values have not been verified analytically"
902
903 def testcollapsingcan2(self):
904 """ solve the collapsing can model using IAPWS-IF97 steam props """
905 M = self._run("collapsingcan2",filename="collapsingcan2.a4c");
906
907 #-------------------------------------------------------------------------------
908 # Testing of the brent-solver EXTERNAL method
909
910 class TestBrent(AscendSelfTester):
911 def testbrent(self):
912 M = self._run('brent1',filename='test/brent.a4c')
913
914 #-------------------------------------------------------------------------------
915 # Testing of IDA's analysis module
916
917 class TestIDA(Ascend):
918 def _run(self,filen,modeln=""):
919 self.L.load('test/ida/%s.a4c' % filen)
920 T = self.L.findType('%s%s' % (filen,modeln))
921 M = T.getSimulation('sim')
922 M.build()
923 I = ascpy.Integrator(M)
924 I.setEngine('IDA')
925 I.analyse()
926 return M;
927
928 def _runfail(self,filen,n,msg="failed"):
929 try:
930 self._run(filen,'fail%d' % n)
931 except Exception,e:
932 print "(EXPECTED) ERROR: %s" % e
933 return
934 self.fail(msg)
935
936 def testsinglederiv(self):
937 self._run('singlederiv')
938
939 def testsinglederivfail1(self):
940 self._runfail('singlederiv',1
941 ,"t.ode_id=-1 did not trigger error")
942
943 def testsinglederivfail2(self):
944 self._runfail('singlederiv',2
945 ,"dy_dt.ode_id=2 did not trigger error")
946
947 def testsinglederivfail3(self):
948 self._runfail('singlederiv',3
949 ,"dy_dt.ode_type=3 did not trigger error")
950
951 def testsinglederivfail4(self):
952 self._runfail('singlederiv',4
953 ,"duplicate ode_type=1 did not trigger error")
954
955 def testsinglederivfail5(self):
956 self._runfail('singlederiv',5
957 ,"duplicate ode_type=1 did not trigger error")
958
959 def testsinglederivfail6(self):
960 self._runfail('singlederiv',6
961 ,"duplicate ode_type=1 did not trigger error")
962
963 def testtwoderiv(self):
964 self._run('twoderiv')
965
966 def testtwoderivfail1(self):
967 self._runfail('twoderiv',1)
968
969 def testtwoderivfail2(self):
970 self._runfail('twoderiv',2)
971
972 def testtwoderivfail3(self):
973 self._runfail('twoderiv',3)
974 def testtwoderivfail4(self):
975 self._runfail('twoderiv',4)
976 def testtwoderivfail5(self):
977 self._runfail('twoderiv',5)
978
979 def testnoderivs(self):
980 self._runfail('noderivs',1)
981
982 def testnoindeps(self):
983 self._runfail('indeps',1)
984
985 def testtwoindeps(self):
986 self._runfail('indeps',2)
987
988 def testfixedvars(self):
989 self._run('fixedvars')
990
991 def testfixedvars1(self):
992 self._run('fixedvars',1)
993
994 # fails the index check
995 # def testfixedvars2(self):
996 # self._run('fixedvars',2)
997
998 # fails the index check
999 # def testfixedvars3(self):
1000 # self._run('fixedvars',3)
1001
1002 def testincidence(self):
1003 self._run('incidence')
1004
1005 def testincidence1(self):
1006 self._run('incidence',1)
1007 def testincidence2(self):
1008 self._run('incidence',2)
1009 def testincidence3(self):
1010 M = self._run('incidence',3)
1011
1012 def testincidence4(self):
1013 self._run('incidence',4)
1014 def testincidencefail5(self):
1015 self._runfail('incidence',5)
1016
1017 def testwritematrix(self):
1018 self.L.load('test/ida/writematrix.a4c')
1019 T = self.L.findType('writematrix')
1020 M = T.getSimulation('sim')
1021 M.build()
1022 I = ascpy.Integrator(M)
1023 I.setEngine('IDA')
1024 I.analyse()
1025 # this stuff fails on Windows because FILE structure is different python vs mingw
1026 # F = os.tmpfile()
1027 # I.writeMatrix(F,"dF/dy")
1028 # F.seek(0)
1029 # print F.read()
1030 # F1 = os.tmpfile()
1031 # I.writeMatrix(F1,"dF/dy'")
1032 # F1.seek(0)
1033 # print F1.read()
1034 # F1 = os.tmpfile()
1035 # I.writeMatrix(F1,"dg/dx")
1036 # F1.seek(0)
1037 # print F1.read()
1038 # # for the moment you'll have to check these results manually.
1039
1040 def testwritematrix2(self):
1041 self.L.load('test/ida/writematrix.a4c')
1042 T = self.L.findType('writematrix2')
1043 M = T.getSimulation('sim')
1044 M.build()
1045 I = ascpy.Integrator(M)
1046 I.setEngine('IDA')
1047 I.analyse()
1048 # this stuff fails on Windows because FILE structure is different python vs mingw
1049 # F = os.tmpfile()
1050 # I.writeMatrix(F,"dF/dy")
1051 # F.seek(0)
1052 # print F.read()
1053 # F1 = os.tmpfile()
1054 # I.writeMatrix(F1,"dF/dy'")
1055 # F1.seek(0)
1056 # print F1.read()
1057 # F1 = os.tmpfile()
1058 # I.writeMatrix(F1,"dg/dx")
1059 # F1.seek(0)
1060 # print F1.read()
1061 #F1 = os.tmpfile()
1062 #I.writeMatrix(F1,"dydp/dyd")
1063 #F1.seek(0)
1064 #print F1.read()
1065 # for the moment you'll have to check these results manually.
1066
1067 def testindexproblem(self):
1068 self.L.load('test/ida/indexproblem.a4c')
1069 T = self.L.findType('indexproblem')
1070 M = T.getSimulation('sim')
1071 M.build()
1072 I = ascpy.Integrator(M)
1073 I.setEngine('IDA')
1074 I.analyse()
1075 pass
1076
1077 def testindexproblem2(self):
1078 self.L.load('test/ida/indexproblem.a4c')
1079 T = self.L.findType('indexproblem2')
1080 M = T.getSimulation('sim')
1081 M.build()
1082 I = ascpy.Integrator(M)
1083 I.setEngine('IDA')
1084 try:
1085 I.analyse()
1086 except Exception,e:
1087 return
1088 self.fail("Index problem not detected")
1089
1090 def testboundaries(self):
1091 self.L.load('test/ida/boundaries.a4c')
1092 T = self.L.findType('boundaries')
1093 M = T.getSimulation('sim')
1094 M.build()
1095 I = ascpy.Integrator(M)
1096 I.setEngine('IDA')
1097 I.analyse()
1098 I.setLogTimesteps(ascpy.Units("s"), 0.1, 20, 20)
1099 I.setParameter('linsolver','DENSE')
1100 I.setParameter('calcic','Y')
1101 I.setParameter('linsolver','DENSE')
1102 I.setParameter('safeeval',False)
1103 I.setReporter(ascpy.IntegratorReporterConsole(I))
1104 I.solve()
1105
1106 # doesn't work yet:
1107 # def testincidence5(self):
1108 # self._run('incidence',5)
1109
1110
1111 #-------------------------------------------------------------------------------
1112 # Testing of IDA models using DENSE linear solver
1113
1114 class TestIDADENSE(Ascend):
1115 """IDA DAE integrator, DENSE linear solver"""
1116
1117 def testlotka(self):
1118 self.L.load('johnpye/lotka.a4c')
1119 M = self.L.findType('lotka').getSimulation('sim')
1120 M.setSolver(ascpy.Solver("QRSlv"))
1121 I = ascpy.Integrator(M)
1122 I.setEngine('IDA')
1123 I.setReporter(ascpy.IntegratorReporterConsole(I))
1124 I.setLinearTimesteps(ascpy.Units("s"), 0, 200, 5);
1125 I.setParameter('linsolver','DENSE')
1126 I.setParameter('rtol',1e-8);
1127 I.analyse()
1128 assert I.getNumVars()==2
1129 assert abs(M.R - 1000) < 1e-300
1130 I.solve()
1131 assert I.getNumObservedVars() == 3
1132 assert abs(M.R - 832) < 1.0
1133 assert abs(M.F - 21.36) < 0.1
1134
1135 def testdenx(self):
1136 print "-----------------------------====="
1137 self.L.load('johnpye/idadenx.a4c')
1138 M = self.L.findType('idadenx').getSimulation('sim')
1139 M.setSolver(ascpy.Solver("QRSlv"))
1140 I = ascpy.Integrator(M)
1141 I.setEngine('IDA')
1142 I.setParameter('calcic','YA_YDP')
1143 I.setParameter('linsolver','DENSE')
1144 I.setParameter('safeeval',True)
1145 I.setReporter(ascpy.IntegratorReporterConsole(I))
1146 I.setLogTimesteps(ascpy.Units("s"), 0.4, 4e10, 11)
1147 I.setMaxSubStep(0);
1148 I.setInitialSubStep(0)
1149 I.setMaxSubSteps(0)
1150 I.setParameter('autodiff',True)
1151 I.analyse()
1152 I.solve()
1153 assert abs(float(M.y1) - 5.1091e-08) < 2e-9
1154 assert abs(float(M.y2) - 2.0437e-13) < 2e-14
1155 assert abs(float(M.y3) - 1.0) < 1e-5
1156
1157 def testhires(self):
1158 self.L.load('test/hires.a4c')
1159 T = self.L.findType('hires')
1160 M = T.getSimulation('sim')
1161 M.setSolver(ascpy.Solver('QRSlv'))
1162 I = ascpy.Integrator(M)
1163 I.setEngine('IDA')
1164 I.setParameter('linsolver','DENSE')
1165 I.setParameter('rtol',1.1e-15)
1166 I.setParameter('atolvect',0)
1167 I.setParameter('atol',1.1e-15)
1168 I.setReporter(ascpy.IntegratorReporterConsole(I))
1169 I.setLogTimesteps(ascpy.Units(""), 1, 321.8122, 5)
1170 I.setInitialSubStep(1e-5)
1171 I.setMaxSubSteps(10000)
1172 I.analyse()
1173 I.solve()
1174 for i in range(8):
1175 print "y[%d] = %.20g" % (i+1, M.y[i+1])
1176 M.run(T.getMethod('self_test'))
1177
1178 def testchemakzo(self):
1179 self.L.load('test/chemakzo.a4c')
1180 T = self.L.findType('chemakzo')
1181 M = T.getSimulation('sim')
1182 M.setSolver(ascpy.Solver('QRSlv'))
1183 I = ascpy.Integrator(M)
1184 I.setEngine('IDA')
1185 I.setParameter('linsolver','DENSE')
1186 I.setParameter('rtol',1e-15)
1187 I.setParameter('atolvect',0)
1188 I.setParameter('atol',1e-15)
1189 I.setReporter(ascpy.IntegratorReporterConsole(I))
1190 I.setLinearTimesteps(ascpy.Units("s"), 1, 180, 5)
1191 I.setInitialSubStep(1e-13)
1192 I.setMaxSubSteps(10000)
1193 I.analyse()
1194 I.solve()
1195 for i in range(6):
1196 print "y[%d] = %.20g" % (i+1, M.y[i+1])
1197 M.run(T.getMethod('self_test'))
1198
1199 def testtransamp(self):
1200 self.L.load('test/transamp.a4c')
1201 T = self.L.findType('transamp')
1202 M = T.getSimulation('sim')
1203 M.setSolver(ascpy.Solver('QRSlv'))
1204 I = ascpy.Integrator(M)
1205 I.setEngine('IDA')
1206 I.setParameter('linsolver','DENSE')
1207 I.setParameter('rtol',1e-7)
1208 I.setParameter('atolvect',0)
1209 I.setParameter('atol',1e-7)
1210 I.setReporter(ascpy.IntegratorReporterConsole(I))
1211 I.setLinearTimesteps(ascpy.Units("s"), 0.05, 0.2, 20)
1212 I.setInitialSubStep(0.00001)
1213 I.setMaxSubSteps(10000)
1214 I.analyse()
1215 I.solve()
1216 for i in range(6):
1217 print "y[%d] = %.20g" % (i+1, M.y[i+1])
1218 M.run(T.getMethod('self_test'))
1219
1220 # MODEL FAILS ANALYSIS: we need to add support for non-incident differential vars
1221 # def testpollution(self):
1222 # self.L.load('test/pollution.a4c')
1223 # T = self.L.findType('pollution')
1224 # M = T.getSimulation('sim')
1225 # M.setSolver(ascpy.Solver('QRSlv'))
1226 # I = ascpy.Integrator(M)
1227 # I.setEngine('IDA')
1228 # I.setParameter('linsolver','DENSE')
1229 # I.setParameter('rtol',1.1e-15)
1230 # I.setParameter('atolvect',0)
1231 # I.setParameter('atol',1.1e-15)
1232 # I.setReporter(ascpy.IntegratorReporterConsole(I))
1233 # I.setLogTimesteps(ascpy.Units("s"), 1, 60, 5)
1234 # I.setInitialSubStep(1e-5)
1235 # I.setMaxSubSteps(10000)
1236 # I.analyse()
1237 # I.solve()
1238 # for i in range(20):
1239 # print "y[%d] = %.20g" % (i+1, M.y[i+1])
1240 # M.run(T.getMethod('self_test'))
1241
1242 ## @TODO fails during IDACalcIC (model too big?)
1243 # def testkryx(self):
1244 # self.L.load('johnpye/idakryx.a4c')
1245 # ascpy.getCompiler().setUseRelationSharing(False)
1246 # M = self.L.findType('idakryx').getSimulation('sim')
1247 # M.setSolver(ascpy.Solver('QRSlv'))
1248 # M.build()
1249 # I = ascpy.Integrator(M)
1250 # I.setEngine('IDA')
1251 # I.setReporter(ascpy.IntegratorReporterConsole(I))
1252 # I.setParameter('linsolver','DENSE')
1253 # I.setParameter('maxl',8)
1254 # I.setParameter('gsmodified',False)
1255 # I.setParameter('autodiff',True)
1256 # I.setParameter('rtol',0)
1257 # I.setParameter('atol',1e-3);
1258 # I.setParameter('atolvect',False)
1259 # I.setParameter('calcic','YA_YDP')
1260 # I.analyse()
1261 # I.setLogTimesteps(ascpy.Units("s"), 0.01, 10.24, 11)
1262 # I.solve()
1263 # assert abs(M.u[2][2].getValue()) < 1e-5
1264
1265 #-------------------------------------------------------------------------------
1266 # Testing of IDA models using SPGMR linear solver (Krylov)
1267
1268 # these tests are disabled until SPGMR preconditioning has been implemented
1269 class TestIDASPGMR:#(Ascend):
1270 def testlotka(self):
1271 self.L.load('johnpye/lotka.a4c')
1272 M = self.L.findType('lotka').getSimulation('sim')
1273 M.setSolver(ascpy.Solver("QRSlv"))
1274 I = ascpy.Integrator(M)
1275 I.setEngine('IDA')
1276 I.setReporter(ascpy.IntegratorReporterConsole(I))
1277 I.setLinearTimesteps(ascpy.Units("s"), 0, 200, 5)
1278 I.setParameter('rtol',1e-8)
1279 I.analyse()
1280 assert I.getNumVars()==2
1281 assert abs(M.R - 1000) < 1e-300
1282 I.solve()
1283 assert I.getNumObservedVars() == 3
1284 assert abs(M.R - 832) < 1.0
1285 assert abs(M.F - 21.36) < 0.1
1286
1287
1288 def testkryx(self):
1289 self.L.load('johnpye/idakryx.a4c')
1290 M = self.L.findType('idakryx').getSimulation('sim')
1291 M.build()
1292 I = ascpy.Integrator(M)
1293 I.setEngine('IDA')
1294 I.setReporter(ascpy.IntegratorReporterConsole(I))
1295 I.setParameter('linsolver','SPGMR')
1296 I.setParameter('prec','JACOBI')
1297 I.setParameter('maxl',8)
1298 I.setParameter('gsmodified',False)
1299 I.setParameter('autodiff',True)
1300 I.setParameter('gsmodified',True)
1301 I.setParameter('rtol',0)
1302 I.setParameter('atol',1e-3);
1303 I.setParameter('atolvect',False)
1304 I.setParameter('calcic','Y')
1305 I.analyse()
1306 I.setLogTimesteps(ascpy.Units("s"), 0.01, 10.24, 10);
1307 print M.udot[1][3]
1308 I.solve()
1309 assert 0
1310
1311 def testzill(self):
1312 self.L.load('johnpye/zill.a4c')
1313 T = self.L.findType('zill')
1314 M = T.getSimulation('sim')
1315 M.setSolver(ascpy.Solver('QRSlv'))
1316 I = ascpy.Integrator(M)
1317 I.setEngine('IDA')
1318 I.setParameter('safeeval',False)
1319 I.setMinSubStep(1e-7)
1320 I.setMaxSubStep(0.001)
1321 I.setMaxSubSteps(10000)
1322 I.setReporter(ascpy.IntegratorReporterConsole(I))
1323 I.setLinearTimesteps(ascpy.Units(), 1.0, 1.5, 5)
1324 I.analyse()
1325 I.solve()
1326 M.run(T.getMethod('self_test'))
1327
1328 def testdenxSPGMR(self):
1329 self.L.load('johnpye/idadenx.a4c')
1330 M = self.L.findType('idadenx').getSimulation('sim')
1331 M.setSolver(ascpy.Solver('QRSlv'))
1332 I = ascpy.Integrator(M)
1333 I.setEngine('IDA')
1334 I.setReporter(ascpy.IntegratorReporterConsole(I))
1335 I.setLogTimesteps(ascpy.Units("s"), 0.4, 4e10, 11)
1336 I.setMaxSubStep(0);
1337 I.setInitialSubStep(0);
1338 I.setMaxSubSteps(0);
1339 I.setParameter('autodiff',True)
1340 I.setParameter('linsolver','SPGMR')
1341 I.setParameter('gsmodified',False)
1342 I.setParameter('maxncf',10)
1343 I.analyse()
1344 I.solve()
1345 assert abs(float(M.y1) - 5.1091e-08) < 1e-10
1346 assert abs(float(M.y2) - 2.0437e-13) < 1e-15
1347 assert abs(float(M.y3) - 1.0) < 1e-5
1348
1349 class TestDOPRI5(Ascend):
1350 def testlotka(self):
1351 self.L.load('test/dopri5/dopri5test.a4c')
1352 M = self.L.findType('dopri5test').getSimulation('sim')
1353 M.setSolver(ascpy.Solver("QRSlv"))
1354 M.solve(ascpy.Solver("QRSlv"),ascpy.SolverReporter())
1355 I = ascpy.Integrator(M)
1356 I.setEngine('DOPRI5')
1357 I.setReporter(ascpy.IntegratorReporterConsole(I))
1358 I.setLinearTimesteps(ascpy.Units("s"), 0, 200, 20)
1359 I.setParameter('rtol',1e-8)
1360 I.analyse()
1361 assert I.getNumVars()==1
1362 I.solve()
1363 def testaren(self):
1364 self.L.load('test/dopri5/aren.a4c')
1365 M = self.L.findType('aren').getSimulation('sim')
1366 M.setSolver(ascpy.Solver("QRSlv"))
1367 M.solve(ascpy.Solver("QRSlv"),ascpy.SolverReporter())
1368 I = ascpy.Integrator(M)
1369 I.setEngine('DOPRI5')
1370 I.setReporter(ascpy.IntegratorReporterConsole(I))
1371 #xend = 17.0652165601579625588917206249
1372 I.setLinearTimesteps(ascpy.Units("s"), 0, 17.0652165601579625588917206249, 10)
1373 I.setParameter('rtol',1e-7)
1374 I.setParameter('atol',1e-7)
1375 I.setParameter('tolvect',False)
1376 I.setMinSubStep(0);
1377 I.setMaxSubStep(0);
1378 I.setInitialSubStep(0);
1379 I.analyse()
1380 I.solve()
1381 print "y[0] = %f" % float(M.y[0])
1382 assert abs(float(M.y[0]) - 0.994) < 1e-5
1383 assert abs(float(M.y[1]) - 0.0) < 1e-5
1384
1385 class TestIPOPT(Ascend):
1386
1387 def ipopt_tester(self,testname,hessian_approx='limited-memory',linear_solver='mumps'):
1388 self.L.load('test/ipopt/%s.a4c' % testname)
1389 T = self.L.findType(testname)
1390 M = T.getSimulation('sim')
1391 M.setSolver(ascpy.Solver("IPOPT"))
1392 M.setParameter('linear_solver',linear_solver)
1393 M.setParameter('hessian_approximation',hessian_approx)
1394 M.solve(ascpy.Solver("IPOPT"),ascpy.SolverReporter())
1395 M.run(T.getMethod('self_test'))
1396
1397 def test2(self):
1398 self.ipopt_tester('test2')
1399
1400 def test3(self):
1401 self.ipopt_tester('test3')
1402
1403 def test4(self):
1404 self.ipopt_tester('test4')
1405
1406 def test5(self):
1407 self.ipopt_tester('test5')
1408
1409 def test6(self):
1410 self.ipopt_tester('test6')
1411
1412 def test7(self):
1413 self.ipopt_tester('test7')
1414
1415 def test7_hsl(self):
1416 self.ipopt_tester('test7',linear_solver='ma27')
1417
1418 def test8(self):
1419 self.ipopt_tester('test8')
1420
1421 def test9(self):
1422 self.ipopt_tester('test9')
1423
1424 def test10(self):
1425 self.ipopt_tester('test10')
1426
1427 def test11(self):
1428 self.ipopt_tester('test11')
1429
1430 def test12(self):
1431 self.ipopt_tester('test12')
1432
1433 def test13(self):
1434 self.ipopt_tester('test13')
1435
1436 def test14(self):
1437 self.ipopt_tester('test14')
1438
1439 # tests with exact hessian routines...
1440
1441 def test2_exact(self):
1442 self.ipopt_tester('test2',hessian_approx='exact')
1443
1444 def test3_exact(self):
1445 self.ipopt_tester('test3',hessian_approx='exact')
1446
1447 def test4_exact(self):
1448 self.ipopt_tester('test4',hessian_approx='exact')
1449
1450 def test5_exact(self):
1451 self.ipopt_tester('test5',hessian_approx='exact')
1452
1453 def test6_exact(self):
1454 self.ipopt_tester('test6',hessian_approx='exact')
1455
1456 def test7_exact(self):
1457 self.ipopt_tester('test7',hessian_approx='exact')
1458
1459 def test8_exact(self):
1460 self.ipopt_tester('test8',hessian_approx='exact')
1461
1462 def test9_exact(self):
1463 self.ipopt_tester('test9',hessian_approx='exact')
1464
1465 def test10_exact(self):
1466 self.ipopt_tester('test10',hessian_approx='exact')
1467
1468 def test11_exact(self):
1469 self.ipopt_tester('test11',hessian_approx='exact')
1470
1471 def test12_exact(self):
1472 self.ipopt_tester('test12',hessian_approx='exact')
1473
1474 def test13_exact(self):
1475 self.ipopt_tester('test13',hessian_approx='exact')
1476
1477 def test14_exact(self):
1478 self.ipopt_tester('test14',hessian_approx='exact')
1479
1480 def testformula(self):
1481 self.ipopt_tester('formula')
1482
1483 class TestCSV(Ascend):
1484 def test1(self):
1485 self.L.load('johnpye/datareader/testcsv.a4c')
1486 M = self.L.findType('testcsv').getSimulation('sim')
1487 M.solve(ascpy.Solver("QRSlv"),ascpy.SolverReporter())
1488
1489
1490 class TestSlvReq(Ascend):
1491 def test1(self):
1492 self.L.load('test/slvreq/test1.a4c')
1493 H = ascpy.SolverHooks(ascpy.SolverReporter())
1494 ascpy.SolverHooksManager_Instance().setHooks(H)
1495 T = self.L.findType('test1')
1496 M = T.getSimulation('sim')
1497 print "\n\n\nRUNNING ON_LOAD EXPLICITLY NOW..."
1498 M.run(T.getMethod('on_load'))
1499
1500 def test2(self):
1501 self.L.load('test/slvreq/test1.a4c')
1502 R = ascpy.SolverReporter()
1503 class SolverHooksPython(ascpy.SolverHooks):
1504 def __init__(self):
1505 print "PYTHON SOLVER HOOKS"
1506 ascpy.SolverHooks.__init__(self,None)
1507 def setSolver(self,solvername,sim):
1508 sim.setSolver(ascpy.Solver(solvername))
1509 print "PYTHON: SOLVER is now %s" % sim.getSolver().getName()
1510 return 0
1511 def setOption(self,optionname,val,sim):
1512 try:
1513 PP = sim.getParameters()
1514 except Exception,e:
1515 print "PYTHON ERROR: ",str(e)
1516 return ascpy.SLVREQ_OPTIONS_UNAVAILABLE
1517 try:
1518 for P in PP:
1519 if P.getName()==optionname:
1520 try:
1521 P.setValueValue(val)
1522 sim.setParameters(PP)
1523 print "PYTHON: SET",optionname,"to",repr(val)
1524 return 0
1525 except Exception,e:
1526 print "PYTHON ERROR: ",str(e)
1527 return ascpy.SLVREQ_WRONG_OPTION_VALUE_TYPE
1528 return ascpy.SLVREQ_INVALID_OPTION_NAME
1529 except Exception,e:
1530 print "PYTHON ERROR: ",str(e)
1531 return ascpy.SLVREQ_INVALID_OPTION_NAME
1532 def doSolve(self,inst,sim):
1533 try:
1534 print "PYTHON: SOLVING",sim.getName(),"WITH",sim.getSolver().getName()
1535 sim.solve(sim.getSolver(),R)
1536 except Exception,e:
1537 print "PYTHON ERROR:",str(e)
1538 return 3
1539 return 0
1540 H = SolverHooksPython()
1541 ascpy.SolverHooksManager_Instance().setHooks(H)
1542 T = self.L.findType('test1')
1543 M = T.getSimulation('sim')
1544
1545 # test some stuff for beam calculations
1546 class TestSection(Ascend):
1547 def test_compound3(self):
1548 self.L.load('johnpye/section.a4c')
1549 T = self.L.findType('compound_section_test3')
1550 M = T.getSimulation('sim')
1551 M.solve(ascpy.Solver("QRSlv"),ascpy.SolverReporter())
1552 M.run(T.getMethod('self_test'))
1553 def test_compound4(self):
1554 self.L.load('johnpye/section.a4c')
1555 T = self.L.findType('compound_section_test4')
1556 M = T.getSimulation('sim')
1557 M.solve(ascpy.Solver("QRSlv"),ascpy.SolverReporter())
1558 M.run(T.getMethod('self_test'))
1559 def test_compound2(self):
1560 self.L.load('johnpye/section.a4c')
1561 T = self.L.findType('compound_section_test2')
1562 M = T.getSimulation('sim')
1563 M.solve(ascpy.Solver("QRSlv"),ascpy.SolverReporter())
1564 M.run(T.getMethod('self_test'))
1565
1566 # move code above down here if you want to temporarily avoid testing it
1567 class NotToBeTested:
1568 def nothing(self):
1569 pass
1570
1571 def testnewton(self):
1572 sys.stderr.write("STARTING TESTNEWTON\n")
1573 self.L.load('johnpye/newton.a4c')
1574 T = self.L.findType('newton')
1575 M = T.getSimulation('sim')
1576 M.solve(ascpy.Solver("QRSlv"),ascpy.SolverReporter())
1577 I = ascpy.Integrator(M)
1578 I.setEngine('IDA')
1579 I.setParameter('linsolver','DENSE')
1580 I.setParameter('safeeval',True)
1581 I.setParameter('rtol',1e-8)
1582 I.setMaxSubStep(0.001)
1583 I.setMaxSubSteps(10000)
1584
1585 I.setReporter(ascpy.IntegratorReporterConsole(I))
1586 I.setLinearTimesteps(ascpy.Units("s"), 0, 2*float(M.v)/float(M.g), 2)
1587 I.analyse()
1588 I.solve()
1589 print "At end of simulation,"
1590 print "x = %f" % M.x
1591 print "v = %f" % M.v
1592 M.run(T.getMethod('self_test'))
1593
1594 def patchpath(VAR,SEP,addvals):
1595 restart = 0
1596 envpath = [os.path.abspath(i) for i in os.environ[VAR].split(SEP)]
1597 for l in addvals:
1598 if l in envpath[len(addvals):]:
1599 envpath.remove(l)
1600 restart = 1
1601 for l in addvals:
1602 if l not in envpath:
1603 envpath.insert(0,l)
1604 restart = 1
1605 os.environ[VAR] = SEP.join(envpath)
1606 return restart
1607
1608 if __name__=='__main__':
1609 # a whole bag of tricks to make sure we get the necessary dirs in our ascend, python and ld path vars
1610 restart = 0
1611
1612 if platform.system()=="Windows":
1613 LD_LIBRARY_PATH="PATH"
1614 SEP = ";"
1615 else:
1616 LD_LIBRARY_PATH="LD_LIBRARY_PATH"
1617 SEP = ":"
1618
1619 solverdir = os.path.abspath(os.path.join(sys.path[0],"solvers"))
1620 solverdirs = [os.path.join(solverdir,s) for s in "qrslv","cmslv","lrslv","conopt","ida","lsode","ipopt","dopri5"]
1621
1622 if not os.environ.get('ASCENDSOLVERS'):
1623 os.environ['ASCENDSOLVERS'] = SEP.join(solverdirs)
1624 restart = 1
1625 else:
1626 if patchpath('ASCENDSOLVERS',SEP,solverdirs):
1627 restart = 1
1628
1629 freesteamdir = os.path.expanduser("~/freesteam/ascend")
1630 modeldirs = [os.path.abspath(os.path.join(sys.path[0],"models")),os.path.abspath(freesteamdir)]
1631
1632 if not os.environ.get('ASCENDLIBRARY'):
1633 os.environ['ASCENDLIBRARY'] = SEP.join(modeldirs)
1634 restart = 1
1635 else:
1636 if patchpath('ASCENDLIBRARY',SEP,modeldirs):
1637 restart = 1
1638
1639 libdirs = ["ascxx","."]
1640 libdirs = [os.path.normpath(os.path.join(sys.path[0],l)) for l in libdirs]
1641 if not os.environ.get(LD_LIBRARY_PATH):
1642 os.environ[LD_LIBRARY_PATH]=SEP.join(libdirs)
1643 restart = 1
1644 else:
1645 envlibdirs = [os.path.normpath(i) for i in os.environ[LD_LIBRARY_PATH].split(SEP)]
1646 for l in libdirs:
1647 if l in envlibdirs[len(libdirs):]:
1648 envlibdirs.remove(l)
1649 restart = 1
1650 for l in libdirs:
1651 if l not in envlibdirs:
1652 envlibdirs.insert(0,l)
1653 restart = 1
1654 os.environ[LD_LIBRARY_PATH] = SEP.join(envlibdirs)
1655
1656 pypath = [os.path.normpath(os.path.join(sys.path[0],i)) for i in ['ascxx','pygtk']]
1657 if not os.environ.get('PYTHONPATH'):
1658 os.environ['PYTHONPATH']= SEP.join(pypath)
1659 else:
1660 envpypath = os.environ['PYTHONPATH'].split(SEP)
1661 for pypath1 in pypath:
1662 if pypath1 not in envpypath:
1663 envpypath.insert(0,pypath1)
1664 os.environ['PYTHONPATH']=SEP.join(envpypath)
1665 restart = 1
1666
1667 if restart and platform.system()!="Windows":
1668 script = os.path.join(sys.path[0],"test.py")
1669 sys.stderr.write("Restarting with...\n")
1670 sys.stderr.write(" export LD_LIBRARY_PATH=%s\n" % os.environ.get(LD_LIBRARY_PATH))
1671 sys.stderr.write(" export PYTHONPATH=%s\n" % os.environ.get('PYTHONPATH'))
1672 sys.stderr.write(" export ASCENDLIBRARY=%s\n" % os.environ.get('ASCENDLIBRARY'))
1673 sys.stderr.write(" export ASCENDSOLVERS=%s\n" % os.environ.get('ASCENDSOLVERS'))
1674 sys.stderr.flush()
1675 os.execvp("python",[script] + sys.argv)
1676 exit(1)
1677 else:
1678 sys.stderr.write("Got...\n")
1679 sys.stderr.write(" LD_LIBRARY_PATH=%s\n" % os.environ.get(LD_LIBRARY_PATH))
1680 sys.stderr.write(" PYTHONPATH=%s\n" % os.environ.get('PYTHONPATH'))
1681 sys.stderr.write(" ASCENDLIBRARY=%s\n" % os.environ.get('ASCENDLIBRARY'))
1682 sys.stderr.write(" ASCENDSOLVERS=%s\n" % os.environ.get('ASCENDSOLVERS'))
1683 sys.stderr.flush()
1684
1685 import ascpy
1686
1687 try:
1688 import cunit
1689 except:
1690 pass
1691
1692 atexit.register(ascpy.shutdown)
1693 #suite = unittest.TestSuite()
1694 #suite = unittest.defaultTestLoader.loadTestsFromName('__main__')
1695 #unittest.TextTestRunner(verbosity=2).run(suite)
1696 unittest.main()
1697
1698 # ex: set ts=4 :
1699

Properties

Name Value
svn:executable *

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