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

Annotation of /trunk/test.py

Parent Directory Parent Directory | Revision Log Revision Log


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

Properties

Name Value
svn:executable *

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