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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 637 - (show annotations)
Thu Mar 23 10:55:31 2006 UTC (17 years ago) by gross
File MIME type: text/x-python
File size: 29276 byte(s)
more copyright statements added
1 # $Id$
2
3 """
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 __copyright__=""" Copyright (c) 2006 by ACcESS MNRF
16 http://www.access.edu.au
17 Primary Business: Queensland, Australia"""
18 __license__="""Licensed under the Open Software License version 3.0
19 http://www.opensource.org/licenses/osl-3.0.php"""
20 __url__="http://www.iservo.edu.au/esys"
21 __version__="$Revision$"
22 __date__="$Date$"
23
24
25 from types import StringType,IntType,FloatType,BooleanType,ListType,DictType
26 from sys import stdout
27 import itertools
28 # import modellib temporarily removed!!!
29
30 # 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 from xml.dom import minidom
37
38 def dataNode(document, tagName, data):
39 """
40 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 """
47 t = document.createTextNode(str(data))
48 n = document.createElement(tagName)
49 n.appendChild(t)
50 return n
51
52 def esysDoc():
53 """
54 Global method for creating an instance of an EsysXML document.
55 """
56 doc = minidom.Document()
57 esys = doc.createElement('ESys')
58 doc.appendChild(esys)
59 return doc, esys
60
61 def all(seq):
62 for x in seq:
63 if not x:
64 return False
65 return True
66
67 def any(seq):
68 for x in seq:
69 if x:
70 return True
71 return False
72
73 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 """
85 Generic parse method for EsysXML. Without this, Links don't work.
86 """
87 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 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 def getComponent(doc):
114 """
115 Used to get components of Simualtions, Models.
116 """
117 for node in doc.childNodes:
118
119 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 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 if node.tagName == 'ParameterSet':
135 parameter_type = node.getAttribute("type")
136 return ParameterSet.fromDom(node)
137 raise "Invalid simulation type, %r" % node.getAttribute("type")
138
139
140 raise ValueError("No Simulation Found")
141
142
143 class Link:
144 """
145 A Link makes an attribute of an object callable::
146
147 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 """
155 Creates a link to the object target. If attribute is given, the link is
156 establised to this attribute of the target. Otherwise the attribute is
157 undefined.
158 """
159 self.target = target
160 self.attribute = None
161 self.setAttributeName(attribute)
162
163 def setAttributeName(self,attribute):
164 """
165 Set a new attribute name to be collected from the target object. The
166 target object must have the attribute with name attribute.
167 """
168 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 self.attribute = attribute
176
177 def hasDefinedAttributeName(self):
178 """
179 Returns true if an attribute name is set.
180 """
181 return self.attribute != None
182
183 def __repr__(self):
184 """
185 Returns a string representation of the link.
186 """
187 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 """
194 Returns the value of the attribute of the target object. If the
195 atrribute is callable then the return value of the call is returned.
196 """
197 if name:
198 out=getattr(self.target, name)
199 else:
200 out=getattr(self.target, self.attribute)
201
202 if callable(out):
203 return out()
204 else:
205 return out
206
207 def toDom(self, document, node):
208 """
209 C{toDom} method of Link. Creates a Link node and appends it to the
210 current XML document.
211 """
212 link = document.createElement('Link')
213 assert (self.target != None), ("Target was none, name was %r" % self.attribute)
214 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
221 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
228 fromDom = classmethod(fromDom)
229
230 def writeXML(self,ostream=stdout):
231 """
232 Writes an XML representation of self to the output stream ostream.
233 If ostream is nor present the standart output stream is used. If
234 esysheader==True the esys XML header is written
235 """
236 print 'I got to the Link writeXML method'
237 document, rootnode = esysDoc()
238 self.toDom(document, rootnode)
239
240 ostream.write(document.toprettyxml())
241
242 class LinkableObject(object):
243 """
244 An object that allows to link its attributes to attributes of other objects
245 via a Link object. For instance::
246
247 p = LinkableObject()
248 p.x = Link(o,"name")
249 print p.x
250
251 links attribute C{x} of C{p} to the attribute name of object C{o}.
252
253 C{p.x} will contain the current value of attribute C{name} of object
254 C{o}.
255
256 If the value of C{getattr(o, "name")} is callable, C{p.x} will return
257 the return value of the call.
258 """
259
260 number_sequence = itertools.count(100)
261
262 def __init__(self, debug=False):
263 """
264 Initializes LinkableObject so that we can operate on Links
265 """
266 self.debug = debug
267 self.__linked_attributes={}
268 self.id = self.number_sequence.next()
269 registerLinkableObject(self.id, self)
270
271 def trace(self, msg):
272 """
273 If debugging is on, print the message, otherwise do nothing
274 """
275 if self.debug:
276 print "%s: %s"%(str(self),msg)
277
278 def __getattr__(self,name):
279 """
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 out = self.getAttributeObject(name)
284 if isinstance(out,Link):
285 return out()
286 else:
287 return out
288
289 def getAttributeObject(self,name):
290 """
291 Return the object stored for attribute name.
292 """
293
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 if self.__class__.__dict__.has_key(name):
301 return self.__class.__dict__[name]
302
303 raise AttributeError,"No attribute %s."%name
304
305 def hasAttribute(self,name):
306 """
307 Returns True if self as attribute name.
308 """
309 return self.__dict__.has_key(name) or self.__linked_attributes.has_key(name) or self.__class__.__dict__.has_key(name)
310
311 def __setattr__(self,name,value):
312 """
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
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 self.trace("attribute %s is now linked by %s."%(name,value))
326 else:
327 self.__dict__[name] = value
328
329 def __delattr__(self,name):
330 """
331 Removes the attribute name.
332 """
333
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 class _ParameterIterator:
342 def __init__(self,parameterset):
343
344 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 class ParameterSet(LinkableObject):
355 """
356 A class which allows to emphazise attributes to be written and read to XML
357
358 Leaves of an ESySParameters object can be:
359
360 - 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
369 Example how to create an ESySParameters object::
370
371 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
375 This can be accessed as::
376
377 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 """
387 def __init__(self, parameters=[], **kwargs):
388 """
389 Creates a ParameterSet with parameters parameters.
390 """
391 LinkableObject.__init__(self, **kwargs)
392 self.parameters = set()
393 self.declareParameters(parameters)
394
395 def __repr__(self):
396 return "<%s %r>" % (self.__class__.__name__,
397 [(p, getattr(self, p, None)) for p in self.parameters])
398
399 def declareParameter(self,**parameters):
400 """
401 Declares a new parameter(s) and its (their) initial value.
402 """
403 self.declareParameters(parameters)
404
405 def declareParameters(self,parameters):
406 """
407 Declares a set of parameters. parameters can be a list, a dictionary
408 or a ParameterSet.
409 """
410 if isinstance(parameters,ListType):
411 parameters = zip(parameters, itertools.repeat(None))
412 if isinstance(parameters,DictType):
413 parameters = parameters.iteritems()
414
415 for prm, value in parameters:
416 setattr(self,prm,value)
417 self.parameters.add(prm)
418
419 self.trace("parameter %s has been declared."%prm)
420
421 def releaseParameters(self,name):
422 """
423 Removes parameter name from the paramameters.
424 """
425 if self.isParameter(name):
426 self.parameters.remove(name)
427 self.trace("parameter %s has been removed."%name)
428
429 def __iter__(self):
430 """
431 Creates an iterator over the parameter and their values.
432 """
433 return _ParameterIterator(self)
434
435 def showParameters(self):
436 """
437 Returns a descrition of the parameters.
438 """
439 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 """
452 Removes the attribute name.
453 """
454 LinkableObject.__delattr__(self,name)
455 try:
456 self.releaseParameter(name)
457 except:
458 pass
459
460 def toDom(self, document, node):
461 """
462 C{toDom} method of ParameterSet class.
463 """
464 pset = document.createElement('ParameterSet')
465 node.appendChild(pset)
466 self._parametersToDom(document, pset)
467
468 def _parametersToDom(self, document, node):
469 node.setAttribute('id', str(self.id))
470 node.setIdAttribute("id")
471 for name,value in self:
472 param = document.createElement('Parameter')
473 param.setAttribute('type', value.__class__.__name__)
474
475 param.appendChild(dataNode(document, 'Name', name))
476
477 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 Remove the empty nodes from the children of this node.
498 """
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
508 def _intfromValue(doc):
509 return int(doc.nodeValue.strip())
510
511 def _boolfromValue(doc):
512 return bool(doc.nodeValue.strip())
513
514 def _nonefromValue(doc):
515 return None
516
517 # 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 "int":_intfromValue,
524 "str":_stringfromValue,
525 "bool":_boolfromValue,
526 "NoneType":_nonefromValue
527 }
528
529 # print doc.toxml()
530
531 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 """
557 Writes the object as an XML object into an output stream.
558 """
559 # ParameterSet(d) with d[Name]=Value
560 document, node = esysDoc()
561 self.toDom(document, node)
562 ostream.write(document.toprettyxml())
563
564 class Model(ParameterSet):
565 """
566 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
571 doInitialization()
572 while not finalize():
573 dt=getSafeTimeStepSize(dt)
574 doStepPreprocessing(dt)
575 while not terminateIteration(): doStep(dt)
576 doStepPostprocessing(dt)
577 doFinalization()
578
579 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 """
585
586 UNDEF_DT=1.e300
587
588 def __init__(self,parameters=[],**kwarg):
589 """
590 Creates a model.
591
592 Just calls the parent constructor.
593 """
594 ParameterSet.__init__(self, parameters=parameters,**kwarg)
595
596 def __str__(self):
597 return "<%s %d>"%(self.__class__,id(self))
598
599 def toDom(self, document, node):
600 """
601 C{toDom} method of Model class
602 """
603 pset = document.createElement('Model')
604 pset.setAttribute('type', self.__class__.__name__)
605 if not self.__class__.__module__.startswith('esys.escript'):
606 pset.setAttribute('module', self.__class__.__module__)
607 node.appendChild(pset)
608 self._parametersToDom(document, pset)
609
610 def doInitialization(self):
611 """
612 Initializes the time stepping scheme.
613
614 This function may be overwritten.
615 """
616 pass
617
618 def getSafeTimeStepSize(self,dt):
619 """
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 return self.UNDEF_DT
627
628 def finalize(self):
629 """
630 Returns False if the time stepping is finalized.
631
632 This function may be overwritten.
633 """
634 return False
635
636 def doFinalization(self):
637 """
638 Finalizes the time stepping.
639
640 This function may be overwritten.
641 """
642 pass
643
644 def doStepPreprocessing(self,dt):
645 """
646 Sets up a time step of step size dt.
647
648 This function may be overwritten.
649 """
650 pass
651
652 def doStep(self,dt):
653 """
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 pass
661
662 def terminateIteration(self):
663 """
664 Returns True if iteration on a time step is terminated.
665 """
666 return True
667
668 def doStepPostprocessing(self,dt):
669 """
670 Finalalizes the time step.
671
672 dt is the currently used time step size.
673
674 This function may be overwritten.
675 """
676 pass
677
678 def writeXML(self, ostream=stdout):
679 document, node = esysDoc()
680 self.toDom(document, node)
681 ostream.write(document.toprettyxml())
682
683
684 class Simulation(Model):
685 """
686 A Simulation object is special Model which runs a sequence of Models.
687
688 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 """
693
694 FAILED_TIME_STEPS_MAX=20
695 MAX_ITER_STEPS=50
696
697 def __init__(self, models=[], **kwargs):
698 """
699 Initiates a simulation from a list of models.
700 """
701 Model.__init__(self, **kwargs)
702 self.__models=[]
703
704 for i in range(len(models)):
705 self[i] = models[i]
706
707
708 def __repr__(self):
709 """
710 Returns a string representation of the Simulation.
711 """
712 return "<Simulation %r>" % self.__models
713
714 def __str__(self):
715 """
716 Returning Simulation as a string.
717 """
718 return "<Simulation %d>"%id(self)
719
720 def iterModels(self):
721 """
722 Returns an iterator over the models.
723 """
724 return self.__models
725
726 def __getitem__(self,i):
727 """
728 Returns the i-th model.
729 """
730 return self.__models[i]
731
732 def __setitem__(self,i,value):
733 """
734 Sets the i-th model.
735 """
736 if not isinstance(value,Model):
737 raise ValueError("assigned value is not a Model, Model was:", i, " an instance of: ", value.__class__.__name__)
738 for j in range(max(i-len(self.__models)+1,0)):
739 self.__models.append(None)
740 self.__models[i]=value
741
742 def __len__(self):
743 """
744 Returns the number of models.
745 """
746 return len(self.__models)
747
748 def toDom(self, document, node):
749 """
750 C{toDom} method of Simulation class.
751 """
752 simulation = document.createElement('Simulation')
753 simulation.setAttribute('type', self.__class__.__name__)
754
755 for rank, sim in enumerate(self.iterModels()):
756 component = document.createElement('Component')
757 component.setAttribute('rank', str(rank))
758
759 sim.toDom(document, component)
760
761 simulation.appendChild(component)
762
763 node.appendChild(simulation)
764
765 def writeXML(self,ostream=stdout):
766 """
767 Writes the object as an XML object into an output stream.
768 """
769 document, rootnode = esysDoc()
770 self.toDom(document, rootnode)
771 targetsList = document.getElementsByTagName('Target')
772
773 for element in targetsList:
774 targetId = int(element.firstChild.nodeValue.strip())
775 if document.getElementById(str(targetId)):
776 continue
777 targetObj = LinkableObjectRegistry[targetId]
778 targetObj.toDom(document, rootnode)
779 ostream.write(document.toprettyxml())
780
781 def getSafeTimeStepSize(self,dt):
782 """
783 Returns a time step size which can safely be used by all models.
784
785 This is the minimum over the time step sizes of all models.
786 """
787 out=min([o.getSafeTimeStepSize(dt) for o in self.iterModels()])
788 #print "%s: safe step size is %e."%(str(self),out)
789 return out
790
791 def doInitialization(self):
792 """
793 Initializes all models.
794 """
795 self.n=0
796 self.tn=0.
797 for o in self.iterModels():
798 o.doInitialization()
799
800 def finalize(self):
801 """
802 Returns True if any of the models is to be finalized.
803 """
804 return any([o.finalize() for o in self.iterModels()])
805
806 def doFinalization(self):
807 """
808 Finalalizes the time stepping for all models.
809 """
810 for i in self.iterModels(): i.doFinalization()
811 self.trace("end of time integation.")
812
813 def doStepPreprocessing(self,dt):
814 """
815 Initializes the time step for all models.
816 """
817 for o in self.iterModels():
818 o.doStepPreprocessing(dt)
819
820 def terminateIteration(self):
821 """
822 Returns True if all iterations for all models are terminated.
823 """
824 out=all([o.terminateIteration() for o in self.iterModels()])
825 return out
826
827 def doStepPostprocessing(self,dt):
828 """
829 Finalalizes the iteration process for all models.
830 """
831 for o in self.iterModels():
832 o.doStepPostprocessing(dt)
833 self.n+=1
834 self.tn+=dt
835
836 def doStep(self,dt):
837 """
838 Executes the iteration step at a time step for all model::
839
840 self.doStepPreprocessing(dt)
841 while not self.terminateIteration():
842 for all models:
843 self.doStep(dt)
844 self.doStepPostprocessing(dt)
845 """
846 self.iter=0
847 while not self.terminateIteration():
848 if self.iter==0: self.trace("iteration at %d-th time step %e starts"%(self.n+1,self.tn+dt))
849 self.iter+=1
850 self.trace("iteration step %d"%(self.iter))
851 for o in self.iterModels():
852 o.doStep(dt)
853 if self.iter>0: self.trace("iteration at %d-th time step %e finalized."%(self.n+1,self.tn+dt))
854
855 def run(self,check_point=None):
856 """
857 Run the simulation by performing essentially::
858
859 self.doInitialization()
860 while not self.finalize():
861 dt=self.getSafeTimeStepSize()
862 self.doStep(dt)
863 if n%check_point==0:
864 self.writeXML()
865 self.doFinalization()
866
867 If one of the models in throws a C{FailedTimeStepError} exception a
868 new time step size is computed through getSafeTimeStepSize() and the
869 time step is repeated.
870
871 If one of the models in throws a C{IterationDivergenceError}
872 exception the time step size is halved and the time step is repeated.
873
874 In both cases the time integration is given up after
875 C{Simulation.FAILED_TIME_STEPS_MAX} attempts.
876 """
877 dt=self.UNDEF_DT
878 self.doInitialization()
879 while not self.finalize():
880 step_fail_counter=0
881 iteration_fail_counter=0
882 dt_new=self.getSafeTimeStepSize(dt)
883 self.trace("%d. time step %e (step size %e.)" % (self.n+1,self.tn+dt_new,dt_new))
884 end_of_step=False
885 while not end_of_step:
886 end_of_step=True
887 if not dt_new>0:
888 raise NonPositiveStepSizeError("non-positive step size in step %d",self.n+1)
889 try:
890 self.doStepPreprocessing(dt_new)
891 self.doStep(dt_new)
892 self.doStepPostprocessing(dt_new)
893 except IterationDivergenceError:
894 dt_new*=0.5
895 end_of_step=False
896 iteration_fail_counter+=1
897 if iteration_fail_counter>self.FAILED_TIME_STEPS_MAX:
898 raise SimulationBreakDownError("reduction of time step to achieve convergence failed.")
899 self.trace("iteration fails. time step is repeated with new step size.")
900 except FailedTimeStepError:
901 dt_new=self.getSafeTimeStepSize(dt)
902 end_of_step=False
903 step_fail_counter+=1
904 self.trace("time step is repeated.")
905 if step_fail_counter>self.FAILED_TIME_STEPS_MAX:
906 raise SimulationBreakDownError("time integration is given up after %d attempts."%step_fail_counter)
907 dt=dt_new
908 if not check_point==None:
909 if n%check_point==0:
910 self.trace("check point is created.")
911 self.writeXML()
912 self.doFinalization()
913
914 def fromDom(cls, doc):
915 sims = []
916 for node in doc.childNodes:
917 if isinstance(node, minidom.Text):
918 continue
919
920 sims.append(getComponent(node))
921
922 return cls(sims)
923
924 fromDom = classmethod(fromDom)
925
926
927 class IterationDivergenceError(Exception):
928 """
929 Exception which is thrown if there is no convergence of the iteration
930 process at a time step.
931
932 But there is a chance that a smaller step could help to reach convergence.
933 """
934 pass
935
936 class FailedTimeStepError(Exception):
937 """
938 Exception which is thrown if the time step fails because of a step
939 size that have been choosen to be too large.
940 """
941 pass
942
943 class SimulationBreakDownError(Exception):
944 """
945 Exception which is thrown if the simulation does not manage to
946 progress in time.
947 """
948 pass
949
950 class NonPositiveStepSizeError(Exception):
951 """
952 Exception which is thrown if the step size is not positive.
953 """
954 pass
955
956 # 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