/[escript]/trunk/escript/py_src/modelframe.py
ViewVC logotype

Contents of /trunk/escript/py_src/modelframe.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 122 - (show annotations)
Thu Jun 9 05:38:05 2005 UTC (13 years, 10 months ago) by jgs
Original Path: trunk/esys2/escript/py_src/modelframe.py
File MIME type: text/x-python
File size: 29908 byte(s)
Merge of development branch back to main trunk on 2005-06-09

1 # $Id$
2
3 from types import StringType,IntType,FloatType,BooleanType,ListType,DictType
4 from sys import stdout
5
6 from xml.dom import minidom
7
8 _HEADER_START="<?xml version=\"1.0\"?><ESyS>\n"
9 _HEADER_END="<\ESyS>\n"
10
11 def dataNode(document, tagName, data):
12 t = document.createTextNode(str(data))
13 n = document.createElement(tagName)
14 n.appendChild(t)
15 return n
16
17 def esysDoc():
18 doc = minidom.Document()
19 esys = doc.createElement('ESys')
20 doc.appendChild(esys)
21 return doc, esys
22
23 class Link:
24 """
25 a Link makes an attribute of an object callable:
26 o.object()
27 o.a=8
28 l=Link(o,"a")
29 assert l()==8
30 """
31
32 def __init__(self,target=None,attribute=None):
33 """creates a link to the object target. If attribute is given, the link is
34 establised to this attribute of the target. Otherwise the attribute is
35 undefined."""
36 self.__target=target
37 self.setAttributeName(attribute)
38
39 def setAttributeName(self,name):
40 """set a new attribute name to be collected from the target object. The
41 target object must have the attribute with name name."""
42 if not name==None:
43 if not hasattr(self.__target,name):
44 raise AttributeError("%s: target %s has no attribute %s."%(self,self.__target,name))
45 self.__attribute=name
46
47 def hasDefinedAttributeName(self):
48 """returns true if an attribute name is set"""
49 if self.__attribute==None:
50 return False
51 else:
52 return True
53
54 def __str__(self):
55 """returns a string representation of the link"""
56 if self.hasDefinedAttributeName():
57 return "reference to attribute %s of %s"%(self.__attribute,self.__target)
58 else:
59 return "reference to target %s"%self.__target
60
61 def __call__(self,name=None):
62 """returns the value of the attribute of the target object. If the
63 atrribute is callable then the return value of the call is returned."""
64 if not name == None:
65 out=getattr(self.__target, name)
66 else:
67 out=getattr(self.__target, self.__attribute)
68 if callable(out):
69 return out()
70 else:
71 return out
72
73 def toDom(self, document, node):
74 """ toDom method of Link. Creates a Link node and appends it to the current XML
75 document """
76 link = document.createElement('link')
77 link.appendChild(dataNode(document, 'target', id(self.__target)))
78 link.appendChild(dataNode(document, 'attribute', self.__attribute))
79 node.appendChild(link)
80
81 def writeXML(self,ostream=stdout):
82 """writes an XML representation of self to the output stream ostream.
83 If ostream is nor present the standart output stream is used. If
84 esysheader==True the esys XML header is written"""
85
86 document, rootnode = esysDoc()
87 self.toDom(document, rootnode)
88
89 ostream.write(document.toprettyxml())
90
91 class LinkableObject(object):
92 """an object that allows to link its attributes to attributes of other
93 objects via a Link object. For instance
94
95 p=LinkableObject()
96 p.x=Link(o,"name")
97 print p.x
98
99 links attribute x of p to the attribute name of object o. p.x will print
100 the current value of attribute name of object o. if the value is
101 callable p.x will rturn the return value of the call.
102 """
103
104 def __init__(self,debug=False):
105 """initiates a model from a list of subsimulations. """
106 self.setDebug(debug)
107 self.__linked_attributes={}
108 #
109 # some basic fuctions:
110 #
111 def debugOn(self):
112 """sets debugging to on"""
113 self.__debug=True
114 def debugOff(self):
115 """sets debugging to off"""
116 self.__debug=False
117 def debug(self):
118 """returns True if debug mode is set to on"""
119 return self.__debug
120 def setDebug(self,flag=False):
121 """sets debugging to flag"""
122 if flag:
123 self.debugOn()
124 else:
125 self.debugOff()
126
127 def trace(self, msg):
128 if self.__debug:
129 print msg
130
131 def __getattr__(self,name):
132 """returns the value of attribute name. If the value is a Link object the
133 object is called and the return value is returned."""
134 out=self.getAttributeObject(name)
135 if isinstance(out,Link):
136 return out()
137 else:
138 return out
139
140 def getAttributeObject(self,name):
141 """return the object stored for attribute name."""
142 if self.__dict__.has_key(name):
143 out=self.__dict__[name]
144 elif self.__linked_attributes.has_key(name):
145 out=self.__linked_attributes[name]
146 else:
147 raise AttributeError,"No attribute %s."%name
148 return out
149
150 def __setattr__(self,name,value):
151 """sets the value for attribute name. If value is a Link the target
152 attribute is set to name if no attribute has been specified."""
153 if self.__dict__.has_key(name):
154 del self.__dict__[name]
155 if isinstance(value,Link):
156 if not value.hasDefinedAttributeName():
157 value.setAttributeName(name)
158 self.__linked_attributes[name]=value
159 if self.debug():
160 print "DEBUG: %s: attribute %s is now linked by %s."%(self,name,value)
161 else:
162 self.__dict__[name]=value
163
164 def __delattr__(self,name):
165 """removes the attribute name."""
166 if self.__linked_attributes.has_key[name]:
167 del self.__linked_attributes[name]
168 elif self.__dict__.has_key(name):
169 del self.__dict__[name]
170 else:
171 raise AttributeError,"No attribute %s."%name
172
173 class SimulationFrame(LinkableObject):
174 """A SimulationFrame objects represents a processess marching over time
175 until a finalizing condition is fullfilled. At each time step an iterative
176 process can be performed and the time step size can be controlled
177 """
178 UNDEF_DT=1.e300
179 MAX_TIME_STEP_REDUCTION=20
180 MAX_ITER_STEPS=50
181
182 def __init__(self,debug=False):
183 """initiates a simulation"""
184 LinkableObject.__init__(self,debug)
185
186 def doInitialization(self,t):
187 """initializes the time stepping scheme. This function may be
188 overwritten."""
189 pass
190
191 def getSafeTimeStepSize(self,dt):
192 """returns a time step size which can savely be used. This function may
193 be overwritten."""
194 return self.UNDEF_DT
195
196 def finalize(self):
197 """returns True if the time stepping is finalized. This function may be
198 overwritten."""
199 return True
200
201 def doFinalization(self):
202 """finalizes the time stepping. This function may be overwritten."""
203 pass
204
205 def doIterationInitialization(self,dt):
206 """initializes the iteration at time step t. This function may be
207 overwritten. (only called if doStep is not overwritten)"""
208 pass
209
210 def doIterationStep(self,dt):
211 """executes the iteration step. This function may be overwritten. (only
212 called if doStep is not overwritten)"""
213 pass
214
215 def terminate(self):
216 """returns True if iteration on a time step is terminated. (only called
217 if doStep is not overwritten)"""
218 return True
219
220 def doIterationFinalization(self,dt):
221 """finalalizes the iteration process. (only called if doStep is not
222 overwritten)"""
223 pass
224
225 def run(self,check_point=None):
226 """run the simulation by performing essentially
227
228 self.doInitialization()
229 while not self.finalize():
230 dt=self.getSafeTimeStepSize()
231 self.doStep(dt)
232 if n%check_point==0: self.writeXML()
233 self.doFinalization()
234
235 """
236 self.__tn=0.
237 self.__n=0
238 self.__dt=None
239 self.doInitialization(self.__tn)
240 while not self.finalize():
241 self.__n+=1
242 self.__dt=self.getSafeTimeStepSize(self.__dt)
243 if self.__dt==None: self.__dt=self.UNDEF_DT
244 if not self.__dt>0:
245 raise NonPositiveStepSizeError("non-positive step size in step %d",self.__n)
246 if self.debug(): print "%s: %d. time step %e (step size %e.)"%(self,self.__n,self.__tn+self.__dt,self.__dt)
247 endoftimestep=False
248 failcounter=0
249 while not endoftimestep:
250 endoftimestep=True
251 try:
252 self.doStep(self.__dt)
253 except FailedTimeStepError:
254 self.__dt=self.getSafeTimeStepSize(self.__dt)
255 if self.__dt==None: self.__dt=self.UNDEF_DT
256 endoftimestep=False
257 if self.debug(): print "%s: time step is repeated with new step size %e."%(self,self.__dt)
258 except IterationDivergenceError:
259 self.__dt*=0.5
260 endoftimestep=False
261 failcounter+=1
262 if failcounter>self.MAX_TIME_STEP_REDUCTION:
263 raise IterationBreakDownError("reduction of time step to achieve convergence failed.")
264
265 self.trace("%s: iteration failes. time step is repeated with new step size %e."
266 % (self,self.__dt))
267 self.__tn+=self.__dt
268 if not check_point==None:
269 if self.__n%check_point==0:
270 self.trace("%s: check point is created."%self)
271 self.writeXML()
272 self.doFinalization()
273
274 def writeXML(self):
275 raise RuntimeError, "Not implemented"
276
277 def doStep(self,dt):
278 """executes a time step by iteration. This function may be overwritten.
279
280 basicly it performs :
281
282 self.doIterationInitialization(dt)
283 while not self.terminate(): self.doIterationStep(dt)
284 self.doIterationFinalization(dt)
285
286 """
287 self.doIterationInitialization(dt)
288 self.__iter=0
289 while not self.terminate():
290 if self.debug(): print "%s: iteration step %d"%(self,self.__iter)
291 self.doIterationStep(dt)
292 self.__iter+=1
293 if self.__iter>self.MAX_ITER_STEPS:
294 raise IterationDivergenceError("%s: divergence in step %d"%(self,self.__iter))
295 self.doIterationFinalization(dt)
296
297 class Simulation(SimulationFrame):
298 """A Simulation object is comprised by SimulationFrame(s) called subsimulations."""
299
300 def __init__(self,subsimulations=[],debug=False):
301 """initiates a simulation from a list of subsimulations. """
302 SimulationFrame.__init__(self,debug)
303 self.__subsimulations=[]
304 for i in range(len(subsimulations)): self[i]=subsimulations[i]
305
306 def iterSubsimulations(self):
307 """returns an iterator over the subsimulations"""
308 return self.__subsimulations
309
310 def __getitem__(self,i):
311 """returns the i-th subsimulation"""
312 return self.__subsimulations[i]
313
314 def __setitem__(self,i,value):
315 """sets the i-th subsimulation"""
316 if not isinstance(value,SimulationFrame):
317 raise ValueError("assigned value is not a Simulation")
318 for j in range(max(i-len(self.__subsimulations)+1,0)): self.__subsimulations.append(None)
319 self.__subsimulations[i]=value
320
321 def __len__(self):
322 """returns the number of subsimulations"""
323 return len(self.__subsimulations)
324
325 def toDom(self, document, node):
326 """ toDom method of Simualtion class """
327 simulation = document.createElement('Simulation')
328 simulation.setAttribute('type', self.__class__.__name__)
329
330 for rank, sim in enumerate(self.iterSubsimulations()):
331 component = document.createElement('Component')
332 component.setAttribute('rank', str(rank))
333
334 sim.toDom(document, component)
335
336 simulation.appendChild(component)
337
338 node.appendChild(simulation)
339
340 def writeXML(self,ostream=stdout):
341 """writes the object as an XML object into an output stream"""
342 document, rootnode = esysDoc()
343 self.toDom(document, rootnode)
344 ostream.write(document.toprettyxml())
345
346 def getSafeTimeStepSize(self,dt):
347 """returns a time step size which can savely be used by all subsimulations"""
348 out=self.UNDEF_DT
349 for i in self.iterSubsimulations():
350 dt_new=i.getSafeTimeStepSize(dt)
351 if dt_new!=None: out=min(out,dt_new)
352 return out
353
354 def doInitialization(self,dt):
355 """initializes all subsimulations """
356 for i in self.iterSubsimulations(): i.doInitialization(dt)
357
358 def finalize(self):
359 """returns True if all subsimulations are finalized"""
360 out=True
361 for i in self.iterSubsimulations(): out = out and i.finalize()
362 return out
363
364 def doFinalization(self):
365 """finalalizes the time stepping for all subsimulations."""
366 for i in self.iterSubsimulations(): i.doFinalization()
367
368 def doIterationInitialization(self,dt):
369 """initializes the iteration at time t for all subsimulations."""
370 self.__iter=0
371 if self.debug(): print "%s: iteration starts"%self
372 for i in self.iterSubsimulations(): i.doIterationInitialization(dt)
373
374 def terminate(self):
375 """returns True if all iterations for all subsimulations are terminated."""
376 out=True
377 for i in self.iterSubsimulations(): out = out and i.terminate()
378 return out
379
380 def doIterationFinalization(self,dt):
381 """finalalizes the iteration process for each of the subsimulations."""
382 for i in self.iterSubsimulations(): i.doIterationFinalization(dt)
383 if self.debug(): print "%s: iteration finalized after %s steps"%(self,self.__iter+1)
384
385 def doIterationStep(self,dt):
386 """executes the iteration step at time step for each subsimulation"""
387 self.__iter+=1
388 if self.debug(): print "%s: iteration step %d"%(self,self.__iter)
389 for i in self.iterSubsimulations(): i.doIterationStep(dt)
390
391 class ExplicitSimulation(Simulation):
392 """This is a modified form of the Simulation class. In fact it overwrites
393 the doStep method by executing the doStep method of all subsimulations
394 rather then iterating over all subsimulations."""
395
396 def doStep(self,dt):
397 """executes the time step for all subsimulation"""
398 for i in self.iterSubsimulations(): i.doStep(dt)
399
400 class _ParameterIterator:
401 def __init__(self,parameterset):
402 self.__set=parameterset
403 self.__iter=iter(parameterset.getParameterList())
404 def next(self):
405 o=self.__iter.next()
406 return (o,self.__set.getAttributeObject(o))
407 def __iter__(self):
408 return self
409
410 class ParameterSet(LinkableObject):
411 """a class which allows to emphazise attributes to be written and read to XML
412
413 Leaves of an ESySParameters objects can be
414
415 a real number
416 a integer number
417 a string
418 a boolean value
419 a ParameterSet object
420 a Simulation object
421 a Model object
422 any other object (not considered by writeESySXML and writeXML)
423
424 Example how to create an ESySParameters object:
425
426 p11=ParameterSet(gamma1=1.,gamma2=2.,gamma3=3.)
427 p1=ParameterSet(dim=2,tol_v=0.001,output_file="/tmp/u.%3.3d.dx",runFlag=True,parm11=p11)
428 parm=ParameterSet(parm1=p1,parm2=ParameterSet(alpha=Link(p11,"gamma1")))
429
430 This can be accessed as
431
432 parm.parm1.gamma=0.
433 parm.parm1.dim=2
434 parm.parm1.tol_v=0.001
435 parm.parm1.output_file="/tmp/u.%3.3d.dx"
436 parm.parm1.runFlag=True
437 parm.parm1.parm11.gamma1=1.
438 parm.parm1.parm11.gamma2=2.
439 parm.parm1.parm11.gamma3=3.
440 parm.parm2.alpha=1. (value of parm.parm1.parm11.gamma1)
441
442 """
443 def __init__(self,parameters=[]):
444 """creates a ParameterSet with parameters parameters"""
445 LinkableObject.__init__(self,debug=False)
446 self.__parameters=[]
447 self.declareParameters(parameters)
448
449 def getParameterList(self):
450 """returns the list of parameters"""
451 return self.__parameters
452
453 def isParameter(self,name):
454 """returns true if name is a parameter"""
455 return name in self.getParameterList()
456
457 def declareParameter(self,**parameters):
458 """declares a new parameter(s) and its (their) inital value."""
459 self.declareParameters(parameters)
460
461 def declareParameters(self,parameters):
462 """declares a set of parameters. parameters can be a list, a dictonary or a ParameterSet."""
463 if isinstance(parameters,ListType):
464 for prm in parameters:
465 setattr(self,prm,None)
466 if not self.isParameter(prm): self.getParameterList().append(prm)
467 if self.debug(): print "%s: parameter %s has been declared."%(self,prm)
468 elif isinstance(parameters,DictType):
469 for prm,value in parameters.iteritems():
470 setattr(self,prm,value)
471 if not self.isParameter(prm): self.getParameterList().append(prm)
472 if self.debug(): print "%s: parameter %s has been declared."%(self,prm)
473 else:
474 for prm,value in parameters:
475 setattr(self,prm,value)
476 if not self.isParameter(prm): self.getParameterList().append(prm)
477 if self.debug(): print "%s: parameter %s has been declared."%(self,prm)
478
479 def releaseParameters(self,name):
480 """removes parameter name from the paramameters"""
481 if self.isParameter(name):
482 self.getParameterList().remove(name)
483 if self.debug(): print "%s: parameter %s has been removed."%(self,prm)
484
485 def __iter__(self):
486 """creates an iterator over the parameter and their values"""
487 return _ParameterIterator(self)
488
489 def showParameters(self):
490 """returns a descrition of the parameters"""
491 out="{"
492 notfirst=False
493 for i,v in self:
494 if notfirst: out=out+","
495 notfirst=True
496 if isinstance(v,ParameterSet):
497 out="%s\"%s\" : %s"%(out,i,v.showParameters())
498 else:
499 out="%s\"%s\" : %s"%(out,i,v)
500 return out+"}"
501
502 def __delattr__(self,name):
503 """removes the attribute name."""
504 LinkableObject.__delattr__(self,name)
505 try:
506 self.releaseParameter(name)
507 except:
508 pass
509
510 def toDom(self, document, node):
511 """ toDom method of ParameterSet class """
512 pset = document.createElement('ParameterSet')
513 node.appendChild(pset)
514 self._parametersToDom(document, pset)
515
516 def _parametersToDom(self, document, node):
517 for name,value in self:
518 param = document.createElement('Parameter')
519 param.setAttribute('type', value.__class__.__name__)
520
521 param.appendChild(dataNode(document, 'Name', name))
522
523 val = document.createElement('Value')
524
525 if isinstance(value,ParameterSet):
526 value.toDom(document, val)
527 elif isinstance(value, Link):
528 value.toDom(document, val)
529 param.appendChild(val)
530 elif isinstance(value,StringType):
531 param.appendChild(dataNode(document, 'Value', value))
532 else:
533 param.appendChild(dataNode(document, 'Value', str(value)))
534
535 node.appendChild(param)
536
537 def writeXML(self,ostream=stdout):
538 """writes the object as an XML object into an output stream"""
539 # ParameterSet(d) with d[Name]=Value
540 document, node = esysDoc()
541 self.toDom(document, node)
542 ostream.write(document.toprettyxml())
543
544 class Model(ParameterSet,SimulationFrame):
545 """a Model is a SimulationFrame which is also a ParameterSet."""
546
547 def __init__(self,debug=False):
548 """creates a model"""
549 ParameterSet.__init__(self)
550 SimulationFrame.__init__(self,debug=debug)
551
552 def writeXML(self,ostream=stdout,esysheader=True,modelheader=True):
553 """writes the object as an XML object into an output stream"""
554 # x=type() with x.Name=Value where type is a Model from the model library.
555 if esysheader: ostream.write(_HEADER_START)
556 ostream.write("<Model type=\"%s\" module=\"%s\">\n"%(self.__class__.__name__,self.__class__.__module__))
557 self. writeParameterXML(ostream)
558 ostream.write("</Model>\n")
559 if esysheader: ostream.write(_HEADER_END)
560
561 class IterationDivergenceError(Exception):
562 """excpetion which should be thrown if there is no convergence of the iteration process at a time step but there is a chance taht a smaller step could help
563 to reach convergence."""
564 pass
565
566 class IterationBreakDownError(Exception):
567 """excpetion which should be thrown if there is no conevregence and there is no chance that a time step reduction would help"""
568 pass
569
570 class FailedTimeStepError(Exception):
571 """excpetion which should be thrown if the time step fails because of a step size that have been choosen to be to large"""
572 pass
573
574 class NonPositiveStepSizeError(Exception):
575 """excpetion which is thrown if the step size is not positive"""
576 pass
577
578 #
579 # ignore this text:
580 #
581 """ the Model class provides a framework to run a time-dependent simulation. A
582 Model has a set of parameter which may be fixed or altered by the Model itself
583 or other Models over time.
584
585 The parameters of a models are declared at instantion, e.g.
586
587 m=Model({"message" : "none" })
588
589 creates a Model with parameters p1 and p2 with inital values 1 and 2.
590 Typically a particular model is defined as a subclass of Model:
591
592 class Messenger(Model):
593 def __init__(self):
594 Model.__init__(self,parameters={"message" : "none" })
595
596 m=MyModel()
597
598 There are various ways how model parameters can be changed:
599
600 1) use object attributes:
601
602 m.message="Hello World!"
603
604 2) use setParamter method
605
606
607 m.setParameters(message="Hello World!")
608
609 3) or dictonaries
610
611 d={ message : "Hello World!" }
612 m.setParameters(**d)
613
614
615 A model executed buy staring the run method of the model:
616
617 m=Messenger()
618 m.run()
619
620 The run methods marches through time. It first calls the
621 doInitialization() method of the Model to set up the process. In each
622 time step the doStep() method is called to get from the current to the
623 next time step. The step size is defined by calling the
624 getSafeTimeStepSize() method. the time integration process is
625 terminated when the finalize() methods return true. Final the
626 doFinalization() method is called to finalize the process. To implement
627 a particular model a subclass of the Model class is defined. The
628 subclass overwrites the default methods of Model.
629
630 The following class defines a messenger printing in the doStep method
631 what ever the current value of its parameter message is:
632
633 class Messenger(Model):
634 def __init__(self):
635 Model.__init__(self,parameters={"message" : "none" })
636
637 def doInitialization(self):
638 print "I start talking now!"
639
640 def doStep(self,t):
641 print "Message (time %e) : %s "%(t,self.message)
642
643 def doFinalization(self):
644 print "I have no more to say!"
645
646 If a instance of the Messenger class is run, it will print the
647 initialization and finalization message only. This is because the
648 default method for finalize() does always returns True and therefore the
649 transition is terminated startcht away.
650
651 Following example for solving the ODE using a forward euler scheme:
652
653 u(t=0)=u0
654 u_t=a*u**2 for all 0<t<=ten
655
656 exact solution is given by u(t)=1/(1/u0-a*t)
657
658 class Ode1(Model):
659 def __init__(self,**args):
660 Model.__init__(self,parameters={"tend" : 1., "dt" : 0.0001 ,"a" : 0.1 ,"u" : 1. },name="test",debug=True)
661
662 def doInitialization(self):
663 self._tn=0
664
665 def doStep(self,t):
666 self.u=self.u+(t-self._tn)*self.a*self.u**2
667 self._tn=t
668
669 def doFinalization(self):
670 print "all done final error = ",abs(self.u-1./(1./3.-self.a*self._tn))
671
672 def getSafeTimeStepSize(self):
673 return self.dt
674
675 def finalize(self):
676 return self._tn>=self.tend
677
678 In some cases at a given time step an iteration process has to be
679 performed to get the state of the Model for the next time step. ` In
680 this case the doStep() method is replaced by a sequance of methods which
681 implements this iterative process. The method then will control the
682 iteration process by initializing the iteration through calling the
683 doIterationInitialization() method. The iteration is preformed by
684 calling the doIterationStep() method until the terminate() method
685 returns True. The doIterationFinalization() method is called to end the
686 iteration.
687 For a particular model these methods have to overwritten by a suitable
688 subclass without touching the doStep() method.
689
690 following example is a modification of the example above. Here an
691 implicit euler scheme is used. in each time step the problem
692
693 0= u_{n+1}-u_{n}+a*dt*u_{n+1}**2
694
695 has to be solved for u_{n+1}. The Newton scheme is used to solve this non-linear problem.
696
697
698 class Ode2(Model):
699
700 def __init__(self,**args):
701 Model.__init__(self,{"tend" : 1., "dt" : 0.1 ,"a" : 10. ,"u" : 1. , "tol " : 1.e-8},"test","bla",None,True)
702
703 def doInitialization(self):
704 self.__tn=0
705
706 def doIterationInitialization(self,t):
707 self.__iter=0
708 self.u_last=self.u
709 self.current_dt=t-self.tn
710 self.__tn=t
711
712 def doIterationStep(self):
713 self.__iter+=1
714 self.u_old=self.u
715 self.u=(self.current_dt*self.a*self.u**2-self.u_last)/(2*self.current_dt*self.a*self.u-1.)
716
717 def terminate(self):
718 return abs(self.u_old-self.u)<self.tol*abs(self.u)
719
720 def doIterationFinalization(self)
721 print "all done"
722
723 def getSafeTimeStepSize(self):
724 return self.dt
725
726 def finalize(self):
727 return self.__tn>self.tend
728
729 A model can be composed from subsimulations. Subsimulations are treated
730 as model parameters. If a model parameter is set or a value of a model
731 parameter is requested, the model will search for this parameter its
732 subsimulations in the case the model does not have this parameter
733 itself. The order in which the subsimulations are searched is critical.
734 By default a Model initializes all its subsimulations, is finalized when
735 all its subsimulations are finalized and finalizes all its
736 subsimulations. In the case an iterative process is applied on a
737 particular time step the iteration is initialized for all
738 subsimulations, then the iteration step is performed for each
739 subsimulation until all subsimulations indicate termination. Then the
740 iteration is finalized for all subsimulations. Finally teh doStop()
741 method for all submethods is called.
742
743 Here we are creating a model which groups ab instantiation of the Ode2 and the Messenger Model
744
745 o=Ode2()
746 m=Messenger()
747 om=Model(subsimulations=[o,m],debug=True)
748 om.dt=0.01
749 om.u=1.
750 m.message="it's me!"
751 om.run()
752
753 Notice that dt and u are parameters of class Ode2 and message is a
754 parameter of the Messenger class. The Model formed from these models
755 automatically hand the assignment of new values down to the
756 subsimulation. om.run() starts this combined model where now the
757 soStep() method of the Messenger object printing the value of its
758 parameter message together with a time stamp is executed in each time
759 step introduced by the Ode2 model.
760
761 A parameter of a Model can be linked to an attribute of onother object,
762 typically an parameter of another Model object.
763
764
765 which is comprised by a set of subsimulations.
766 The simulation is run through its run method which in the simplest case has the form:
767
768 s=Model()
769 s.run()
770
771 The run has an initializion and finalization phase. The latter is called
772 if all subsimulations are to be finalized. The simulation is processing
773 in time through calling the stepForward methods which updates the
774 observables of each subsimulation. A time steps size which is save for
775 all subsimulation is choosen.
776
777 At given time step an iterative process may be performed to make sure
778 that all observables are consistent across all subsimulations. In this
779 case, similar the time dependence, an initialization and finalization of
780 the iteration is performed.
781
782 A Model has input and output parameters where each input parameter can
783 be constant, time dependent or may depend on an output parameter of
784 another model or the model itself. To create a parameter name of a model
785 and to assign a value to it one can use the statement
786
787 model.name=object
788
789
790 At any time the current value of the parameter name can be obtained by
791
792 value=model.name
793
794 If the object that has been assigned to the paramter/attribute name has
795 the attribute/parameter name isself the current value of this attribute
796 of the object is returned (e.g. for model.name=object where object has
797 an attribute name, the statement value=model.name whould assign the
798 value object.name to value.). If the name of the parameters of a model
799 and an object don't match the setParameter method of model can be used.
800 So
801
802 model.setParameter(name,object,name_for_object)
803
804 links the parameter name of model with the parameter name_for_object of
805 object.
806
807 The run method initiates checkpointing (it is not clear how to do this
808 yet)
809 =====
810
811 """
812
813
814
815 if __name__=="__main__":
816 import math
817 #
818 # test for parameter set
819 #
820 p11=ParameterSet()
821 p11.declareParameter(gamma1=1.,gamma2=2.,gamma3=3.)
822 p1=ParameterSet()
823 p1.declareParameter(dim=2,tol_v=0.001,output_file="/tmp/u.%3.3d.dx",runFlag=True,parm11=p11)
824 parm=ParameterSet({ "parm1" : p1 , "parm2" : ParameterSet(["alpha"])})
825 parm.parm2.alpha=Link(p11,"gamma1")
826 parm.x="that should not be here!"
827 print parm.showParameters()
828 # should be something like: {"parm2" : {"alpha" : reference to attribute
829 # gamma1 of <__main__.ParameterSet instance at 0xf6db51cc>},"parm1" : {"dim"
830 # : 2,"runFlag" : True,"tol_v": 0.001,"parm11" : {"gamma3" : 3.0,"gamma2" :
831 # 2.0,"gamma1" : 1.0},"output_file" : /tmp/u.%3.3d.dx}}
832 assert parm.parm2.alpha==1.
833 parm.writeXML()
834
835 #=======================
836 class Messenger(Model):
837 def __init__(self):
838 Model.__init__(self)
839 self.declareParameter(message="none")
840
841 def doInitialization(self,t):
842 self.__t=t
843 print "I start talking now!"
844
845 def doStep(self,dt):
846 self.__t+=dt
847 print "Message (time %e) : %s "%(self.__t,self.message)
848
849 def doFinalization(self):
850 print "I have no more to say!"
851
852 class ODETEST(Model):
853 """ implements a solver for the ODE
854
855 du/dt=a*u+f(t)
856
857 we use a implicit euler scheme :
858
859 u_n-u_{n-1}= dt*a u_n + st*f(t_n)
860
861 to get u_n we run an iterative process
862
863 u_{n.k}=u_{n-1}+dt*(a u_{n.i-1} + f(t_n))
864
865
866 input for this model are step size dt, end time tend and a value for
867 a, f and initial value for u. we need also a tolerance tol for a
868 stopping criterion.
869
870 """
871
872 def __init__(self):
873 Model.__init__(self,debug=True)
874 self.declareParameter(tend=1.,dt=0.1,a=0.9,u=10.,f=0.,message="",tol=1.e-8)
875
876 def doInitialization(self,t):
877 self.__tn=t
878 self.__iter=0
879
880 def doIterationInitialization(self,dt):
881 self.__iter=0
882 self.__u_last=self.u
883
884 def doIterationStep(self,dt):
885 self.__iter+=1
886 self.__u_old=self.u
887 self.u=self.__u_last+dt*(self.a*self.__u_old+self.f)
888
889 def terminate(self):
890 if self.__iter<1:
891 return False
892 else:
893 return abs(self.__u_old-self.u)<self.tol*abs(self.u)
894
895 def doIterationFinalization(self,dt):
896 self.__tn+=dt
897 self.message="current error = %e"%abs(self.u-10.*math.exp((self.a-1.)*self.__tn))
898
899 def getSafeTimeStepSize(self,dt):
900 return min(self.dt,1./(abs(self.a)+1.))
901
902 def finalize(self):
903 return self.__tn>=self.tend
904
905 #
906 # s solves the coupled ODE:
907 #
908 # du/dt=a*u+ v
909 # dv/dt= u+a*v
910 #
911 # each equation is treated through the ODETEST class. The equations are
912 # linked and iteration over each time step is performed. the current
913 # error of v is reported by the Messenger class.
914 #
915 o1=ODETEST()
916 o1.u=10
917 o2=ODETEST()
918 o2.u=-10.
919 o1.f=Link(o2,"u")
920 o2.f=Link(o1,"u")
921 m=Messenger()
922 o1.dt=0.01
923 m.message=Link(o1)
924 s=ExplicitSimulation([Simulation([o1,o2],debug=True),m],debug=True)
925 s.run()
926 s.writeXML()

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26