/[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 871 - (hide annotations)
Sat Oct 14 08:25:54 2006 UTC (12 years, 6 months ago) by elspeth
File MIME type: text/x-python
File size: 31201 byte(s)
Numarray support added to modelframe and tests for numarray added to run_xml.
Numarrays will look like this in ESysXML:

A 2 x 3 numarray will look like this:

<Parameter type="NumArray">
                        <Name>
                                numtest
                        </Name>
                        <Value>
                                <ArrayType>
                                        Float64
                                </ArrayType>
                                <Shape>
                                        2 3
                                </Shape>
                                <Data>
                                        1.0 2.0 3.0 3.0 4.0 5.0
                                </Data>
                        </Value>
                </Parameter>

A numarray vector will look like this:

<Parameter type="NumArray">
                        <Name>
                                numtest
                        </Name>
                        <Value>
                                <ArrayType>
                                        Float64
                                </ArrayType>
                                <Shape>
                                        3
                                </Shape>
                                <Data>
                                        3.0 4.0 5.0
                                </Data>
                        </Value>
                </Parameter>



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     - any other object (not considered by writeESySXML and writeXML)
370 jgs 123
371 jgs 149 Example how to create an ESySParameters object::
372 jgs 123
373 jgs 149 p11=ParameterSet(gamma1=1.,gamma2=2.,gamma3=3.)
374     p1=ParameterSet(dim=2,tol_v=0.001,output_file="/tmp/u.%3.3d.dx",runFlag=True,parm11=p11)
375     parm=ParameterSet(parm1=p1,parm2=ParameterSet(alpha=Link(p11,"gamma1")))
376 jgs 123
377 jgs 149 This can be accessed as::
378 jgs 123
379 jgs 149 parm.parm1.gamma=0.
380     parm.parm1.dim=2
381     parm.parm1.tol_v=0.001
382     parm.parm1.output_file="/tmp/u.%3.3d.dx"
383     parm.parm1.runFlag=True
384     parm.parm1.parm11.gamma1=1.
385     parm.parm1.parm11.gamma2=2.
386     parm.parm1.parm11.gamma3=3.
387     parm.parm2.alpha=1. (value of parm.parm1.parm11.gamma1)
388 jgs 123 """
389     def __init__(self, parameters=[], **kwargs):
390 jgs 149 """
391     Creates a ParameterSet with parameters parameters.
392     """
393 jgs 123 LinkableObject.__init__(self, **kwargs)
394     self.parameters = set()
395     self.declareParameters(parameters)
396 jgs 147
397     def __repr__(self):
398     return "<%s %r>" % (self.__class__.__name__,
399     [(p, getattr(self, p, None)) for p in self.parameters])
400 jgs 123
401     def declareParameter(self,**parameters):
402 jgs 149 """
403     Declares a new parameter(s) and its (their) initial value.
404     """
405 jgs 123 self.declareParameters(parameters)
406    
407     def declareParameters(self,parameters):
408 jgs 149 """
409     Declares a set of parameters. parameters can be a list, a dictionary
410     or a ParameterSet.
411     """
412 jgs 123 if isinstance(parameters,ListType):
413     parameters = zip(parameters, itertools.repeat(None))
414     if isinstance(parameters,DictType):
415     parameters = parameters.iteritems()
416 jgs 121
417 jgs 123 for prm, value in parameters:
418     setattr(self,prm,value)
419     self.parameters.add(prm)
420 jgs 121
421 jgs 147 self.trace("parameter %s has been declared."%prm)
422 jgs 121
423 jgs 123 def releaseParameters(self,name):
424 jgs 149 """
425     Removes parameter name from the paramameters.
426     """
427 jgs 123 if self.isParameter(name):
428     self.parameters.remove(name)
429 jgs 147 self.trace("parameter %s has been removed."%name)
430 jgs 123
431     def __iter__(self):
432 jgs 149 """
433     Creates an iterator over the parameter and their values.
434     """
435 jgs 123 return _ParameterIterator(self)
436    
437     def showParameters(self):
438 jgs 149 """
439     Returns a descrition of the parameters.
440     """
441 jgs 123 out="{"
442     notfirst=False
443     for i,v in self:
444     if notfirst: out=out+","
445     notfirst=True
446     if isinstance(v,ParameterSet):
447     out="%s\"%s\" : %s"%(out,i,v.showParameters())
448     else:
449     out="%s\"%s\" : %s"%(out,i,v)
450     return out+"}"
451    
452     def __delattr__(self,name):
453 jgs 149 """
454     Removes the attribute name.
455     """
456 jgs 123 LinkableObject.__delattr__(self,name)
457     try:
458     self.releaseParameter(name)
459     except:
460     pass
461 jgs 121
462 jgs 123 def toDom(self, document, node):
463 jgs 149 """
464     C{toDom} method of ParameterSet class.
465     """
466 jgs 123 pset = document.createElement('ParameterSet')
467     node.appendChild(pset)
468     self._parametersToDom(document, pset)
469 jgs 121
470 jgs 123 def _parametersToDom(self, document, node):
471 elspeth 282 node.setAttribute('id', str(self.id))
472     node.setIdAttribute("id")
473 jgs 123 for name,value in self:
474     param = document.createElement('Parameter')
475     param.setAttribute('type', value.__class__.__name__)
476 jgs 121
477 jgs 123 param.appendChild(dataNode(document, 'Name', name))
478 jgs 121
479 jgs 123 val = document.createElement('Value')
480    
481     if isinstance(value,ParameterSet):
482     value.toDom(document, val)
483     param.appendChild(val)
484     elif isinstance(value, Link):
485     value.toDom(document, val)
486     param.appendChild(val)
487 elspeth 871 elif isinstance(value, numarray.NumArray):
488     shape = value.getshape()
489     if isinstance(shape, tuple):
490     size = reduce(operator.mul, shape)
491     shape = ' '.join(map(str, shape))
492     else:
493     size = shape
494     shape = str(shape)
495    
496     arraytype = value.type()
497     val.appendChild(dataNode(document, 'ArrayType', str(arraytype)))
498     val.appendChild(dataNode(document, 'Shape', shape))
499     val.appendChild(dataNode(document, 'Data', ' '.join(
500     [str(x) for x in numarray.reshape(value, size)])))
501     param.appendChild(val)
502 jgs 123 elif isinstance(value,StringType):
503     param.appendChild(dataNode(document, 'Value', value))
504     else:
505     param.appendChild(dataNode(document, 'Value', str(value)))
506    
507     node.appendChild(param)
508    
509     def fromDom(cls, doc):
510    
511     # Define a host of helper functions to assist us.
512     def _children(node):
513     """
514 jgs 149 Remove the empty nodes from the children of this node.
515 jgs 123 """
516 elspeth 871 ret = []
517     for x in node.childNodes:
518     if isinstance(x, minidom.Text):
519     if x.nodeValue.strip():
520     ret.append(x)
521     else:
522     ret.append(x)
523     return ret
524 jgs 123
525     def _floatfromValue(doc):
526     return float(doc.nodeValue.strip())
527    
528     def _stringfromValue(doc):
529     return str(doc.nodeValue.strip())
530 jgs 126
531     def _intfromValue(doc):
532     return int(doc.nodeValue.strip())
533    
534     def _boolfromValue(doc):
535     return bool(doc.nodeValue.strip())
536 elspeth 266
537     def _nonefromValue(doc):
538     return None
539 elspeth 871
540     def _numarrayfromValue(doc):
541     for node in doc:
542     if node.tagName == 'ArrayType':
543     arraytype = node.firstChild.nodeValue.strip()
544     if node.tagName == 'Shape':
545     shape = node.firstChild.nodeValue.strip()
546     shape = [int(x) for x in shape.split()]
547     if node.tagName == 'Data':
548     data = node.firstChild.nodeValue.strip()
549     data = [float(x) for x in data.split()]
550     return numarray.reshape(numarray.array(data, type=getattr(numarray, arraytype)),
551     shape)
552 jgs 126
553 jgs 123 # Mapping from text types in the xml to methods used to process trees of that type
554     ptypemap = {"Simulation": Simulation.fromDom,
555     "Model":Model.fromDom,
556     "ParameterSet":ParameterSet.fromDom,
557     "Link":Link.fromDom,
558     "float":_floatfromValue,
559 jgs 126 "int":_intfromValue,
560 jgs 123 "str":_stringfromValue,
561 elspeth 269 "bool":_boolfromValue,
562 elspeth 871 "NoneType":_nonefromValue,
563 jgs 123 }
564    
565 jgs 147 # print doc.toxml()
566    
567 jgs 123 parameters = {}
568     for node in _children(doc):
569     ptype = node.getAttribute("type")
570    
571     pname = pvalue = None
572     for childnode in _children(node):
573    
574     if childnode.tagName == "Name":
575     pname = childnode.firstChild.nodeValue.strip()
576    
577     if childnode.tagName == "Value":
578     nodes = _children(childnode)
579 elspeth 871 if ptype == 'NumArray':
580     pvalue = _numarrayfromValue(nodes)
581     else:
582     pvalue = ptypemap[ptype](nodes[0])
583 jgs 123
584     parameters[pname] = pvalue
585    
586     # Create the instance of ParameterSet
587     o = cls()
588     o.declareParameters(parameters)
589     registerLinkableObject(doc.getAttribute("id"), o)
590     return o
591    
592     fromDom = classmethod(fromDom)
593    
594     def writeXML(self,ostream=stdout):
595 jgs 149 """
596     Writes the object as an XML object into an output stream.
597     """
598 jgs 123 # ParameterSet(d) with d[Name]=Value
599     document, node = esysDoc()
600     self.toDom(document, node)
601     ostream.write(document.toprettyxml())
602    
603 jgs 147 class Model(ParameterSet):
604     """
605 jgs 149 A Model object represents a processess marching over time until a
606     finalizing condition is fullfilled. At each time step an iterative
607     process can be performed and the time step size can be controlled. A
608     Model has the following work flow::
609 jgs 121
610 jgs 147 doInitialization()
611     while not finalize():
612     dt=getSafeTimeStepSize(dt)
613     doStepPreprocessing(dt)
614     while not terminateIteration(): doStep(dt)
615     doStepPostprocessing(dt)
616     doFinalization()
617 jgs 121
618 jgs 149 where C{doInitialization}, C{finalize}, C{getSafeTimeStepSize},
619     C{doStepPreprocessing}, C{terminateIteration}, C{doStepPostprocessing},
620     C{doFinalization} are methods of the particular instance of a Model. The
621     default implementations of these methods have to be overwritten by the
622     subclass implementing a Model.
623 jgs 147 """
624 jgs 121
625 jgs 147 UNDEF_DT=1.e300
626 jgs 121
627 jgs 147 def __init__(self,parameters=[],**kwarg):
628 jgs 149 """
629     Creates a model.
630 jgs 121
631 jgs 149 Just calls the parent constructor.
632 jgs 147 """
633     ParameterSet.__init__(self, parameters=parameters,**kwarg)
634 jgs 121
635 jgs 147 def __str__(self):
636     return "<%s %d>"%(self.__class__,id(self))
637 jgs 121
638 jgs 147 def toDom(self, document, node):
639 jgs 149 """
640     C{toDom} method of Model class
641     """
642 jgs 147 pset = document.createElement('Model')
643     pset.setAttribute('type', self.__class__.__name__)
644 jgs 150 if not self.__class__.__module__.startswith('esys.escript'):
645     pset.setAttribute('module', self.__class__.__module__)
646 jgs 147 node.appendChild(pset)
647     self._parametersToDom(document, pset)
648 jgs 121
649 jgs 147 def doInitialization(self):
650 jgs 149 """
651     Initializes the time stepping scheme.
652    
653     This function may be overwritten.
654     """
655 jgs 147 pass
656    
657     def getSafeTimeStepSize(self,dt):
658 jgs 149 """
659     Returns a time step size which can safely be used.
660    
661     C{dt} gives the previously used step size.
662    
663     This function may be overwritten.
664     """
665 jgs 147 return self.UNDEF_DT
666    
667     def finalize(self):
668 jgs 149 """
669     Returns False if the time stepping is finalized.
670    
671     This function may be overwritten.
672     """
673 jgs 147 return False
674 jgs 123
675 jgs 147 def doFinalization(self):
676 jgs 149 """
677     Finalizes the time stepping.
678    
679     This function may be overwritten.
680     """
681 jgs 147 pass
682    
683     def doStepPreprocessing(self,dt):
684 jgs 149 """
685     Sets up a time step of step size dt.
686    
687     This function may be overwritten.
688     """
689 jgs 147 pass
690    
691     def doStep(self,dt):
692 jgs 149 """
693     Executes an iteration step at a time step.
694    
695     C{dt} is the currently used time step size.
696    
697     This function may be overwritten.
698     """
699 jgs 147 pass
700    
701     def terminateIteration(self):
702 jgs 149 """
703     Returns True if iteration on a time step is terminated.
704     """
705 jgs 147 return True
706    
707     def doStepPostprocessing(self,dt):
708 jgs 149 """
709     Finalalizes the time step.
710    
711     dt is the currently used time step size.
712    
713     This function may be overwritten.
714     """
715 jgs 147 pass
716    
717     def writeXML(self, ostream=stdout):
718     document, node = esysDoc()
719     self.toDom(document, node)
720     ostream.write(document.toprettyxml())
721    
722 jgs 121
723 jgs 147 class Simulation(Model):
724     """
725 jgs 149 A Simulation object is special Model which runs a sequence of Models.
726 jgs 121
727 jgs 149 The methods C{doInitialization}, C{finalize}, C{getSafeTimeStepSize},
728     C{doStepPreprocessing}, C{terminateIteration}, C{doStepPostprocessing},
729     C{doFinalization} are executing the corresponding methods of the models in
730     the simulation.
731 jgs 147 """
732    
733     FAILED_TIME_STEPS_MAX=20
734     MAX_ITER_STEPS=50
735 gross 836 MAX_CHANGE_OF_DT=2.
736 jgs 147
737     def __init__(self, models=[], **kwargs):
738 jgs 149 """
739     Initiates a simulation from a list of models.
740     """
741 jgs 147 Model.__init__(self, **kwargs)
742     self.__models=[]
743 jgs 149
744 jgs 147 for i in range(len(models)):
745     self[i] = models[i]
746 jgs 149
747 jgs 121
748 jgs 147 def __repr__(self):
749     """
750 jgs 149 Returns a string representation of the Simulation.
751 jgs 147 """
752     return "<Simulation %r>" % self.__models
753 jgs 121
754 jgs 147 def __str__(self):
755     """
756 jgs 149 Returning Simulation as a string.
757 jgs 147 """
758     return "<Simulation %d>"%id(self)
759    
760     def iterModels(self):
761 jgs 149 """
762     Returns an iterator over the models.
763     """
764 jgs 147 return self.__models
765    
766     def __getitem__(self,i):
767 jgs 149 """
768     Returns the i-th model.
769     """
770 jgs 147 return self.__models[i]
771    
772     def __setitem__(self,i,value):
773 jgs 149 """
774     Sets the i-th model.
775     """
776 jgs 147 if not isinstance(value,Model):
777 gross 728 raise ValueError,"assigned value is not a Model but instance of %s"%(value.__class__.__name__,)
778 jgs 149 for j in range(max(i-len(self.__models)+1,0)):
779     self.__models.append(None)
780 jgs 147 self.__models[i]=value
781    
782     def __len__(self):
783 jgs 149 """
784     Returns the number of models.
785     """
786 jgs 147 return len(self.__models)
787 jgs 121
788 jgs 147 def toDom(self, document, node):
789 jgs 149 """
790     C{toDom} method of Simulation class.
791     """
792 jgs 147 simulation = document.createElement('Simulation')
793     simulation.setAttribute('type', self.__class__.__name__)
794 jgs 121
795 jgs 147 for rank, sim in enumerate(self.iterModels()):
796     component = document.createElement('Component')
797     component.setAttribute('rank', str(rank))
798 jgs 121
799 jgs 147 sim.toDom(document, component)
800 jgs 121
801 jgs 147 simulation.appendChild(component)
802 jgs 121
803 jgs 147 node.appendChild(simulation)
804 jgs 121
805 jgs 147 def writeXML(self,ostream=stdout):
806 jgs 149 """
807     Writes the object as an XML object into an output stream.
808     """
809 jgs 147 document, rootnode = esysDoc()
810     self.toDom(document, rootnode)
811 jgs 149 targetsList = document.getElementsByTagName('Target')
812 elspeth 282
813     for element in targetsList:
814     targetId = int(element.firstChild.nodeValue.strip())
815     if document.getElementById(str(targetId)):
816     continue
817 jgs 149 targetObj = LinkableObjectRegistry[targetId]
818     targetObj.toDom(document, rootnode)
819 jgs 147 ostream.write(document.toprettyxml())
820    
821     def getSafeTimeStepSize(self,dt):
822 jgs 149 """
823     Returns a time step size which can safely be used by all models.
824    
825     This is the minimum over the time step sizes of all models.
826     """
827 jgs 147 out=min([o.getSafeTimeStepSize(dt) for o in self.iterModels()])
828 jgs 150 #print "%s: safe step size is %e."%(str(self),out)
829 jgs 147 return out
830    
831     def doInitialization(self):
832 jgs 149 """
833     Initializes all models.
834     """
835 jgs 147 self.n=0
836     self.tn=0.
837     for o in self.iterModels():
838     o.doInitialization()
839    
840     def finalize(self):
841 jgs 149 """
842     Returns True if any of the models is to be finalized.
843     """
844 jgs 147 return any([o.finalize() for o in self.iterModels()])
845 jgs 123
846 jgs 147 def doFinalization(self):
847 jgs 149 """
848     Finalalizes the time stepping for all models.
849     """
850 jgs 147 for i in self.iterModels(): i.doFinalization()
851     self.trace("end of time integation.")
852    
853     def doStepPreprocessing(self,dt):
854 jgs 149 """
855     Initializes the time step for all models.
856     """
857 jgs 147 for o in self.iterModels():
858     o.doStepPreprocessing(dt)
859    
860     def terminateIteration(self):
861 jgs 149 """
862     Returns True if all iterations for all models are terminated.
863     """
864 jgs 147 out=all([o.terminateIteration() for o in self.iterModels()])
865     return out
866 jgs 123
867 jgs 147 def doStepPostprocessing(self,dt):
868 jgs 149 """
869     Finalalizes the iteration process for all models.
870     """
871 jgs 147 for o in self.iterModels():
872     o.doStepPostprocessing(dt)
873     self.n+=1
874     self.tn+=dt
875    
876     def doStep(self,dt):
877     """
878 jgs 149 Executes the iteration step at a time step for all model::
879 jgs 147
880 jgs 149 self.doStepPreprocessing(dt)
881     while not self.terminateIteration():
882     for all models:
883     self.doStep(dt)
884     self.doStepPostprocessing(dt)
885 jgs 147 """
886     self.iter=0
887     while not self.terminateIteration():
888     if self.iter==0: self.trace("iteration at %d-th time step %e starts"%(self.n+1,self.tn+dt))
889     self.iter+=1
890     self.trace("iteration step %d"%(self.iter))
891     for o in self.iterModels():
892     o.doStep(dt)
893     if self.iter>0: self.trace("iteration at %d-th time step %e finalized."%(self.n+1,self.tn+dt))
894 jgs 121
895 jgs 147 def run(self,check_point=None):
896     """
897 jgs 149 Run the simulation by performing essentially::
898 jgs 123
899 jgs 149 self.doInitialization()
900     while not self.finalize():
901     dt=self.getSafeTimeStepSize()
902     self.doStep(dt)
903     if n%check_point==0:
904     self.writeXML()
905     self.doFinalization()
906 jgs 121
907 jgs 149 If one of the models in throws a C{FailedTimeStepError} exception a
908     new time step size is computed through getSafeTimeStepSize() and the
909     time step is repeated.
910 jgs 121
911 jgs 149 If one of the models in throws a C{IterationDivergenceError}
912     exception the time step size is halved and the time step is repeated.
913 jgs 121
914 jgs 149 In both cases the time integration is given up after
915     C{Simulation.FAILED_TIME_STEPS_MAX} attempts.
916 jgs 147 """
917     dt=self.UNDEF_DT
918     self.doInitialization()
919     while not self.finalize():
920     step_fail_counter=0
921     iteration_fail_counter=0
922 gross 836 if self.n==0:
923     dt_new=self.getSafeTimeStepSize(dt)
924     else:
925 gross 838 dt_new=min(max(self.getSafeTimeStepSize(dt),dt/self.MAX_CHANGE_OF_DT),dt*self.MAX_CHANGE_OF_DT)
926 jgs 147 self.trace("%d. time step %e (step size %e.)" % (self.n+1,self.tn+dt_new,dt_new))
927     end_of_step=False
928     while not end_of_step:
929     end_of_step=True
930     if not dt_new>0:
931 gross 829 raise NonPositiveStepSizeError("non-positive step size in step %d"%(self.n+1))
932 jgs 147 try:
933     self.doStepPreprocessing(dt_new)
934     self.doStep(dt_new)
935     self.doStepPostprocessing(dt_new)
936     except IterationDivergenceError:
937     dt_new*=0.5
938     end_of_step=False
939     iteration_fail_counter+=1
940     if iteration_fail_counter>self.FAILED_TIME_STEPS_MAX:
941 gross 836 raise SimulationBreakDownError("reduction of time step to achieve convergence failed after %s steps."%self.FAILED_TIME_STEPS_MAX)
942     self.trace("Iteration failed. Time step is repeated with new step size %s."%dt_new)
943 jgs 147 except FailedTimeStepError:
944     dt_new=self.getSafeTimeStepSize(dt)
945     end_of_step=False
946     step_fail_counter+=1
947 gross 836 self.trace("Time step is repeated with new time step size %s."%dt_new)
948 jgs 147 if step_fail_counter>self.FAILED_TIME_STEPS_MAX:
949 gross 836 raise SimulationBreakDownError("Time integration is given up after %d attempts."%step_fail_counter)
950 jgs 147 dt=dt_new
951     if not check_point==None:
952     if n%check_point==0:
953     self.trace("check point is created.")
954     self.writeXML()
955     self.doFinalization()
956 jgs 121
957 jgs 147 def fromDom(cls, doc):
958     sims = []
959     for node in doc.childNodes:
960     if isinstance(node, minidom.Text):
961     continue
962 jgs 121
963 jgs 147 sims.append(getComponent(node))
964 jgs 121
965 jgs 147 return cls(sims)
966 jgs 121
967 jgs 147 fromDom = classmethod(fromDom)
968 jgs 121
969    
970 jgs 147 class IterationDivergenceError(Exception):
971     """
972 jgs 149 Exception which is thrown if there is no convergence of the iteration
973     process at a time step.
974 jgs 121
975 jgs 149 But there is a chance that a smaller step could help to reach convergence.
976 jgs 147 """
977     pass
978 jgs 121
979 jgs 147 class FailedTimeStepError(Exception):
980 jgs 149 """
981     Exception which is thrown if the time step fails because of a step
982     size that have been choosen to be too large.
983     """
984 jgs 147 pass
985 jgs 121
986 jgs 147 class SimulationBreakDownError(Exception):
987 jgs 149 """
988     Exception which is thrown if the simulation does not manage to
989     progress in time.
990     """
991 jgs 147 pass
992 jgs 121
993 jgs 147 class NonPositiveStepSizeError(Exception):
994 jgs 149 """
995     Exception which is thrown if the step size is not positive.
996     """
997 jgs 147 pass
998 jgs 149
999     # 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