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

Annotation of /trunk/test.py

Parent Directory Parent Directory | Revision Log Revision Log


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

Properties

Name Value
svn:executable *

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