/[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 3892 - (show annotations)
Tue Apr 10 08:57:23 2012 UTC (6 years, 11 months ago) by jfenwick
File MIME type: text/x-python
File size: 47784 byte(s)
Merged changes across from the attempt2 branch.
This version builds and passes python2 tests.
It also passes most python3 tests.



1 # -*- coding: utf-8 -*-
2
3 ########################################################
4 #
5 # Copyright (c) 2003-2010 by University of Queensland
6 # Earth Systems Science Computational Center (ESSCC)
7 # http://www.uq.edu.au/esscc
8 #
9 # Primary Business: Queensland, Australia
10 # Licensed under the Open Software License version 3.0
11 # http://www.opensource.org/licenses/osl-3.0.php
12 #
13 ########################################################
14
15 __copyright__="""Copyright (c) 2003-2010 by University of Queensland
16 Earth Systems Science Computational Center (ESSCC)
17 http://www.uq.edu.au/esscc
18 Primary Business: Queensland, Australia"""
19 __license__="""Licensed under the Open Software License version 3.0
20 http://www.opensource.org/licenses/osl-3.0.php"""
21 __url__="https://launchpad.net/escript-finley"
22
23 """
24 Environment for implementing models in escript
25
26 :var __author__: name of author
27 :var __copyright__: copyrights
28 :var __license__: licence agreement
29 :var __url__: url entry point on documentation
30 :var __version__: version
31 :var __date__: date of the version
32 """
33 __author__="Lutz Gross, l.gross@uq.edu.au"
34
35
36 import sys
37 import numpy
38 import operator
39 import itertools
40 import time
41 import os
42 import collections
43 from functools import reduce
44
45
46 # import the 'set' module if it's not defined (python2.3/2.4 difference)
47 try:
48 set
49 except NameError:
50 from sets import Set as set
51
52 from xml.dom import minidom
53
54
55 def all(seq):
56 """
57 Returns True if no element in ``seq`` is ``None``, False otherwise.
58 """
59 for x in seq:
60 if not x:
61 return False
62 return True
63
64 def any(seq):
65 """
66 Returns True if not all elements in ``seq`` are ``None``, False otherwise.
67 """
68 for x in seq:
69 if x:
70 return True
71 return False
72
73 def importName(modulename, name):
74 """
75 Imports a named object from a module in the context of this function,
76 which means you should use fully qualified module paths.
77 Returns None on failure.
78
79 This function is from U{http://aspn.activestate.com/ASPN/Cookbook/Python/Recipe/52241}
80 """
81 module = __import__(modulename, globals(), locals(), [name])
82
83 try:
84 return vars(module)[name]
85 except KeyError:
86 raise ImportError("Could not import %s from %s" % (name, modulename))
87
88 class ESySXMLParser(object):
89 """
90 Parser for an ESysXML file.
91 """
92 def __init__(self,xml, debug=False):
93 if sys.version_info[0]<3:
94 xml=str(xml) # xml might be unicode
95 #print("\n")
96 #print(type(xml))
97 #print(xml)
98 #print("\n")
99 self.__dom = minidom.parseString(xml)
100 self.__linkable_object_registry= {}
101 self.__link_registry= []
102 self.__esys=self.__dom.getElementsByTagName('ESys')[0]
103 self.debug=debug
104
105 def getClassPath(self, node):
106 type = node.getAttribute("type")
107 if (node.getAttribute("module")):
108 module = node.getAttribute("module")
109 return importName(module, type)
110 else:
111 return importName("__main__", type)
112
113 def setLinks(self):
114 for obj_id, link in self.__link_registry:
115 link.target = self.__linkable_object_registry[obj_id]
116
117 def parse(self):
118 """
119 Parses EsysXML and returns the list of generating ParameterSets.
120 """
121 found=[]
122 for node in self.__esys.childNodes:
123 if isinstance(node, minidom.Element):
124 if node.tagName == 'Simulation':
125 found.append(Simulation.fromDom(self, node))
126 elif node.tagName == 'Model':
127 found.append(self.getClassPath(node).fromDom(self, node))
128 elif node.tagName == 'ParameterSet':
129 found.append(self.getClassPath(node).fromDom(self, node))
130 else:
131 raise "Invalid type, %r" % node.getAttribute("type")
132 self.setLinks()
133 return found
134
135 def registerLink(self,obj_id, link):
136 self.__link_registry.append((int(obj_id),link))
137
138 def registerLinkableObject(self,obj, node):
139 id_str=node.getAttribute('id').strip()
140 if len(id_str)>0:
141 id=int(id_str)
142 if id in self.__linkable_object_registry:
143 raise ValueError("Object id %s already exists."%id)
144 else:
145 self.__linkable_object_registry[id]=obj
146
147 def getComponent(self, node):
148 """
149 Returns a single component + rank from a simulation.
150 """
151 rank = int(node.getAttribute("rank"))
152 for n in node.childNodes:
153 if isinstance(n, minidom.Element):
154 if n.tagName == 'Simulation':
155 return (rank, Simulation.fromDom(self, n))
156 elif n.tagName == 'Model':
157 return (rank, self.getClassPath(n).fromDom(self, n))
158 elif n.tagName == 'ParameterSet':
159 return (rank, self.getClassPath(n).fromDom(self, n))
160 else:
161 raise ValueError("illegal component type %s"%n.tagName)
162 raise ValueError("cannot resolve Component")
163
164 class ESySXMLCreator(object):
165 """
166 Creates an XML Dom representation.
167 """
168 def __init__(self):
169 self.__dom=minidom.Document()
170 self.__esys =self.__dom.createElement('ESys')
171 self.__dom.appendChild(self.__esys)
172 self.__linkable_object_registry={}
173 self.__number_sequence = itertools.count(100)
174
175 def getRoot(self):
176 return self.__esys
177
178 def createElement(self,name):
179 return self.__dom.createElement(name)
180
181 def createTextNode(self,name):
182 return self.__dom.createTextNode(name)
183
184 def getElementById(self,name):
185 return self.__dom.getElementById(name)
186
187 def createDataNode(self, tagName, data):
188 """
189 ``createDataNode`` s are the building blocks of the XML documents
190 constructed in this module.
191
192 :param tagName: the associated XML tag
193 :param data: the values in the tag
194 """
195 n = self.createElement(tagName)
196 n.appendChild(self.createTextNode(str(data)))
197 return n
198
199 def getLinkableObjectId(self, obj):
200 for id, o in list(self.__linkable_object_registry.items()):
201 if o == obj: return id
202 id =next(self.__number_sequence)
203 self.__linkable_object_registry[id]=obj
204 return id
205
206 def registerLinkableObject(self, obj, node):
207 """
208 Returns a unique object id for object ``obj``.
209 """
210 id=self.getLinkableObjectId(obj)
211 node.setAttribute('id',str(id))
212 node.setIdAttribute("id")
213
214 def includeTargets(self):
215 target_written=True
216 while target_written:
217 targetsList =self.__dom.getElementsByTagName('Target')
218 target_written=False
219 for element in targetsList:
220 targetId = int(element.firstChild.nodeValue.strip())
221 if self.getElementById(str(targetId)): continue
222 targetObj = self.__linkable_object_registry[targetId]
223 targetObj.toDom(self, self.__esys)
224 target_written=True
225
226 def toprettyxml(self):
227 self.includeTargets()
228 return self.__dom.toprettyxml()
229
230 class Link:
231 """
232 A Link makes an attribute of an object callable::
233
234 o.object()
235 o.a=8
236 l=Link(o,"a")
237 assert l()==8
238 """
239
240 def __init__(self,target,attribute=None):
241 """
242 Creates a link to the object target. If attribute is given, the link is
243 established to this attribute of the target. Otherwise the attribute is
244 undefined.
245 """
246 self.target = target
247 self.attribute = None
248 self.setAttributeName(attribute)
249
250 def getTarget(self):
251 """
252 Returns the target.
253 """
254 return self.target
255
256 def getAttributeName(self):
257 """
258 Returns the name of the attribute the link is pointing to.
259 """
260 return self.attribute
261
262 def setAttributeName(self,attribute):
263 """
264 Sets a new attribute name to be collected from the target object. The
265 target object must have the attribute with name attribute.
266 """
267 if attribute and self.target:
268 if isinstance(self.target,LinkableObject):
269 if not self.target.hasAttribute(attribute):
270 raise AttributeError("%s: target %s has no attribute %s."%(self, self.target, attribute))
271 else:
272 if not hasattr(self.target,attribute):
273 raise AttributeError("%s: target %s has no attribute %s."%(self, self.target, attribute))
274 self.attribute = attribute
275
276 def hasDefinedAttributeName(self):
277 """
278 Returns true if an attribute name is set.
279 """
280 return self.attribute != None
281
282 def __repr__(self):
283 """
284 Returns a string representation of the link.
285 """
286 if self.hasDefinedAttributeName():
287 return "<Link to attribute %s of %s>" % (self.attribute, self.target)
288 else:
289 return "<Link to target %s>" % self.target
290
291 def __call__(self,name=None):
292 """
293 Returns the value of the attribute of the target object. If the
294 attribute is callable then the return value of the call is returned.
295 """
296 if name:
297 out=getattr(self.target, name)
298 else:
299 out=getattr(self.target, self.attribute)
300
301 if isinstance(out, collections.Callable):
302 return out()
303 else:
304 return out
305
306 def toDom(self, esysxml, node):
307 """
308 ``toDom`` method of Link. Creates a Link node and appends it to the
309 current XML esysxml.
310 """
311 link = esysxml.createElement('Link')
312 assert (self.target != None), ("Target was none, name was %r" % self.attribute)
313 link.appendChild(esysxml.createDataNode('Target', esysxml.getLinkableObjectId(self.target)))
314 # this use of id will not work for purposes of being able to retrieve the intended
315 # target from the xml later. I need a better unique identifier.
316 assert self.attribute, "You can't xmlify a Link without a target attribute"
317 link.appendChild(esysxml.createDataNode('Attribute', self.attribute))
318 node.appendChild(link)
319
320 def fromDom(cls, esysxml, node):
321 targetid = int(node.getElementsByTagName("Target")[0].firstChild.nodeValue.strip())
322 attribute =str(node.getElementsByTagName("Attribute")[0].firstChild.nodeValue.strip())
323 l = cls(None, attribute)
324 esysxml.registerLink(targetid, l)
325 return l
326
327 fromDom = classmethod(fromDom)
328
329 class LinkableObject(object):
330 """
331 An object that allows to link its attributes to attributes of other objects
332 via a Link object. For instance::
333
334 p = LinkableObject()
335 p.x = Link(o,"name")
336 print p.x
337
338 links attribute ``x`` of ``p`` to the attribute name of object ``o``.
339
340 ``p.x`` will contain the current value of attribute ``name`` of object ``o``.
341
342 If the value of ``getattr(o, "name")`` is callable, ``p.x`` will return
343 the return value of the call.
344 """
345
346 def __init__(self, id = None, debug=False):
347 """
348 Initializes LinkableObject so that we can operate on Links.
349 """
350 self.debug = debug
351 self.__linked_attributes={}
352
353 def trace(self, msg):
354 """
355 If debugging is on, prints the message, otherwise does nothing.
356 """
357 if self.debug:
358 print("%s: %s"%(str(self),msg))
359
360 def __getattr__(self,name):
361 """
362 Returns the value of attribute name. If the value is a Link object the
363 object is called and the return value is returned.
364 """
365 out = self.getAttributeObject(name)
366 if isinstance(out,Link):
367 return out()
368 else:
369 return out
370
371 def getAttributeObject(self,name):
372 """
373 Returns the object stored for attribute ``name``.
374 """
375
376 if name in self.__dict__:
377 return self.__dict__[name]
378
379 if name in self.__linked_attributes:
380 return self.__linked_attributes[name]
381
382 if name in self.__class__.__dict__:
383 return self.__class.__dict__[name]
384
385 raise AttributeError("No attribute %s."%name)
386
387 def hasAttribute(self,name):
388 """
389 Returns True if self has attribute ``name``.
390 """
391 return name in self.__dict__ or name in self.__linked_attributes or name in self.__class__.__dict__
392
393 def __setattr__(self,name,value):
394 """
395 Sets the value for attribute name. If value is a Link the target
396 attribute is set to name if no attribute has been specified.
397 """
398
399 if name in self.__dict__:
400 del self.__dict__[name]
401
402 if isinstance(value,Link):
403 if not value.hasDefinedAttributeName():
404 value.setAttributeName(name)
405 self.__linked_attributes[name] = value
406
407 self.trace("attribute %s is now linked by %s."%(name,value))
408 else:
409 self.__dict__[name] = value
410
411 def __delattr__(self,name):
412 """
413 Removes the attribute ``name``.
414 """
415
416 if self.__linked_attributes.has_key[name]:
417 del self.__linked_attributes[name]
418 elif name in self.__dict__:
419 del self.__dict__[name]
420 else:
421 raise AttributeError("No attribute %s."%name)
422
423 class _ParameterIterator:
424 def __init__(self,parameterset):
425
426 self.__set=parameterset
427 self.__iter=iter(parameterset.parameters)
428
429 def __next__(self):
430 o=next(self.__iter)
431 return (o,self.__set.getAttributeObject(o))
432
433 def next(self): #Still needed by py2.6
434 return self.__next__()
435
436 def __iter__(self):
437 return self
438
439 class ParameterSet(LinkableObject):
440 """
441 A class which allows to emphasize attributes to be written and read to XML.
442
443 Leaves of an ESySParameters object can be:
444
445 - a real number
446 - an integer number
447 - a string
448 - a boolean value
449 - a ParameterSet object
450 - a Simulation object
451 - a Model object
452 - a numpy object
453 - a list of booleans
454 - any other object (not considered by writeESySXML and writeXML)
455
456 Example for how to create an ESySParameters object::
457
458 p11=ParameterSet(gamma1=1.,gamma2=2.,gamma3=3.)
459 p1=ParameterSet(dim=2,tol_v=0.001,output_file="/tmp/u.%3.3d.dx",runFlag=True,parm11=p11)
460 parm=ParameterSet(parm1=p1,parm2=ParameterSet(alpha=Link(p11,"gamma1")))
461
462 This can be accessed as::
463
464 parm.parm1.gamma=0.
465 parm.parm1.dim=2
466 parm.parm1.tol_v=0.001
467 parm.parm1.output_file="/tmp/u.%3.3d.dx"
468 parm.parm1.runFlag=True
469 parm.parm1.parm11.gamma1=1.
470 parm.parm1.parm11.gamma2=2.
471 parm.parm1.parm11.gamma3=3.
472 parm.parm2.alpha=1. (value of parm.parm1.parm11.gamma1)
473 """
474 def __init__(self, parameters=[], **kwargs):
475 """
476 Creates a ParameterSet with given parameters.
477 """
478 LinkableObject.__init__(self, **kwargs)
479 self.parameters = set()
480 self.declareParameters(parameters)
481
482 def __repr__(self):
483 return "<%s %d>"%(self.__class__.__name__,id(self))
484
485 def declareParameter(self,**parameters):
486 """
487 Declares one or more new parameters and their initial value.
488 """
489 self.declareParameters(parameters)
490
491 def declareParameters(self,parameters):
492 """
493 Declares a set of parameters. parameters can be a list, a dictionary
494 or a ParameterSet.
495 """
496 if isinstance(parameters,type([])):
497 parameters = list(zip(parameters, itertools.repeat(None)))
498 if isinstance(parameters,type(dict())):
499 parameters = iter(list(parameters.items()))
500
501 for prm, value in parameters:
502 setattr(self,prm,value)
503 self.parameters.add(prm)
504
505 def releaseParameters(self,name):
506 """
507 Removes parameter name from the parameters.
508 """
509 if self.isParameter(name):
510 self.parameters.remove(name)
511 self.trace("parameter %s has been removed."%name)
512
513 def checkLinkTargets(self, models, hash):
514 """
515 Returns a set of tuples ("<self>(<name>)", <target model>) if the
516 parameter <name> is linked to model <target model> but <target model>
517 is not in the list of models. If a parameter is linked to another
518 parameter set which is not in the hash list the parameter set is
519 checked for its models. hash gives the call history.
520 """
521 out=set()
522 for name, value in self:
523 if isinstance(value, Link):
524 m=value.getTarget()
525 if isinstance(m, Model):
526 if not m in models: out.add( (str(self)+"("+name+")",m) )
527 elif isinstance(m, ParameterSet) and not m in hash:
528 out|=set( [ (str(self)+"("+name+")."+f[0],f[1]) for f in m.checkLinkTargets(models, hash+[ self ] ) ] )
529 return out
530
531 def __iter__(self):
532 """
533 Creates an iterator over the parameter and their values.
534 """
535 return _ParameterIterator(self)
536
537 def showParameters(self):
538 """
539 Returns a description of the parameters.
540 """
541 out="{"
542 notfirst=False
543 for i,v in self:
544 if notfirst: out=out+","
545 notfirst=True
546 if isinstance(v,ParameterSet):
547 out="%s\"%s\" : %s"%(out,i,v.showParameters())
548 else:
549 out="%s\"%s\" : %s"%(out,i,v)
550 return out+"}"
551
552 def __delattr__(self,name):
553 """
554 Removes the attribute ``name``.
555 """
556 LinkableObject.__delattr__(self,name)
557 try:
558 self.releaseParameter(name)
559 except:
560 pass
561
562 def toDom(self, esysxml, node):
563 """
564 ``toDom`` method of Model class.
565 """
566 pset = esysxml.createElement('ParameterSet')
567 pset.setAttribute('type', self.__class__.__name__)
568 pset.setAttribute('module', self.__class__.__module__)
569 esysxml.registerLinkableObject(self, pset)
570 self._parametersToDom(esysxml, pset)
571 node.appendChild(pset)
572
573 def _parametersToDom(self, esysxml, node):
574 for name,value in self:
575 # convert list to numpy when possible:
576 if isinstance (value, list):
577 elem_type=-1
578 for i in value:
579 if isinstance(i, bool):
580 elem_type = max(elem_type,0)
581 elif isinstance(i, int):
582 elem_type = max(elem_type,1)
583 elif isinstance(i, float):
584 elem_type = max(elem_type,2)
585 if elem_type == 0: value = numpy.array(value,numpy.bool_)
586 if elem_type == 1: value = numpy.array(value,numpy.int_)
587 if elem_type == 2: value = numpy.array(value,numpy.float_)
588
589 param = esysxml.createElement('Parameter')
590 param.setAttribute('type', value.__class__.__name__)
591
592 param.appendChild(esysxml.createDataNode('Name', name))
593
594 val = esysxml.createElement('Value')
595 if isinstance(value,(ParameterSet,Link,DataSource)):
596 value.toDom(esysxml, val)
597 param.appendChild(val)
598 elif isinstance(value, numpy.ndarray):
599 shape = value.shape
600 if isinstance(shape, tuple):
601 size = reduce(operator.mul, shape)
602 shape = ' '.join(map(str, shape))
603 else:
604 size = shape
605 shape = str(shape)
606
607 arraytype = value.dtype.kind
608 if arraytype=='b':
609 arraytype_str="bool_"
610 elif arraytype=='i':
611 arraytype_str="int_"
612 elif arraytype=='f':
613 arraytype_str="float_"
614 elif arraytype=='c':
615 arraytype_str="complex_"
616 else:
617 arraytype_str=str(arraytype)
618 ndarrayElement = esysxml.createElement('ndarray')
619 ndarrayElement.appendChild(esysxml.createDataNode('ArrayType', arraytype_str))
620 ndarrayElement.appendChild(esysxml.createDataNode('Shape', shape))
621 ndarrayElement.appendChild(esysxml.createDataNode('Data', ' '.join(
622 [str(x) for x in numpy.reshape(value, size)])))
623 val.appendChild(ndarrayElement)
624 param.appendChild(val)
625 elif isinstance(value, list):
626 param.appendChild(esysxml.createDataNode('Value', ' '.join([str(x) for x in value]) ))
627 elif isinstance(value, (str, bool, int, float, type(None))):
628 param.appendChild(esysxml.createDataNode('Value', str(value)))
629 elif isinstance(value, dict):
630 dic = esysxml.createElement('dictionary')
631 if len(list(value.keys()))>0:
632 dic.setAttribute('key_type', list(value.keys())[0].__class__.__name__)
633 dic.setAttribute('value_type', value[list(value.keys())[0]].__class__.__name__)
634 for k,v in list(value.items()):
635 i=esysxml.createElement('item')
636 i.appendChild(esysxml.createDataNode('key', k))
637 i.appendChild(esysxml.createDataNode('value', v))
638 dic.appendChild(i)
639 param.appendChild(dic)
640 else:
641 raise ValueError("cannot serialize %s type to XML."%str(value.__class__))
642
643 node.appendChild(param)
644
645 def fromDom(cls, esysxml, node):
646 # Define a host of helper functions to assist us.
647 def _children(node):
648 """
649 Remove the empty nodes from the children of this node.
650 """
651 ret = []
652 for x in node.childNodes:
653 if isinstance(x, minidom.Text):
654 if x.nodeValue.strip():
655 ret.append(x)
656 else:
657 ret.append(x)
658 return ret
659
660 def _floatfromValue(esysxml, node):
661 return float(node.nodeValue.strip())
662
663 def _stringfromValue(esysxml, node):
664 return str(node.nodeValue.strip())
665
666 def _intfromValue(esysxml, node):
667 return int(node.nodeValue.strip())
668
669 def _boolfromValue(esysxml, node):
670 return _boolfromstring(node.nodeValue.strip())
671
672 def _nonefromValue(esysxml, node):
673 return None
674
675 def _ndarrayfromValue(esysxml, node):
676 for node in _children(node):
677 if node.tagName == 'ArrayType':
678 arraytype = node.firstChild.nodeValue.strip()
679 if node.tagName == 'Shape':
680 shape = node.firstChild.nodeValue.strip()
681 shape = [int(x) for x in shape.split()]
682 if node.tagName == 'Data':
683 data = node.firstChild.nodeValue.strip()
684 ndtype=getattr(numpy,arraytype)
685 if ndtype==numpy.bool_:
686 data=[(x=="True") for x in data.split()]
687 else:
688 data=[ndtype(x) for x in data.split()]
689 return numpy.reshape(numpy.array(data, dtype=ndtype),
690 shape)
691
692 def _listfromValue(esysxml, node):
693 return [x for x in node.nodeValue.split()]
694
695 def _boolfromstring(s):
696 if s == 'True':
697 return True
698 else:
699 return False
700 # Mapping from text types in the xml to methods used to process trees of that type
701 ptypemap = {"Simulation": Simulation.fromDom,
702 "Model":Model.fromDom,
703 "ParameterSet":ParameterSet.fromDom,
704 "Link":Link.fromDom,
705 "DataSource":DataSource.fromDom,
706 "float":_floatfromValue,
707 "int":_intfromValue,
708 "str":_stringfromValue,
709 "bool":_boolfromValue,
710 "list":_listfromValue,
711 "ndarray" : _ndarrayfromValue,
712 "NoneType":_nonefromValue,
713 }
714
715 parameters = {}
716 for n in _children(node):
717 ptype = n.getAttribute("type")
718 if ptype not in ptypemap:
719 raise KeyError("cannot handle parameter type %s."%ptype)
720
721 pname = pvalue = None
722 for childnode in _children(n):
723 if childnode.tagName == "Name":
724 pname = childnode.firstChild.nodeValue.strip()
725
726 if childnode.tagName == "Value":
727 nodes = _children(childnode)
728 pvalue = ptypemap[ptype](esysxml, nodes[0])
729
730 parameters[pname] = pvalue
731
732 # Create the instance of ParameterSet
733 try:
734 o = cls(debug=esysxml.debug)
735 except TypeError as inst:
736 print(inst.args[0])
737 if inst.args[0]=="__init__() got an unexpected keyword argument 'debug'":
738 raise TypeError("The Model class %s __init__ needs to have argument 'debug'.")
739 else:
740 raise inst
741 o.declareParameters(parameters)
742 esysxml.registerLinkableObject(o, node)
743 return o
744
745 fromDom = classmethod(fromDom)
746
747 def writeXML(self,ostream=sys.stdout):
748 """
749 Writes the object as an XML object into an output stream.
750 """
751 esysxml=ESySXMLCreator()
752 self.toDom(esysxml, esysxml.getRoot())
753 if sys.version_info[0]<3:
754 ostream.write(unicode(esysxml.toprettyxml()))
755 else:
756 ostream.write(esysxml.toprettyxml())
757
758 class Model(ParameterSet):
759 """
760 A Model object represents a process marching over time until a
761 finalizing condition is fulfilled. At each time step an iterative
762 process can be performed and the time step size can be controlled. A
763 Model has the following work flow::
764
765 doInitialization()
766 while not terminateInitialIteration(): doInitialStep()
767 doInitialPostprocessing()
768 while not finalize():
769 dt=getSafeTimeStepSize(dt)
770 doStepPreprocessing(dt)
771 while not terminateIteration(): doStep(dt)
772 doStepPostprocessing(dt)
773 doFinalization()
774
775 where ``doInitialization``, ``finalize``, ``getSafeTimeStepSize``,
776 ``doStepPreprocessing``, ``terminateIteration``, ``doStepPostprocessing``,
777 ``doFinalization`` are methods of the particular instance of a Model. The
778 default implementations of these methods have to be overwritten by the
779 subclass implementing a Model.
780 """
781
782 UNDEF_DT=1.e300
783
784 def __init__(self,parameters=[],**kwargs):
785 """
786 Creates a model.
787
788 Just calls the parent constructor.
789 """
790 ParameterSet.__init__(self, parameters=parameters,**kwargs)
791
792 def __str__(self):
793 return "<%s %d>"%(self.__class__.__name__,id(self))
794
795
796 def setUp(self):
797 """
798 Sets up the model.
799
800 This function may be overwritten.
801 """
802 pass
803
804 def doInitialization(self):
805 """
806 Initializes the time stepping scheme. This method is not called in
807 case of a restart.
808
809 This function may be overwritten.
810 """
811 pass
812
813 def doInitialStep(self):
814 """
815 Performs an iteration step in the initialization phase. This method
816 is not called in case of a restart.
817
818 This function may be overwritten.
819 """
820 pass
821
822 def terminateInitialIteration(self):
823 """
824 Returns True if iteration at the inital phase is terminated.
825 """
826 return True
827
828 def doInitialPostprocessing(self):
829 """
830 Finalises the initialization iteration process. This method is not
831 called in case of a restart.
832
833 This function may be overwritten.
834 """
835 pass
836
837 def getSafeTimeStepSize(self,dt):
838 """
839 Returns a time step size which can be safely used.
840
841 ``dt`` gives the previously used step size.
842
843 This function may be overwritten.
844 """
845 return self.UNDEF_DT
846
847 def finalize(self):
848 """
849 Returns False if the time stepping is finalized.
850
851 This function may be overwritten.
852 """
853 return False
854
855 def doFinalization(self):
856 """
857 Finalizes the time stepping.
858
859 This function may be overwritten.
860 """
861 pass
862
863 def doStepPreprocessing(self,dt):
864 """
865 Sets up a time step of step size dt.
866
867 This function may be overwritten.
868 """
869 pass
870
871 def doStep(self,dt):
872 """
873 Executes an iteration step at a time step.
874
875 ``dt`` is the currently used time step size.
876
877 This function may be overwritten.
878 """
879 pass
880
881 def terminateIteration(self):
882 """
883 Returns True if iteration on a time step is terminated.
884 """
885 return True
886
887 def doStepPostprocessing(self,dt):
888 """
889 Finalises the time step.
890
891 dt is the currently used time step size.
892
893 This function may be overwritten.
894 """
895 pass
896
897 def toDom(self, esysxml, node):
898 """
899 ``toDom`` method of Model class.
900 """
901 pset = esysxml.createElement('Model')
902 pset.setAttribute('type', self.__class__.__name__)
903 pset.setAttribute('module', self.__class__.__module__)
904 esysxml.registerLinkableObject(self, pset)
905 node.appendChild(pset)
906 self._parametersToDom(esysxml, pset)
907
908 class Simulation(Model):
909 """
910 A Simulation object is a special Model which runs a sequence of Models.
911
912 The methods ``doInitialization``, ``finalize``, ``getSafeTimeStepSize``,
913 ``doStepPreprocessing``, ``terminateIteration``, ``doStepPostprocessing``,
914 ``doFinalization`` execute the corresponding methods of the models in
915 the simulation.
916 """
917
918 FAILED_TIME_STEPS_MAX=20
919 MAX_ITER_STEPS=50
920 MAX_CHANGE_OF_DT=2.
921
922 def __init__(self, models=[], **kwargs):
923 """
924 Initiates a simulation from a list of models.
925 """
926 super(Simulation, self).__init__(**kwargs)
927 self.declareParameter(time=0.,
928 time_step=0,
929 dt = self.UNDEF_DT)
930 for m in models:
931 if not isinstance(m, Model):
932 raise TypeError("%s is not a subclass of Model."%m)
933 self.__models=[]
934 for i in range(len(models)):
935 self[i] = models[i]
936
937
938 def __repr__(self):
939 """
940 Returns a string representation of the Simulation.
941 """
942 return "<Simulation %r>" % self.__models
943
944 def __str__(self):
945 """
946 Returns Simulation as a string.
947 """
948 return "<Simulation %d>" % id(self)
949
950 def iterModels(self):
951 """
952 Returns an iterator over the models.
953 """
954 return self.__models
955
956 def __getitem__(self,i):
957 """
958 Returns the i-th model.
959 """
960 return self.__models[i]
961
962 def __setitem__(self,i,value):
963 """
964 Sets the i-th model.
965 """
966 if not isinstance(value,Model):
967 raise ValueError("assigned value is not a Model but instance of %s"%(value.__class__.__name__,))
968 for j in range(max(i-len(self.__models)+1,0)):
969 self.__models.append(None)
970 self.__models[i]=value
971
972 def __len__(self):
973 """
974 Returns the number of models.
975 """
976 return len(self.__models)
977
978 def getAllModels(self):
979 """
980 Returns a list of all models used in the Simulation including
981 subsimulations.
982 """
983 out=[]
984 for m in self.iterModels():
985 if isinstance(m, Simulation):
986 out+=m.getAllModels()
987 else:
988 out.append(m)
989 return list(set(out))
990
991 def checkModels(self, models, hash):
992 """
993 Returns a list of (model, parameter, target model) if the parameter
994 of model is linking to the target_model which is not in the list of
995 models.
996 """
997 out=self.checkLinkTargets(models, hash + [self])
998 for m in self.iterModels():
999 if isinstance(m, Simulation):
1000 out|=m.checkModels(models, hash)
1001 else:
1002 out|=m.checkLinkTargets(models, hash + [self])
1003 return set( [ (str(self)+"."+f[0],f[1]) for f in out ] )
1004
1005
1006 def getSafeTimeStepSize(self,dt):
1007 """
1008 Returns a time step size which can be safely used by all models.
1009
1010 This is the minimum over the time step sizes of all models.
1011 """
1012 out=min([o.getSafeTimeStepSize(dt) for o in self.iterModels()])
1013 return out
1014
1015 def setUp(self):
1016 """
1017 Performs the setup for all models.
1018 """
1019 for o in self.iterModels():
1020 o.setUp()
1021
1022 def doInitialization(self):
1023 """
1024 Initializes all models.
1025 """
1026 for o in self.iterModels():
1027 o.doInitialization()
1028
1029 def doInitialStep(self):
1030 """
1031 Performs an iteration step in the initialization step for all models.
1032 """
1033 iter=0
1034 while not self.terminateInitialIteration():
1035 if iter==0: self.trace("iteration for initialization starts")
1036 iter+=1
1037 self.trace("iteration step %d"%(iter))
1038 for o in self.iterModels():
1039 o.doInitialStep()
1040 if iter>self.MAX_ITER_STEPS:
1041 raise IterationDivergenceError("initial iteration did not converge after %s steps."%iter)
1042 self.trace("Initialization finalized after %s iteration steps."%iter)
1043
1044 def doInitialPostprocessing(self):
1045 """
1046 Finalises the initialization iteration process for all models.
1047 """
1048 for o in self.iterModels():
1049 o.doInitialPostprocessing()
1050
1051 def finalize(self):
1052 """
1053 Returns True if any of the models is to be finalized.
1054 """
1055 return any([o.finalize() for o in self.iterModels()])
1056
1057 def doFinalization(self):
1058 """
1059 Finalises the time stepping for all models.
1060 """
1061 for i in self.iterModels(): i.doFinalization()
1062 self.trace("end of time integation.")
1063
1064 def doStepPreprocessing(self,dt):
1065 """
1066 Initializes the time step for all models.
1067 """
1068 for o in self.iterModels():
1069 o.doStepPreprocessing(dt)
1070
1071 def terminateIteration(self):
1072 """
1073 Returns True if all iterations for all models are terminated.
1074 """
1075 out=all([o.terminateIteration() for o in self.iterModels()])
1076 return out
1077
1078 def terminateInitialIteration(self):
1079 """
1080 Returns True if all initial iterations for all models are terminated.
1081 """
1082 out=all([o.terminateInitialIteration() for o in self.iterModels()])
1083 return out
1084
1085 def doStepPostprocessing(self,dt):
1086 """
1087 Finalises the iteration process for all models.
1088 """
1089 for o in self.iterModels():
1090 o.doStepPostprocessing(dt)
1091 self.time_step+=1
1092 self.time+=dt
1093 self.dt=dt
1094
1095 def doStep(self,dt):
1096 """
1097 Executes the iteration step at a time step for all models::
1098
1099 self.doStepPreprocessing(dt)
1100 while not self.terminateIteration():
1101 for all models:
1102 self.doStep(dt)
1103 self.doStepPostprocessing(dt)
1104 """
1105 self.iter=0
1106 while not self.terminateIteration():
1107 if self.iter==0: self.trace("iteration at %d-th time step %e starts"%(self.time_step+1,self.time+dt))
1108 self.iter+=1
1109 self.trace("iteration step %d"%(self.iter))
1110 for o in self.iterModels():
1111 o.doStep(dt)
1112 if self.iter>0: self.trace("iteration at %d-th time step %e finalized."%(self.time_step+1,self.time+dt))
1113
1114 def run(self,check_pointing=None):
1115 """
1116 Runs the simulation by performing essentially::
1117
1118 self.setUp()
1119 if not restart:
1120 self.doInitialization()
1121 while not self.terminateInitialIteration(): self.doInitialStep()
1122 self.doInitialPostprocessing()
1123 while not self.finalize():
1124 dt=self.getSafeTimeStepSize()
1125 self.doStepPreprocessing(dt_new)
1126 self.doStep(dt_new)
1127 self.doStepPostprocessing(dt_new)
1128 self.doFinalization()
1129
1130 If one of the models throws a ``FailedTimeStepError`` exception a
1131 new time step size is computed through getSafeTimeStepSize() and the
1132 time step is repeated.
1133
1134 If one of the models throws a ``IterationDivergenceError``
1135 exception the time step size is halved and the time step is repeated.
1136
1137 In both cases the time integration is given up after
1138 ``Simulation.FAILED_TIME_STEPS_MAX`` attempts.
1139 """
1140 # check the completeness of the models:
1141 # first a list of all the models involved in the simulation including
1142 # subsimulations:
1143 #
1144 missing=self.checkModels(self.getAllModels(), [])
1145 if len(missing)>0:
1146 msg=""
1147 for l in missing:
1148 msg+="\n\t"+str(l[1])+" at "+l[0]
1149 raise MissingLink("link targets missing in the Simulation: %s"%msg)
1150 #==============================
1151 self.setUp()
1152 if self.time_step < 1:
1153 self.doInitialization()
1154 self.doInitialStep()
1155 self.doInitialPostprocessing()
1156 while not self.finalize():
1157 step_fail_counter=0
1158 iteration_fail_counter=0
1159 if self.time_step==0:
1160 dt_new=self.getSafeTimeStepSize(self.dt)
1161 else:
1162 dt_new=min(max(self.getSafeTimeStepSize(self.dt),self.dt/self.MAX_CHANGE_OF_DT),self.dt*self.MAX_CHANGE_OF_DT)
1163 self.trace("%d. time step %e (step size %e.)" % (self.time_step+1,self.time+dt_new,dt_new))
1164 end_of_step=False
1165 while not end_of_step:
1166 end_of_step=True
1167 if not dt_new>0:
1168 raise NonPositiveStepSizeError("non-positive step size in step %d"%(self.time_step+1))
1169 try:
1170 self.doStepPreprocessing(dt_new)
1171 self.doStep(dt_new)
1172 self.doStepPostprocessing(dt_new)
1173 except IterationDivergenceError:
1174 dt_new*=0.5
1175 end_of_step=False
1176 iteration_fail_counter+=1
1177 if iteration_fail_counter>self.FAILED_TIME_STEPS_MAX:
1178 raise SimulationBreakDownError("reduction of time step to achieve convergence failed after %s steps."%self.FAILED_TIME_STEPS_MAX)
1179 self.trace("Iteration failed. Time step is repeated with new step size %s."%dt_new)
1180 except FailedTimeStepError:
1181 dt_new=self.getSafeTimeStepSize(self.dt)
1182 end_of_step=False
1183 step_fail_counter+=1
1184 self.trace("Time step is repeated with new time step size %s."%dt_new)
1185 if step_fail_counter>self.FAILED_TIME_STEPS_MAX:
1186 raise SimulationBreakDownError("Time integration is given up after %d attempts."%step_fail_counter)
1187 if not check_pointing==None:
1188 if check_pointing.doDump():
1189 self.trace("check point is created.")
1190 self.writeXML()
1191 self.doFinalization()
1192
1193
1194 def toDom(self, esysxml, node):
1195 """
1196 ``toDom`` method of Simulation class.
1197 """
1198 simulation = esysxml.createElement('Simulation')
1199 esysxml.registerLinkableObject(self, simulation)
1200 for rank, sim in enumerate(self.iterModels()):
1201 component = esysxml.createElement('Component')
1202 component.setAttribute('rank', str(rank))
1203 sim.toDom(esysxml, component)
1204 simulation.appendChild(component)
1205 node.appendChild(simulation)
1206
1207 def fromDom(cls, esysxml, node):
1208 sims = []
1209 for n in node.childNodes:
1210 if isinstance(n, minidom.Text):
1211 continue
1212 sims.append(esysxml.getComponent(n))
1213 sims.sort(key=_cmpkey)
1214 sim=cls([s[1] for s in sims], debug=esysxml.debug)
1215 esysxml.registerLinkableObject(sim, node)
1216 return sim
1217
1218 fromDom = classmethod(fromDom)
1219
1220 def _cmpkey(a):
1221 return a[0]
1222
1223 def _comp(a,b):
1224 if a[0]<a[1]:
1225 return 1
1226 elif a[0]>a[1]:
1227 return -1
1228 else:
1229 return 0
1230
1231 class IterationDivergenceError(Exception):
1232 """
1233 Exception which is thrown if there is no convergence of the iteration
1234 process at a time step.
1235
1236 But there is a chance that a smaller step could help to reach convergence.
1237 """
1238 pass
1239
1240 class FailedTimeStepError(Exception):
1241 """
1242 Exception which is thrown if the time step fails because of a step
1243 size that has been chosen to be too large.
1244 """
1245 pass
1246
1247 class SimulationBreakDownError(Exception):
1248 """
1249 Exception which is thrown if the simulation does not manage to progress
1250 in time.
1251 """
1252 pass
1253
1254 class NonPositiveStepSizeError(Exception):
1255 """
1256 Exception which is thrown if the step size is not positive.
1257 """
1258 pass
1259
1260 class MissingLink(Exception):
1261 """
1262 Exception which is thrown when a link is missing.
1263 """
1264 pass
1265
1266 class DataSource(object):
1267 """
1268 Class for handling data sources, including local and remote files.
1269 This class is under development.
1270 """
1271
1272 def __init__(self, uri="file.ext", fileformat="unknown"):
1273 self.uri = uri
1274 self.fileformat = fileformat
1275
1276 def toDom(self, esysxml, node):
1277 """
1278 ``toDom`` method of DataSource. Creates a DataSource node and appends
1279 it to the current XML esysxml.
1280 """
1281 ds = esysxml.createElement('DataSource')
1282 ds.appendChild(esysxml.createDataNode('URI', self.uri))
1283 ds.appendChild(esysxml.createDataNode('FileFormat', self.fileformat))
1284 node.appendChild(ds)
1285
1286 def fromDom(cls, esysxml, node):
1287 uri= str(node.getElementsByTagName("URI")[0].firstChild.nodeValue.strip())
1288 fileformat= str(node.getElementsByTagName("FileFormat")[0].firstChild.nodeValue.strip())
1289 ds = cls(uri, fileformat)
1290 return ds
1291
1292 def getLocalFileName(self):
1293 return self.uri
1294
1295 fromDom = classmethod(fromDom)
1296
1297 class RestartManager(object):
1298 """
1299 A restart manager which does two things: it decides when restart files
1300 were created (when doDump returns true) and manages directories for
1301 restart files. The method getNewDumper creates a new directory and
1302 returns its name.
1303
1304 This restart manager will decide to dump restart files every dump_step
1305 calls of doDump or if more than dump_time since the last dump has
1306 elapsed. The restart manager controls two directories for dumping
1307 restart data, namely for the current and previous dump. This way the
1308 previous dump can be used for restart in the case the current dump failed.
1309
1310 :cvar SEC: unit of seconds, for instance 5*RestartManager.SEC defines 5 seconds
1311 :cvar MIN: unit of minutes, for instance 5*RestartManager.MIN defines 5 minutes
1312 :cvar H: unit of hours, for instance 5*RestartManager.H defines 5 hours
1313 :cvar D: unit of days, for instance 5*RestartManager.D defines 5 days
1314 """
1315 SEC=1.
1316 MIN=60.
1317 H=360.
1318 D=8640.
1319 def __init__(self,dump_time=1080., dump_step=None, dumper=None):
1320 """
1321 Initializes the RestartManager.
1322
1323 :param dump_time: defines the minimum time interval in SEC between two
1324 dumps. If ``None``, time is not used as criterion.
1325 :param dump_step: defines the number of calls of doDump between two
1326 dump events. If ``None``, the call counter is not
1327 used as criterion.
1328 :param dumper: defines the directory for dumping restart files.
1329 Additionally, the directories dumper+"_bkp" and
1330 dumper+"_bkp2" are used. If the directory does not
1331 exist it is created. If dumper is not present a unique
1332 directory within the current working directory is used.
1333 """
1334 self.__dump_step=dump_time
1335 self.__dump_time=dump_step
1336 self.__counter=0
1337 self.__saveMarker()
1338 if dumper == None:
1339 self.__dumper="restart"+str(os.getpid())
1340 else:
1341 self.__dumper=dumper
1342 self.__dumper_bkp=self.__dumper+"_bkp"
1343 self.__dumper_bkp2=self.__dumper+"_bkp2"
1344 self.__current_dumper=None
1345
1346 def __saveMarker(self):
1347 self.__last_restart_time=time.time()
1348 self.__last_restart_counter=self.__counter
1349
1350 def getCurrentDumper(self):
1351 """
1352 Returns the name of the currently used dumper.
1353 """
1354 return self.__current_dumper
1355
1356 def doDump(self):
1357 """
1358 Returns true if the restart should be dumped. Use ``getNewDumper`` to
1359 retrieve the directory name to be used for dumping.
1360 """
1361 if self.__dump_step == None:
1362 if self.__dump_step == None:
1363 out = False
1364 else:
1365 out = (self.__dump_step + self.__last_restart_counter) <= self.__counter
1366 else:
1367 if dump_step == None:
1368 out = (self.__last_restart_time + self.__dump_time) <= time.time()
1369 else:
1370 out = ( (self.__dump_step + self.__last_restart_counter) <= self.__counter) \
1371 or ( (self.__last_restart_time + self.__dump_time) <= time.time() )
1372 if out: self.__saveMarker()
1373 self__counter+=1
1374
1375 def getNewDumper(self):
1376 """
1377 Creates a new directory to be used for dumping and returns its name.
1378 """
1379 if os.access(self.__dumper_bkp,os.F_OK):
1380 if os.access(self.__dumper_bkp2, os.F_OK):
1381 raise RunTimeError("please remove %s."%self.__dumper_bkp2)
1382 try:
1383 os.rename(self.__dumper_bkp, self.__dumper_bkp2)
1384 except:
1385 self.__current_dumper=self.__dumper
1386 raise RunTimeError("renaming backup directory %s failed. Use %s for restart."%(self.__dumper_bkp,self.__dumper))
1387 if os.access(self.__dumper,os.F_OK):
1388 if os.access(self.__dumper_bkp, os.F_OK):
1389 raise RunTimeError("please remove %s."%self.__dumper_bkp)
1390 try:
1391 os.rename(self.__dumper, self.__dumper_bkp)
1392 except:
1393 self.__current_dumper=self.__dumper_bkp2
1394 raise RunTimeError("moving directory %s to backup failed. Use %s for restart."%(self.__dumper,self.__dumper_bkp2))
1395 try:
1396 os.mkdir(self.__dumper)
1397 except:
1398 self.__current_dumper=self.__dumper_bkp
1399 raise RunTimeError("creating a new restart directory %s failed. Use %s for restart."%(self.__dumper,self.__dumper_bkp))
1400 if os.access(self.__dumper_bkp2, os.F_OK): os.rmdir(self.__dumper_bkp2)
1401 return self.getCurrentDumper()
1402
1403
1404 # 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