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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 874 - (hide annotations)
Tue Oct 17 12:06:11 2006 UTC (13 years ago) by elspeth
File MIME type: text/x-python
File size: 31714 byte(s)
Changed numarray xml structure slightly, to fit better with the rest of 
the XML. 
Added lists of booleans as a data type in modelframe.
Fixed bug with bools in modelframe.

1 jgs 121 # $Id$
2    
3 gross 637 """
4     Environment for implementing models in escript
5    
6     @var __author__: name of author
7     @var __copyright__: copyrights
8     @var __license__: licence agreement
9     @var __url__: url entry point on documentation
10     @var __version__: version
11     @var __date__: date of the version
12     """
13    
14     __author__="Lutz Gross, l.gross@uq.edu.au"
15 elspeth 609 __copyright__=""" Copyright (c) 2006 by ACcESS MNRF
16     http://www.access.edu.au
17     Primary Business: Queensland, Australia"""
18 elspeth 614 __license__="""Licensed under the Open Software License version 3.0
19     http://www.opensource.org/licenses/osl-3.0.php"""
20 gross 637 __url__="http://www.iservo.edu.au/esys"
21     __version__="$Revision$"
22     __date__="$Date$"
23 elspeth 609
24 gross 637
25 jgs 122 from types import StringType,IntType,FloatType,BooleanType,ListType,DictType
26     from sys import stdout
27 elspeth 871 import numarray
28     import operator
29 jgs 123 import itertools
30 jgs 142 # import modellib temporarily removed!!!
31 jgs 121
32 jgs 123 # import the 'set' module if it's not defined (python2.3/2.4 difference)
33     try:
34     set
35     except NameError:
36     from sets import Set as set
37    
38 jgs 122 from xml.dom import minidom
39 jgs 121
40 jgs 122 def dataNode(document, tagName, data):
41 jgs 126 """
42 jgs 149 C{dataNode}s are the building blocks of the xml documents constructed in
43     this module.
44    
45     @param document: the current xml document
46     @param tagName: the associated xml tag
47     @param data: the values in the tag
48 jgs 126 """
49 jgs 123 t = document.createTextNode(str(data))
50     n = document.createElement(tagName)
51     n.appendChild(t)
52     return n
53 jgs 121
54 jgs 122 def esysDoc():
55 jgs 126 """
56     Global method for creating an instance of an EsysXML document.
57     """
58 jgs 123 doc = minidom.Document()
59     esys = doc.createElement('ESys')
60     doc.appendChild(esys)
61     return doc, esys
62 jgs 121
63 jgs 123 def all(seq):
64     for x in seq:
65     if not x:
66     return False
67     return True
68    
69 jgs 147 def any(seq):
70     for x in seq:
71     if x:
72     return True
73     return False
74    
75 jgs 123 LinkableObjectRegistry = {}
76    
77     def registerLinkableObject(obj_id, o):
78     LinkableObjectRegistry[obj_id] = o
79    
80     LinkRegistry = []
81    
82     def registerLink(obj_id, l):
83     LinkRegistry.append((obj_id,l))
84    
85     def parse(xml):
86 jgs 147 """
87 jgs 149 Generic parse method for EsysXML. Without this, Links don't work.
88 jgs 147 """
89 jgs 123 global LinkRegistry, LinkableObjectRegistry
90     LinkRegistry = []
91     LinkableObjectRegistry = {}
92    
93     doc = minidom.parseString(xml)
94     sim = getComponent(doc.firstChild)
95     for obj_id, link in LinkRegistry:
96     link.target = LinkableObjectRegistry[obj_id]
97    
98     return sim
99    
100 jgs 150 def importName(modulename, name):
101     """ Import a named object from a module in the context of this function,
102     which means you should use fully qualified module paths.
103    
104     Return None on failure.
105    
106     This function from: http://aspn.activestate.com/ASPN/Cookbook/Python/Recipe/52241
107     """
108     module = __import__(modulename, globals(), locals(), [name])
109    
110     try:
111     return vars(module)[name]
112     except KeyError:
113     raise ImportError("Could not import %s from %s" % (name, modulename))
114    
115 jgs 123 def getComponent(doc):
116 jgs 147 """
117     Used to get components of Simualtions, Models.
118     """
119 jgs 123 for node in doc.childNodes:
120 jgs 147
121 jgs 123 if isinstance(node, minidom.Element):
122     if node.tagName == 'Simulation':
123     if node.getAttribute("type") == 'Simulation':
124     return Simulation.fromDom(node)
125     if node.tagName == 'Model':
126 jgs 150 if (node.getAttribute("module")):
127     model_module = node.getAttribute("module")
128     model_type = node.getAttribute("type")
129     return importName(model_module, model_type).fromDom(node)
130     else:
131     model_type = node.getAttribute("type")
132     model_subclasses = Model.__subclasses__()
133     for model in model_subclasses:
134     if model_type == model.__name__:
135     return Model.fromDom(node)
136 jgs 147 if node.tagName == 'ParameterSet':
137     parameter_type = node.getAttribute("type")
138     return ParameterSet.fromDom(node)
139 jgs 123 raise "Invalid simulation type, %r" % node.getAttribute("type")
140 jgs 147
141 jgs 123
142     raise ValueError("No Simulation Found")
143    
144    
145 jgs 122 class Link:
146 jgs 123 """
147 jgs 149 A Link makes an attribute of an object callable::
148    
149 jgs 123 o.object()
150     o.a=8
151     l=Link(o,"a")
152     assert l()==8
153     """
154    
155     def __init__(self,target,attribute=None):
156 jgs 126 """
157 jgs 149 Creates a link to the object target. If attribute is given, the link is
158 jgs 123 establised to this attribute of the target. Otherwise the attribute is
159 jgs 126 undefined.
160     """
161 jgs 123 self.target = target
162     self.attribute = None
163     self.setAttributeName(attribute)
164    
165     def setAttributeName(self,attribute):
166 jgs 126 """
167 jgs 149 Set a new attribute name to be collected from the target object. The
168 jgs 126 target object must have the attribute with name attribute.
169     """
170 jgs 147 if attribute and self.target:
171     if isinstance(self.target,LinkableObject):
172     if not self.target.hasAttribute(attribute):
173     raise AttributeError("%s: target %s has no attribute %s."%(self, self.target, attribute))
174     else:
175     if not hasattr(self.target,attribute):
176     raise AttributeError("%s: target %s has no attribute %s."%(self, self.target, attribute))
177 jgs 123 self.attribute = attribute
178    
179     def hasDefinedAttributeName(self):
180 jgs 126 """
181 jgs 149 Returns true if an attribute name is set.
182 jgs 126 """
183 jgs 123 return self.attribute != None
184    
185     def __repr__(self):
186 jgs 126 """
187 jgs 149 Returns a string representation of the link.
188 jgs 126 """
189 jgs 123 if self.hasDefinedAttributeName():
190     return "<Link to attribute %s of %s>" % (self.attribute, self.target)
191     else:
192     return "<Link to target %s>" % self.target
193    
194     def __call__(self,name=None):
195 jgs 126 """
196 jgs 149 Returns the value of the attribute of the target object. If the
197 jgs 126 atrribute is callable then the return value of the call is returned.
198     """
199 jgs 123 if name:
200     out=getattr(self.target, name)
201     else:
202     out=getattr(self.target, self.attribute)
203 jgs 121
204 jgs 123 if callable(out):
205     return out()
206     else:
207     return out
208 jgs 121
209 jgs 123 def toDom(self, document, node):
210 jgs 126 """
211 jgs 149 C{toDom} method of Link. Creates a Link node and appends it to the
212     current XML document.
213 jgs 126 """
214 jgs 123 link = document.createElement('Link')
215 jgs 147 assert (self.target != None), ("Target was none, name was %r" % self.attribute)
216 jgs 123 link.appendChild(dataNode(document, 'Target', self.target.id))
217     # this use of id will not work for purposes of being able to retrieve the intended
218     # target from the xml later. I need a better unique identifier.
219     assert self.attribute, "You can't xmlify a Link without a target attribute"
220     link.appendChild(dataNode(document, 'Attribute', self.attribute))
221     node.appendChild(link)
222 jgs 121
223 jgs 123 def fromDom(cls, doc):
224     targetid = doc.getElementsByTagName("Target")[0].firstChild.nodeValue.strip()
225     attribute = doc.getElementsByTagName("Attribute")[0].firstChild.nodeValue.strip()
226     l = cls(None, attribute)
227     registerLink(targetid, l)
228     return l
229 jgs 121
230 jgs 123 fromDom = classmethod(fromDom)
231    
232     def writeXML(self,ostream=stdout):
233 jgs 126 """
234 jgs 149 Writes an XML representation of self to the output stream ostream.
235 jgs 123 If ostream is nor present the standart output stream is used. If
236 jgs 126 esysheader==True the esys XML header is written
237     """
238 jgs 147 print 'I got to the Link writeXML method'
239 jgs 123 document, rootnode = esysDoc()
240     self.toDom(document, rootnode)
241    
242     ostream.write(document.toprettyxml())
243    
244 jgs 122 class LinkableObject(object):
245 jgs 123 """
246     An object that allows to link its attributes to attributes of other objects
247 jgs 149 via a Link object. For instance::
248 jgs 123
249     p = LinkableObject()
250     p.x = Link(o,"name")
251     print p.x
252    
253 jgs 149 links attribute C{x} of C{p} to the attribute name of object C{o}.
254 jgs 121
255 jgs 149 C{p.x} will contain the current value of attribute C{name} of object
256     C{o}.
257 jgs 121
258 jgs 149 If the value of C{getattr(o, "name")} is callable, C{p.x} will return
259     the return value of the call.
260 jgs 123 """
261    
262     number_sequence = itertools.count(100)
263    
264     def __init__(self, debug=False):
265 jgs 149 """
266     Initializes LinkableObject so that we can operate on Links
267     """
268 jgs 123 self.debug = debug
269     self.__linked_attributes={}
270     self.id = self.number_sequence.next()
271 jgs 149 registerLinkableObject(self.id, self)
272 jgs 123
273     def trace(self, msg):
274     """
275 jgs 149 If debugging is on, print the message, otherwise do nothing
276     """
277 jgs 123 if self.debug:
278 jgs 147 print "%s: %s"%(str(self),msg)
279 jgs 123
280     def __getattr__(self,name):
281 jgs 149 """
282     Returns the value of attribute name. If the value is a Link object the
283     object is called and the return value is returned.
284     """
285 jgs 123 out = self.getAttributeObject(name)
286     if isinstance(out,Link):
287     return out()
288     else:
289     return out
290    
291     def getAttributeObject(self,name):
292 jgs 149 """
293     Return the object stored for attribute name.
294     """
295 jgs 123
296     if self.__dict__.has_key(name):
297     return self.__dict__[name]
298    
299     if self.__linked_attributes.has_key(name):
300     return self.__linked_attributes[name]
301    
302 jgs 147 if self.__class__.__dict__.has_key(name):
303     return self.__class.__dict__[name]
304    
305 jgs 123 raise AttributeError,"No attribute %s."%name
306    
307 jgs 147 def hasAttribute(self,name):
308 jgs 149 """
309     Returns True if self as attribute name.
310     """
311 jgs 147 return self.__dict__.has_key(name) or self.__linked_attributes.has_key(name) or self.__class__.__dict__.has_key(name)
312    
313 jgs 123 def __setattr__(self,name,value):
314 jgs 149 """
315     Sets the value for attribute name. If value is a Link the target
316     attribute is set to name if no attribute has been specified.
317     """
318 jgs 123
319     if self.__dict__.has_key(name):
320     del self.__dict__[name]
321    
322     if isinstance(value,Link):
323     if not value.hasDefinedAttributeName():
324     value.setAttributeName(name)
325     self.__linked_attributes[name] = value
326    
327 jgs 147 self.trace("attribute %s is now linked by %s."%(name,value))
328 jgs 123 else:
329     self.__dict__[name] = value
330    
331     def __delattr__(self,name):
332 jgs 149 """
333     Removes the attribute name.
334     """
335 jgs 123
336     if self.__linked_attributes.has_key[name]:
337     del self.__linked_attributes[name]
338     elif self.__dict__.has_key(name):
339     del self.__dict__[name]
340     else:
341     raise AttributeError,"No attribute %s."%name
342    
343 jgs 122 class _ParameterIterator:
344 jgs 123 def __init__(self,parameterset):
345 jgs 121
346 jgs 123 self.__set=parameterset
347     self.__iter=iter(parameterset.parameters)
348    
349     def next(self):
350     o=self.__iter.next()
351     return (o,self.__set.getAttributeObject(o))
352    
353     def __iter__(self):
354     return self
355    
356 jgs 122 class ParameterSet(LinkableObject):
357 jgs 149 """
358     A class which allows to emphazise attributes to be written and read to XML
359 jgs 123
360 jgs 149 Leaves of an ESySParameters object can be:
361 jgs 123
362 jgs 149 - a real number
363     - a integer number
364     - a string
365     - a boolean value
366     - a ParameterSet object
367     - a Simulation object
368     - a Model object
369 elspeth 874 - a numarray object
370     - a list of booleans
371     - any other object (not considered by writeESySXML and writeXML)
372 jgs 123
373 jgs 149 Example how to create an ESySParameters object::
374 jgs 123
375 jgs 149 p11=ParameterSet(gamma1=1.,gamma2=2.,gamma3=3.)
376     p1=ParameterSet(dim=2,tol_v=0.001,output_file="/tmp/u.%3.3d.dx",runFlag=True,parm11=p11)
377     parm=ParameterSet(parm1=p1,parm2=ParameterSet(alpha=Link(p11,"gamma1")))
378 jgs 123
379 jgs 149 This can be accessed as::
380 jgs 123
381 jgs 149 parm.parm1.gamma=0.
382     parm.parm1.dim=2
383     parm.parm1.tol_v=0.001
384     parm.parm1.output_file="/tmp/u.%3.3d.dx"
385     parm.parm1.runFlag=True
386     parm.parm1.parm11.gamma1=1.
387     parm.parm1.parm11.gamma2=2.
388     parm.parm1.parm11.gamma3=3.
389     parm.parm2.alpha=1. (value of parm.parm1.parm11.gamma1)
390 jgs 123 """
391     def __init__(self, parameters=[], **kwargs):
392 jgs 149 """
393     Creates a ParameterSet with parameters parameters.
394     """
395 jgs 123 LinkableObject.__init__(self, **kwargs)
396     self.parameters = set()
397     self.declareParameters(parameters)
398 jgs 147
399     def __repr__(self):
400     return "<%s %r>" % (self.__class__.__name__,
401     [(p, getattr(self, p, None)) for p in self.parameters])
402 jgs 123
403     def declareParameter(self,**parameters):
404 jgs 149 """
405     Declares a new parameter(s) and its (their) initial value.
406     """
407 jgs 123 self.declareParameters(parameters)
408    
409     def declareParameters(self,parameters):
410 jgs 149 """
411     Declares a set of parameters. parameters can be a list, a dictionary
412     or a ParameterSet.
413     """
414 jgs 123 if isinstance(parameters,ListType):
415     parameters = zip(parameters, itertools.repeat(None))
416     if isinstance(parameters,DictType):
417     parameters = parameters.iteritems()
418 jgs 121
419 jgs 123 for prm, value in parameters:
420     setattr(self,prm,value)
421     self.parameters.add(prm)
422 jgs 121
423 jgs 147 self.trace("parameter %s has been declared."%prm)
424 jgs 121
425 jgs 123 def releaseParameters(self,name):
426 jgs 149 """
427     Removes parameter name from the paramameters.
428     """
429 jgs 123 if self.isParameter(name):
430     self.parameters.remove(name)
431 jgs 147 self.trace("parameter %s has been removed."%name)
432 jgs 123
433     def __iter__(self):
434 jgs 149 """
435     Creates an iterator over the parameter and their values.
436     """
437 jgs 123 return _ParameterIterator(self)
438    
439     def showParameters(self):
440 jgs 149 """
441     Returns a descrition of the parameters.
442     """
443 jgs 123 out="{"
444     notfirst=False
445     for i,v in self:
446     if notfirst: out=out+","
447     notfirst=True
448     if isinstance(v,ParameterSet):
449     out="%s\"%s\" : %s"%(out,i,v.showParameters())
450     else:
451     out="%s\"%s\" : %s"%(out,i,v)
452     return out+"}"
453    
454     def __delattr__(self,name):
455 jgs 149 """
456     Removes the attribute name.
457     """
458 jgs 123 LinkableObject.__delattr__(self,name)
459     try:
460     self.releaseParameter(name)
461     except:
462     pass
463 jgs 121
464 jgs 123 def toDom(self, document, node):
465 jgs 149 """
466     C{toDom} method of ParameterSet class.
467     """
468 jgs 123 pset = document.createElement('ParameterSet')
469     node.appendChild(pset)
470     self._parametersToDom(document, pset)
471 jgs 121
472 jgs 123 def _parametersToDom(self, document, node):
473 elspeth 282 node.setAttribute('id', str(self.id))
474     node.setIdAttribute("id")
475 jgs 123 for name,value in self:
476     param = document.createElement('Parameter')
477     param.setAttribute('type', value.__class__.__name__)
478 jgs 121
479 jgs 123 param.appendChild(dataNode(document, 'Name', name))
480 jgs 121
481 jgs 123 val = document.createElement('Value')
482    
483 elspeth 874 if isinstance(value,(ParameterSet,Link)):
484 jgs 123 value.toDom(document, val)
485     param.appendChild(val)
486 elspeth 871 elif isinstance(value, numarray.NumArray):
487     shape = value.getshape()
488     if isinstance(shape, tuple):
489     size = reduce(operator.mul, shape)
490     shape = ' '.join(map(str, shape))
491     else:
492     size = shape
493     shape = str(shape)
494    
495     arraytype = value.type()
496 elspeth 874 numarrayElement = document.createElement('NumArray')
497     numarrayElement.appendChild(dataNode(document, 'ArrayType', str(arraytype)))
498     numarrayElement.appendChild(dataNode(document, 'Shape', shape))
499     numarrayElement.appendChild(dataNode(document, 'Data', ' '.join(
500 elspeth 871 [str(x) for x in numarray.reshape(value, size)])))
501 elspeth 874 val.appendChild(numarrayElement)
502 elspeth 871 param.appendChild(val)
503 elspeth 874 elif isinstance (value, list):
504     param.appendChild(dataNode(document, 'Value', ' '.join(
505     [str(x) for x in value])
506     ))
507 jgs 123 else:
508     param.appendChild(dataNode(document, 'Value', str(value)))
509    
510     node.appendChild(param)
511    
512     def fromDom(cls, doc):
513    
514     # Define a host of helper functions to assist us.
515     def _children(node):
516     """
517 jgs 149 Remove the empty nodes from the children of this node.
518 jgs 123 """
519 elspeth 871 ret = []
520     for x in node.childNodes:
521     if isinstance(x, minidom.Text):
522     if x.nodeValue.strip():
523     ret.append(x)
524     else:
525     ret.append(x)
526     return ret
527 jgs 123
528     def _floatfromValue(doc):
529     return float(doc.nodeValue.strip())
530    
531     def _stringfromValue(doc):
532     return str(doc.nodeValue.strip())
533 jgs 126
534     def _intfromValue(doc):
535     return int(doc.nodeValue.strip())
536    
537     def _boolfromValue(doc):
538 elspeth 874 return _boolfromstring(doc.nodeValue.strip())
539 elspeth 266
540     def _nonefromValue(doc):
541     return None
542 elspeth 871
543     def _numarrayfromValue(doc):
544 elspeth 874 for node in _children(doc):
545 elspeth 871 if node.tagName == 'ArrayType':
546     arraytype = node.firstChild.nodeValue.strip()
547     if node.tagName == 'Shape':
548     shape = node.firstChild.nodeValue.strip()
549     shape = [int(x) for x in shape.split()]
550     if node.tagName == 'Data':
551     data = node.firstChild.nodeValue.strip()
552     data = [float(x) for x in data.split()]
553     return numarray.reshape(numarray.array(data, type=getattr(numarray, arraytype)),
554     shape)
555 elspeth 874
556     def _listfromValue(doc):
557     return [_boolfromstring(x) for x in doc.nodeValue.split()]
558    
559    
560     def _boolfromstring(s):
561     if s == 'True':
562     return True
563     else:
564     return False
565 jgs 123 # Mapping from text types in the xml to methods used to process trees of that type
566     ptypemap = {"Simulation": Simulation.fromDom,
567     "Model":Model.fromDom,
568     "ParameterSet":ParameterSet.fromDom,
569     "Link":Link.fromDom,
570     "float":_floatfromValue,
571 jgs 126 "int":_intfromValue,
572 jgs 123 "str":_stringfromValue,
573 elspeth 269 "bool":_boolfromValue,
574 elspeth 874 "list":_listfromValue,
575     "NumArray":_numarrayfromValue,
576 elspeth 871 "NoneType":_nonefromValue,
577 jgs 123 }
578    
579 jgs 147 # print doc.toxml()
580    
581 jgs 123 parameters = {}
582     for node in _children(doc):
583     ptype = node.getAttribute("type")
584    
585     pname = pvalue = None
586     for childnode in _children(node):
587    
588     if childnode.tagName == "Name":
589     pname = childnode.firstChild.nodeValue.strip()
590    
591     if childnode.tagName == "Value":
592     nodes = _children(childnode)
593 elspeth 874 # if ptype == 'NumArray':
594     # pvalue = _numarrayfromValue(nodes)
595     # else:
596     pvalue = ptypemap[ptype](nodes[0])
597 jgs 123
598     parameters[pname] = pvalue
599    
600     # Create the instance of ParameterSet
601     o = cls()
602     o.declareParameters(parameters)
603     registerLinkableObject(doc.getAttribute("id"), o)
604     return o
605    
606     fromDom = classmethod(fromDom)
607    
608     def writeXML(self,ostream=stdout):
609 jgs 149 """
610     Writes the object as an XML object into an output stream.
611     """
612 jgs 123 # ParameterSet(d) with d[Name]=Value
613     document, node = esysDoc()
614     self.toDom(document, node)
615     ostream.write(document.toprettyxml())
616    
617 jgs 147 class Model(ParameterSet):
618     """
619 jgs 149 A Model object represents a processess marching over time until a
620     finalizing condition is fullfilled. At each time step an iterative
621     process can be performed and the time step size can be controlled. A
622     Model has the following work flow::
623 jgs 121
624 jgs 147 doInitialization()
625     while not finalize():
626     dt=getSafeTimeStepSize(dt)
627     doStepPreprocessing(dt)
628     while not terminateIteration(): doStep(dt)
629     doStepPostprocessing(dt)
630     doFinalization()
631 jgs 121
632 jgs 149 where C{doInitialization}, C{finalize}, C{getSafeTimeStepSize},
633     C{doStepPreprocessing}, C{terminateIteration}, C{doStepPostprocessing},
634     C{doFinalization} are methods of the particular instance of a Model. The
635     default implementations of these methods have to be overwritten by the
636     subclass implementing a Model.
637 jgs 147 """
638 jgs 121
639 jgs 147 UNDEF_DT=1.e300
640 jgs 121
641 jgs 147 def __init__(self,parameters=[],**kwarg):
642 jgs 149 """
643     Creates a model.
644 jgs 121
645 jgs 149 Just calls the parent constructor.
646 jgs 147 """
647     ParameterSet.__init__(self, parameters=parameters,**kwarg)
648 jgs 121
649 jgs 147 def __str__(self):
650     return "<%s %d>"%(self.__class__,id(self))
651 jgs 121
652 jgs 147 def toDom(self, document, node):
653 jgs 149 """
654     C{toDom} method of Model class
655     """
656 jgs 147 pset = document.createElement('Model')
657     pset.setAttribute('type', self.__class__.__name__)
658 jgs 150 if not self.__class__.__module__.startswith('esys.escript'):
659     pset.setAttribute('module', self.__class__.__module__)
660 jgs 147 node.appendChild(pset)
661     self._parametersToDom(document, pset)
662 jgs 121
663 jgs 147 def doInitialization(self):
664 jgs 149 """
665     Initializes the time stepping scheme.
666    
667     This function may be overwritten.
668     """
669 jgs 147 pass
670    
671     def getSafeTimeStepSize(self,dt):
672 jgs 149 """
673     Returns a time step size which can safely be used.
674    
675     C{dt} gives the previously used step size.
676    
677     This function may be overwritten.
678     """
679 jgs 147 return self.UNDEF_DT
680    
681     def finalize(self):
682 jgs 149 """
683     Returns False if the time stepping is finalized.
684    
685     This function may be overwritten.
686     """
687 jgs 147 return False
688 jgs 123
689 jgs 147 def doFinalization(self):
690 jgs 149 """
691     Finalizes the time stepping.
692    
693     This function may be overwritten.
694     """
695 jgs 147 pass
696    
697     def doStepPreprocessing(self,dt):
698 jgs 149 """
699     Sets up a time step of step size dt.
700    
701     This function may be overwritten.
702     """
703 jgs 147 pass
704    
705     def doStep(self,dt):
706 jgs 149 """
707     Executes an iteration step at a time step.
708    
709     C{dt} is the currently used time step size.
710    
711     This function may be overwritten.
712     """
713 jgs 147 pass
714    
715     def terminateIteration(self):
716 jgs 149 """
717     Returns True if iteration on a time step is terminated.
718     """
719 jgs 147 return True
720    
721     def doStepPostprocessing(self,dt):
722 jgs 149 """
723     Finalalizes the time step.
724    
725     dt is the currently used time step size.
726    
727     This function may be overwritten.
728     """
729 jgs 147 pass
730    
731     def writeXML(self, ostream=stdout):
732     document, node = esysDoc()
733     self.toDom(document, node)
734     ostream.write(document.toprettyxml())
735    
736 jgs 121
737 jgs 147 class Simulation(Model):
738     """
739 jgs 149 A Simulation object is special Model which runs a sequence of Models.
740 jgs 121
741 jgs 149 The methods C{doInitialization}, C{finalize}, C{getSafeTimeStepSize},
742     C{doStepPreprocessing}, C{terminateIteration}, C{doStepPostprocessing},
743     C{doFinalization} are executing the corresponding methods of the models in
744     the simulation.
745 jgs 147 """
746    
747     FAILED_TIME_STEPS_MAX=20
748     MAX_ITER_STEPS=50
749 gross 836 MAX_CHANGE_OF_DT=2.
750 jgs 147
751     def __init__(self, models=[], **kwargs):
752 jgs 149 """
753     Initiates a simulation from a list of models.
754     """
755 jgs 147 Model.__init__(self, **kwargs)
756     self.__models=[]
757 jgs 149
758 jgs 147 for i in range(len(models)):
759     self[i] = models[i]
760 jgs 149
761 jgs 121
762 jgs 147 def __repr__(self):
763     """
764 jgs 149 Returns a string representation of the Simulation.
765 jgs 147 """
766     return "<Simulation %r>" % self.__models
767 jgs 121
768 jgs 147 def __str__(self):
769     """
770 jgs 149 Returning Simulation as a string.
771 jgs 147 """
772     return "<Simulation %d>"%id(self)
773    
774     def iterModels(self):
775 jgs 149 """
776     Returns an iterator over the models.
777     """
778 jgs 147 return self.__models
779    
780     def __getitem__(self,i):
781 jgs 149 """
782     Returns the i-th model.
783     """
784 jgs 147 return self.__models[i]
785    
786     def __setitem__(self,i,value):
787 jgs 149 """
788     Sets the i-th model.
789     """
790 jgs 147 if not isinstance(value,Model):
791 gross 728 raise ValueError,"assigned value is not a Model but instance of %s"%(value.__class__.__name__,)
792 jgs 149 for j in range(max(i-len(self.__models)+1,0)):
793     self.__models.append(None)
794 jgs 147 self.__models[i]=value
795    
796     def __len__(self):
797 jgs 149 """
798     Returns the number of models.
799     """
800 jgs 147 return len(self.__models)
801 jgs 121
802 jgs 147 def toDom(self, document, node):
803 jgs 149 """
804     C{toDom} method of Simulation class.
805     """
806 jgs 147 simulation = document.createElement('Simulation')
807     simulation.setAttribute('type', self.__class__.__name__)
808 jgs 121
809 jgs 147 for rank, sim in enumerate(self.iterModels()):
810     component = document.createElement('Component')
811     component.setAttribute('rank', str(rank))
812 jgs 121
813 jgs 147 sim.toDom(document, component)
814 jgs 121
815 jgs 147 simulation.appendChild(component)
816 jgs 121
817 jgs 147 node.appendChild(simulation)
818 jgs 121
819 jgs 147 def writeXML(self,ostream=stdout):
820 jgs 149 """
821     Writes the object as an XML object into an output stream.
822     """
823 jgs 147 document, rootnode = esysDoc()
824     self.toDom(document, rootnode)
825 jgs 149 targetsList = document.getElementsByTagName('Target')
826 elspeth 282
827     for element in targetsList:
828     targetId = int(element.firstChild.nodeValue.strip())
829     if document.getElementById(str(targetId)):
830     continue
831 jgs 149 targetObj = LinkableObjectRegistry[targetId]
832     targetObj.toDom(document, rootnode)
833 jgs 147 ostream.write(document.toprettyxml())
834    
835     def getSafeTimeStepSize(self,dt):
836 jgs 149 """
837     Returns a time step size which can safely be used by all models.
838    
839     This is the minimum over the time step sizes of all models.
840     """
841 jgs 147 out=min([o.getSafeTimeStepSize(dt) for o in self.iterModels()])
842 jgs 150 #print "%s: safe step size is %e."%(str(self),out)
843 jgs 147 return out
844    
845     def doInitialization(self):
846 jgs 149 """
847     Initializes all models.
848     """
849 jgs 147 self.n=0
850     self.tn=0.
851     for o in self.iterModels():
852     o.doInitialization()
853    
854     def finalize(self):
855 jgs 149 """
856     Returns True if any of the models is to be finalized.
857     """
858 jgs 147 return any([o.finalize() for o in self.iterModels()])
859 jgs 123
860 jgs 147 def doFinalization(self):
861 jgs 149 """
862     Finalalizes the time stepping for all models.
863     """
864 jgs 147 for i in self.iterModels(): i.doFinalization()
865     self.trace("end of time integation.")
866    
867     def doStepPreprocessing(self,dt):
868 jgs 149 """
869     Initializes the time step for all models.
870     """
871 jgs 147 for o in self.iterModels():
872     o.doStepPreprocessing(dt)
873    
874     def terminateIteration(self):
875 jgs 149 """
876     Returns True if all iterations for all models are terminated.
877     """
878 jgs 147 out=all([o.terminateIteration() for o in self.iterModels()])
879     return out
880 jgs 123
881 jgs 147 def doStepPostprocessing(self,dt):
882 jgs 149 """
883     Finalalizes the iteration process for all models.
884     """
885 jgs 147 for o in self.iterModels():
886     o.doStepPostprocessing(dt)
887     self.n+=1
888     self.tn+=dt
889    
890     def doStep(self,dt):
891     """
892 jgs 149 Executes the iteration step at a time step for all model::
893 jgs 147
894 jgs 149 self.doStepPreprocessing(dt)
895     while not self.terminateIteration():
896     for all models:
897     self.doStep(dt)
898     self.doStepPostprocessing(dt)
899 jgs 147 """
900     self.iter=0
901     while not self.terminateIteration():
902     if self.iter==0: self.trace("iteration at %d-th time step %e starts"%(self.n+1,self.tn+dt))
903     self.iter+=1
904     self.trace("iteration step %d"%(self.iter))
905     for o in self.iterModels():
906     o.doStep(dt)
907     if self.iter>0: self.trace("iteration at %d-th time step %e finalized."%(self.n+1,self.tn+dt))
908 jgs 121
909 jgs 147 def run(self,check_point=None):
910     """
911 jgs 149 Run the simulation by performing essentially::
912 jgs 123
913 jgs 149 self.doInitialization()
914     while not self.finalize():
915     dt=self.getSafeTimeStepSize()
916     self.doStep(dt)
917     if n%check_point==0:
918     self.writeXML()
919     self.doFinalization()
920 jgs 121
921 jgs 149 If one of the models in throws a C{FailedTimeStepError} exception a
922     new time step size is computed through getSafeTimeStepSize() and the
923     time step is repeated.
924 jgs 121
925 jgs 149 If one of the models in throws a C{IterationDivergenceError}
926     exception the time step size is halved and the time step is repeated.
927 jgs 121
928 jgs 149 In both cases the time integration is given up after
929     C{Simulation.FAILED_TIME_STEPS_MAX} attempts.
930 jgs 147 """
931     dt=self.UNDEF_DT
932     self.doInitialization()
933     while not self.finalize():
934     step_fail_counter=0
935     iteration_fail_counter=0
936 gross 836 if self.n==0:
937     dt_new=self.getSafeTimeStepSize(dt)
938     else:
939 gross 838 dt_new=min(max(self.getSafeTimeStepSize(dt),dt/self.MAX_CHANGE_OF_DT),dt*self.MAX_CHANGE_OF_DT)
940 jgs 147 self.trace("%d. time step %e (step size %e.)" % (self.n+1,self.tn+dt_new,dt_new))
941     end_of_step=False
942     while not end_of_step:
943     end_of_step=True
944     if not dt_new>0:
945 gross 829 raise NonPositiveStepSizeError("non-positive step size in step %d"%(self.n+1))
946 jgs 147 try:
947     self.doStepPreprocessing(dt_new)
948     self.doStep(dt_new)
949     self.doStepPostprocessing(dt_new)
950     except IterationDivergenceError:
951     dt_new*=0.5
952     end_of_step=False
953     iteration_fail_counter+=1
954     if iteration_fail_counter>self.FAILED_TIME_STEPS_MAX:
955 gross 836 raise SimulationBreakDownError("reduction of time step to achieve convergence failed after %s steps."%self.FAILED_TIME_STEPS_MAX)
956     self.trace("Iteration failed. Time step is repeated with new step size %s."%dt_new)
957 jgs 147 except FailedTimeStepError:
958     dt_new=self.getSafeTimeStepSize(dt)
959     end_of_step=False
960     step_fail_counter+=1
961 gross 836 self.trace("Time step is repeated with new time step size %s."%dt_new)
962 jgs 147 if step_fail_counter>self.FAILED_TIME_STEPS_MAX:
963 gross 836 raise SimulationBreakDownError("Time integration is given up after %d attempts."%step_fail_counter)
964 jgs 147 dt=dt_new
965     if not check_point==None:
966     if n%check_point==0:
967     self.trace("check point is created.")
968     self.writeXML()
969     self.doFinalization()
970 jgs 121
971 jgs 147 def fromDom(cls, doc):
972     sims = []
973     for node in doc.childNodes:
974     if isinstance(node, minidom.Text):
975     continue
976 jgs 121
977 jgs 147 sims.append(getComponent(node))
978 jgs 121
979 jgs 147 return cls(sims)
980 jgs 121
981 jgs 147 fromDom = classmethod(fromDom)
982 jgs 121
983    
984 jgs 147 class IterationDivergenceError(Exception):
985     """
986 jgs 149 Exception which is thrown if there is no convergence of the iteration
987     process at a time step.
988 jgs 121
989 jgs 149 But there is a chance that a smaller step could help to reach convergence.
990 jgs 147 """
991     pass
992 jgs 121
993 jgs 147 class FailedTimeStepError(Exception):
994 jgs 149 """
995     Exception which is thrown if the time step fails because of a step
996     size that have been choosen to be too large.
997     """
998 jgs 147 pass
999 jgs 121
1000 jgs 147 class SimulationBreakDownError(Exception):
1001 jgs 149 """
1002     Exception which is thrown if the simulation does not manage to
1003     progress in time.
1004     """
1005 jgs 147 pass
1006 jgs 121
1007 jgs 147 class NonPositiveStepSizeError(Exception):
1008 jgs 149 """
1009     Exception which is thrown if the step size is not positive.
1010     """
1011 jgs 147 pass
1012 jgs 149
1013     # vim: expandtab shiftwidth=4:

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26