/[escript]/branches/symbolic_from_3470/escript/py_src/modelframe.py
ViewVC logotype

Contents of /branches/symbolic_from_3470/escript/py_src/modelframe.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3789 - (show annotations)
Tue Jan 31 04:55:05 2012 UTC (7 years, 2 months ago) by caltinay
File MIME type: text/x-python
File size: 47713 byte(s)
Fast forward to latest trunk revision 3788.

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