/[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 126 - (show annotations)
Fri Jul 22 03:53:08 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: 38044 byte(s)
Merge of development branch back to main trunk on 2005-07-22

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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26