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