/[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 142 - (show annotations)
Mon Jul 25 05:28:20 2005 UTC (14 years, 3 months ago) by jgs
Original Path: trunk/esys2/escript/py_src/modelframe.py
File MIME type: text/x-python
File size: 38069 byte(s)
Merge of development branch back to main trunk on 2005-07-25

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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26