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

Diff of /trunk/esys2/escript/py_src/modelframe.py

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 122 by jgs, Thu Jun 9 05:38:05 2005 UTC revision 126 by jgs, Fri Jul 22 03:53:08 2005 UTC
# Line 2  Line 2 
2    
3  from types import StringType,IntType,FloatType,BooleanType,ListType,DictType  from types import StringType,IntType,FloatType,BooleanType,ListType,DictType
4  from sys import stdout  from sys import stdout
5    import itertools
6    
7  from xml.dom import minidom  import modellib
8    
9    # import the 'set' module if it's not defined (python2.3/2.4 difference)
10    try:
11        set
12    except NameError:
13        from sets import Set as set
14    
15  _HEADER_START="<?xml version=\"1.0\"?><ESyS>\n"  from xml.dom import minidom
 _HEADER_END="<\ESyS>\n"  
16    
17  def dataNode(document, tagName, data):  def dataNode(document, tagName, data):
18      t = document.createTextNode(str(data))      """
19      n = document.createElement(tagName)      dataNodes are the building blocks of the xml documents constructed in
20      n.appendChild(t)      this module. document is the current xml document, tagName is the
21      return n      associated xml tag, and data is the values in the tag.
22        """
23        t = document.createTextNode(str(data))
24        n = document.createElement(tagName)
25        n.appendChild(t)
26        return n
27    
28  def esysDoc():  def esysDoc():
29      doc = minidom.Document()      """
30      esys = doc.createElement('ESys')      Global method for creating an instance of an EsysXML document.
31      doc.appendChild(esys)      """
32      return doc, esys      doc = minidom.Document()
33        esys = doc.createElement('ESys')
34        doc.appendChild(esys)
35        return doc, esys
36    
37    def all(seq):
38        for x in seq:
39            if not x:
40                return False
41        return True
42    
43    LinkableObjectRegistry = {}
44    
45    def registerLinkableObject(obj_id, o):
46        LinkableObjectRegistry[obj_id] = o
47    
48    LinkRegistry = []
49    
50    def registerLink(obj_id, l):
51        LinkRegistry.append((obj_id,l))
52    
53    def parse(xml):
54        global LinkRegistry, LinkableObjectRegistry
55        LinkRegistry = []
56        LinkableObjectRegistry = {}
57    
58        doc = minidom.parseString(xml)
59    
60        sim = getComponent(doc.firstChild)
61        for obj_id, link in LinkRegistry:
62            link.target = LinkableObjectRegistry[obj_id]
63    
64        return sim
65    
66    def getComponent(doc):
67        for node in doc.childNodes:
68            if isinstance(node, minidom.Element):
69                if node.tagName == 'Simulation':
70                    if node.getAttribute("type") == 'Simulation':
71                        return Simulation.fromDom(node)
72                    elif node.getAttribute("type") == 'ExplicitSimulation':
73                        return ExplicitSimulation.fromDom(node)
74                if node.tagName == 'Model':
75                    model_type = node.getAttribute("type")
76                    model_subclasses = Model.__subclasses__()
77                    for model in model_subclasses:
78                        if model_type == model.__name__:
79                            return model.fromDom(node)
80                        
81                raise "Invalid simulation type, %r" % node.getAttribute("type")
82    
83        raise ValueError("No Simulation Found")
84                
85    
86  class Link:  class Link:
87      """      """
88      a Link makes an attribute of an object callable:      a Link makes an attribute of an object callable:
89            o.object()            o.object()
90            o.a=8            o.a=8
91            l=Link(o,"a")            l=Link(o,"a")
92            assert l()==8            assert l()==8
93       """       """
94            
95      def __init__(self,target=None,attribute=None):      def __init__(self,target,attribute=None):
96          """creates a link to the object target. If attribute is given, the link is          """
97          establised to this attribute of the target.  Otherwise the attribute is          creates a link to the object target. If attribute is given, the link is
98          undefined."""          establised to this attribute of the target.  Otherwise the attribute is
99          self.__target=target          undefined.
100          self.setAttributeName(attribute)          """
101                self.target = target
102      def setAttributeName(self,name):          self.attribute = None
103          """set a new attribute name to be collected from the target object. The          self.setAttributeName(attribute)
104          target object must have the attribute with name name."""      
105          if not name==None:      def setAttributeName(self,attribute):
106              if not hasattr(self.__target,name):          """
107                  raise AttributeError("%s: target %s has no attribute %s."%(self,self.__target,name))          set a new attribute name to be collected from the target object. The
108          self.__attribute=name          target object must have the attribute with name attribute.
109                """
110      def hasDefinedAttributeName(self):          if attribute and self.target and not hasattr(self.target, attribute):
111          """returns true if an attribute name is set"""              raise AttributeError("%s: target %s has no attribute %s."%(self, self.target, attribute))
112          if self.__attribute==None:          self.attribute = attribute
113              return False      
114          else:      def hasDefinedAttributeName(self):
115              return True          """
116                returns true if an attribute name is set
117      def __str__(self):          """
118          """returns a string representation of the link"""          return self.attribute != None
119          if self.hasDefinedAttributeName():      
120              return "reference to attribute %s of %s"%(self.__attribute,self.__target)      def __repr__(self):
121          else:          """
122              return "reference to target %s"%self.__target          returns a string representation of the link
123                """
124      def __call__(self,name=None):          if self.hasDefinedAttributeName():
125          """returns the value of the attribute of the target object. If the              return "<Link to attribute %s of %s>" % (self.attribute, self.target)
126          atrribute is callable then the return value of the call is returned."""          else:
127          if not name == None:              return "<Link to target %s>" % self.target
128              out=getattr(self.__target, name)      
129          else:      def __call__(self,name=None):
130              out=getattr(self.__target, self.__attribute)          """
131          if callable(out):          returns the value of the attribute of the target object. If the
132              return out()          atrribute is callable then the return value of the call is returned.
133          else:          """
134              return out          if name:
135                out=getattr(self.target, name)
136      def toDom(self, document, node):          else:
137          """ toDom method of Link. Creates a Link node and appends it to the current XML              out=getattr(self.target, self.attribute)
138          document """  
139          link = document.createElement('link')          if callable(out):
140          link.appendChild(dataNode(document, 'target', id(self.__target)))              return out()
141          link.appendChild(dataNode(document, 'attribute', self.__attribute))          else:
142          node.appendChild(link)              return out
143        
144      def writeXML(self,ostream=stdout):      def toDom(self, document, node):
145          """writes an XML representation of self to the output stream ostream.          """
146          If ostream is nor present the standart output stream is used.  If          toDom method of Link. Creates a Link node and appends it to the current XML
147          esysheader==True the esys XML header is written"""          document
148            """
149            link = document.createElement('Link')
150            link.appendChild(dataNode(document, 'Target', self.target.id))
151            # this use of id will not work for purposes of being able to retrieve the intended
152            # target from the xml later. I need a better unique identifier.
153            assert self.attribute, "You can't xmlify a Link without a target attribute"
154            link.appendChild(dataNode(document, 'Attribute', self.attribute))
155            node.appendChild(link)
156    
157        def fromDom(cls, doc):
158            targetid = doc.getElementsByTagName("Target")[0].firstChild.nodeValue.strip()
159            attribute = doc.getElementsByTagName("Attribute")[0].firstChild.nodeValue.strip()
160            l = cls(None, attribute)
161            registerLink(targetid, l)
162            return l
163    
164        fromDom = classmethod(fromDom)
165        
166        def writeXML(self,ostream=stdout):
167            """
168            writes an XML representation of self to the output stream ostream.
169            If ostream is nor present the standart output stream is used.  If
170            esysheader==True the esys XML header is written
171            """
172    
173          document, rootnode = esysDoc()          document, rootnode = esysDoc()
174          self.toDom(document, rootnode)          self.toDom(document, rootnode)
175    
176          ostream.write(document.toprettyxml())          ostream.write(document.toprettyxml())
177    
178  class LinkableObject(object):  class LinkableObject(object):
179      """an object that allows to link its attributes to attributes of other      """
180      objects via a Link object. For instance      An object that allows to link its attributes to attributes of other objects
181                  via a Link object. For instance
182             p=LinkableObject()            
183             p.x=Link(o,"name")             p = LinkableObject()
184             print p.x             p.x = Link(o,"name")
185                   print p.x
186          links attribute x of p to the attribute name of object o. p.x will print      
187          the current value of attribute name of object o.  if the value is      links attribute x of p to the attribute name of object o.
188          callable p.x will rturn the return value of the call.  
189      """      p.x will contain the current value of attribute name of object o.  
190        
191      def __init__(self,debug=False):      If the value of getattr(o, "name") is callable, p.x will rturn the return
192          """initiates a model from a list of subsimulations. """      value of the call.
193          self.setDebug(debug)      """
194          self.__linked_attributes={}    
195      #      number_sequence = itertools.count(100)
196      # some basic fuctions:      
197      #      def __init__(self, debug=False):
198      def debugOn(self):          """ initializes LinkableObject so that we can operate on Links """
199          """sets debugging to on"""          self.debug = debug
200          self.__debug=True          self.__linked_attributes={}
201      def debugOff(self):          self.id = self.number_sequence.next()
202          """sets debugging to off"""  
203          self.__debug=False      def trace(self, msg):
204      def debug(self):          """ If debugging is on, print the message, otherwise do nothing
205          """returns True if debug mode is set to on"""          """
206          return self.__debug          if self.debug:
207      def setDebug(self,flag=False):              print msg
208          """sets debugging to flag"""      
209          if flag:      def __getattr__(self,name):
210              self.debugOn()          """returns the value of attribute name. If the value is a Link object the
211          else:          object is called and the return value is returned."""
212              self.debugOff()          out = self.getAttributeObject(name)
213            if isinstance(out,Link):
214      def trace(self, msg):              return out()
215          if self.__debug:          else:
216              print msg              return out
217            
218      def __getattr__(self,name):      def getAttributeObject(self,name):
219          """returns the value of attribute name. If the value is a Link object the          """return the object stored for attribute name."""
220          object is called and the return value is returned."""  
221          out=self.getAttributeObject(name)          if self.__dict__.has_key(name):
222          if isinstance(out,Link):              return self.__dict__[name]
223              return out()  
224          else:          if self.__linked_attributes.has_key(name):
225              return out              return self.__linked_attributes[name]
226        
227      def getAttributeObject(self,name):          raise AttributeError,"No attribute %s."%name
228          """return the object stored for attribute name."""      
229          if self.__dict__.has_key(name):      def __setattr__(self,name,value):
230              out=self.__dict__[name]          """sets the value for attribute name. If value is a Link the target
231          elif self.__linked_attributes.has_key(name):          attribute is set to name if no attribute has been specified."""
232              out=self.__linked_attributes[name]  
233          else:  
234              raise AttributeError,"No attribute %s."%name          if self.__dict__.has_key(name):
235          return out              del self.__dict__[name]
236        
237      def __setattr__(self,name,value):          if isinstance(value,Link):
238          """sets the value for attribute name. If value is a Link the target              if not value.hasDefinedAttributeName():
239          attribute is set to name if no attribute has been specified."""                  value.setAttributeName(name)
240          if self.__dict__.has_key(name):              self.__linked_attributes[name] = value
241              del self.__dict__[name]  
242          if isinstance(value,Link):              self.trace("DEBUG: %s: attribute %s is now linked by %s."%(self,name,value))
243             if not value.hasDefinedAttributeName():          else:
244                 value.setAttributeName(name)              self.__dict__[name] = value
245             self.__linked_attributes[name]=value      
246             if self.debug():      def __delattr__(self,name):
247                 print "DEBUG: %s: attribute %s is now linked by %s."%(self,name,value)          """removes the attribute name."""
248          else:  
249             self.__dict__[name]=value          if self.__linked_attributes.has_key[name]:
250                    del self.__linked_attributes[name]
251      def __delattr__(self,name):          elif self.__dict__.has_key(name):
252          """removes the attribute name."""              del self.__dict__[name]
253          if self.__linked_attributes.has_key[name]:          else:
254              del self.__linked_attributes[name]              raise AttributeError,"No attribute %s."%name
         elif self.__dict__.has_key(name):  
             del self.__dict__[name]  
         else:  
             raise AttributeError,"No attribute %s."%name  
255    
256  class SimulationFrame(LinkableObject):  class SimulationFrame(LinkableObject):
257      """A SimulationFrame objects represents a processess marching over time      """A SimulationFrame objects represents a processess marching over time
258      until a finalizing condition is fullfilled.  At each time step an iterative      until a finalizing condition is fullfilled.  At each time step an iterative
259      process can be performed and the time step size can be controlled      process can be performed and the time step size can be controlled
260      """      """
261      UNDEF_DT=1.e300      UNDEF_DT=1.e300
262      MAX_TIME_STEP_REDUCTION=20      MAX_TIME_STEP_REDUCTION=20
263      MAX_ITER_STEPS=50      MAX_ITER_STEPS=50
264            
265      def __init__(self,debug=False):      def __init__(self,**kwargs):
266          """initiates a simulation"""          """
267          LinkableObject.__init__(self,debug)          Initialises a simulation
268                
269      def doInitialization(self,t):          Just calls the parent constructor.
270          """initializes the time stepping scheme. This function may be          """
271          overwritten."""          LinkableObject.__init__(self,**kwargs)
272          pass      
273            def doInitialization(self,t):
274      def getSafeTimeStepSize(self,dt):          """initializes the time stepping scheme. This function may be
275          """returns a time step size which can savely be used. This function may          overwritten."""
276          be overwritten."""          pass
277          return self.UNDEF_DT      
278            def getSafeTimeStepSize(self,dt):
279      def finalize(self):          """returns a time step size which can savely be used. This function may
280          """returns True if the time stepping is finalized. This function may be          be overwritten."""
281          overwritten."""          return self.UNDEF_DT
282          return True      
283              def finalize(self):
284      def doFinalization(self):          """returns True if the time stepping is finalized. This function may be
285          """finalizes the time stepping.  This function may be overwritten."""          overwritten."""
286          pass          return True
287              
288      def doIterationInitialization(self,dt):      def doFinalization(self):
289          """initializes the iteration at time step t. This function may be          """finalizes the time stepping.  This function may be overwritten."""
290          overwritten. (only called if doStep is not overwritten)"""          pass
291          pass      
292            def doIterationInitialization(self,dt):
293      def doIterationStep(self,dt):          """initializes the iteration at time step t. This function may be
294          """executes the iteration step. This function may be overwritten. (only          overwritten. (only called if doStep is not overwritten)"""
295          called if doStep is not overwritten)"""          pass
296          pass      
297            def doIterationStep(self,dt):
298      def terminate(self):          """executes the iteration step. This function may be overwritten. (only
299          """returns True if iteration on a time step is terminated. (only called          called if doStep is not overwritten)"""
300          if doStep is not overwritten)"""          pass
301          return True      
302              def terminate(self):
303      def doIterationFinalization(self,dt):          """returns True if iteration on a time step is terminated. (only called
304          """finalalizes the iteration process. (only called if doStep is not          if doStep is not overwritten)"""
305          overwritten)"""          return True
306          pass        
307            def doIterationFinalization(self,dt):
308      def run(self,check_point=None):          """finalalizes the iteration process. (only called if doStep is not
309          """run the simulation by performing essentially          overwritten)"""
310                pass
311              self.doInitialization()      
312              while not self.finalize():      def run(self,check_point=None):
313                 dt=self.getSafeTimeStepSize()          """run the simulation by performing essentially
314                 self.doStep(dt)      
315                 if n%check_point==0: self.writeXML()              self.doInitialization()
316              self.doFinalization()              while not self.finalize():
317                       dt=self.getSafeTimeStepSize()
318          """                 self.doStep(dt)
319          self.__tn=0.                 if n%check_point==0: self.writeXML()
320          self.__n=0              self.doFinalization()
321          self.__dt=None      
322          self.doInitialization(self.__tn)          """
323          while not self.finalize():          self.__tn=0.
324              self.__n+=1          self.__n=0
325              self.__dt=self.getSafeTimeStepSize(self.__dt)          self.__dt=None
326              if self.__dt==None: self.__dt=self.UNDEF_DT          self.doInitialization(self.__tn)
327              if not self.__dt>0:          while not self.finalize():
328                   raise NonPositiveStepSizeError("non-positive step size in step %d",self.__n)              self.__n+=1
329              if self.debug(): print "%s: %d. time step %e (step size %e.)"%(self,self.__n,self.__tn+self.__dt,self.__dt)              self.__dt=self.getSafeTimeStepSize(self.__dt)
330              endoftimestep=False              if self.__dt==None: self.__dt=self.UNDEF_DT
331              failcounter=0              if not self.__dt>0:
332              while not endoftimestep:                   raise NonPositiveStepSizeError("non-positive step size in step %d",self.__n)
333                  endoftimestep=True              self.trace("%s: %d. time step %e (step size %e.)" %
334                  try:                            (self,self.__n,self.__tn+self.__dt,self.__dt))
335                      self.doStep(self.__dt)              endoftimestep=False
336                  except FailedTimeStepError:              failcounter=0
337                      self.__dt=self.getSafeTimeStepSize(self.__dt)              while not endoftimestep:
338                      if self.__dt==None: self.__dt=self.UNDEF_DT                  endoftimestep=True
339                      endoftimestep=False                  try:
340                      if self.debug(): print "%s: time step is repeated with new step size %e."%(self,self.__dt)                      self.doStep(self.__dt)
341                  except IterationDivergenceError:                  except FailedTimeStepError:
342                      self.__dt*=0.5                      self.__dt=self.getSafeTimeStepSize(self.__dt)
343                      endoftimestep=False                      if self.__dt==None: self.__dt=self.UNDEF_DT
344                      failcounter+=1                      endoftimestep=False
345                      if failcounter>self.MAX_TIME_STEP_REDUCTION:                      self.trace("%s: time step is repeated with new step size %e."%(self,self.__dt))
346                          raise IterationBreakDownError("reduction of time step to achieve convergence failed.")                  except IterationDivergenceError:
347                        self.__dt*=0.5
348                      self.trace("%s: iteration failes. time step is repeated with new step size %e."                      endoftimestep=False
349                          % (self,self.__dt))                      failcounter+=1
350              self.__tn+=self.__dt                      if failcounter>self.MAX_TIME_STEP_REDUCTION:
351              if not check_point==None:                          raise IterationBreakDownError("reduction of time step to achieve convergence failed.")
352                  if self.__n%check_point==0:  
353                      self.trace("%s: check point is created."%self)                      self.trace("%s: iteration failes. time step is repeated with new step size %e."
354                      self.writeXML()                          % (self,self.__dt))
355          self.doFinalization()              self.__tn+=self.__dt
356                if not check_point==None:
357      def writeXML(self):                  if self.__n%check_point==0:
358          raise RuntimeError, "Not implemented"                      self.trace("%s: check point is created."%self)
359                            self.writeXML()
360      def doStep(self,dt):          self.doFinalization()
361          """executes a time step by iteration. This function may be overwritten.  
362                  def writeXML(self):
363             basicly it performs :          raise RuntimeError, "Not implemented"
364            
365              self.doIterationInitialization(dt)      def doStep(self,dt):
366              while not self.terminate(): self.doIterationStep(dt)          """executes a time step by iteration. This function may be overwritten.
367              self.doIterationFinalization(dt)            
368                   basicly it performs :
369          """      
370          self.doIterationInitialization(dt)              self.doIterationInitialization(dt)
371          self.__iter=0              while not self.terminate(): self.doIterationStep(dt)
372          while not self.terminate():              self.doIterationFinalization(dt)
373              if self.debug(): print "%s: iteration step %d"%(self,self.__iter)      
374              self.doIterationStep(dt)          """
375              self.__iter+=1          self.doIterationInitialization(dt)
376              if self.__iter>self.MAX_ITER_STEPS:          self.__iter=0
377                  raise IterationDivergenceError("%s: divergence in step %d"%(self,self.__iter))          while not self.terminate():
378          self.doIterationFinalization(dt)              self.trace("%s: iteration step %d"%(self,self.__iter))
379                self.doIterationStep(dt)
380                self.__iter+=1
381                if self.__iter>self.MAX_ITER_STEPS:
382                    raise IterationDivergenceError("%s: divergence in step %d"%(self,self.__iter))
383            self.doIterationFinalization(dt)
384    
385  class Simulation(SimulationFrame):  class Simulation(SimulationFrame):
386      """A Simulation object is comprised by SimulationFrame(s) called subsimulations."""      """A Simulation object is comprised by SimulationFrame(s) called subsimulations."""
387            
388      def __init__(self,subsimulations=[],debug=False):      def __init__(self, subsimulations=[], **kwargs):
389          """initiates a simulation from a list of subsimulations. """          """initiates a simulation from a list of subsimulations. """
390          SimulationFrame.__init__(self,debug)          SimulationFrame.__init__(self, **kwargs)
391          self.__subsimulations=[]          self.__subsimulations=[]
392          for i in range(len(subsimulations)): self[i]=subsimulations[i]  
393                for i in range(len(subsimulations)):
394      def iterSubsimulations(self):              self[i] = subsimulations[i]
395          """returns an iterator over the subsimulations"""      
396          return self.__subsimulations      def iterSubsimulations(self):
397                """returns an iterator over the subsimulations"""
398      def __getitem__(self,i):          return self.__subsimulations
399          """returns the i-th subsimulation"""      
400          return self.__subsimulations[i]      def __getitem__(self,i):
401                """returns the i-th subsimulation"""
402      def __setitem__(self,i,value):          return self.__subsimulations[i]
403          """sets the i-th subsimulation"""      
404          if not isinstance(value,SimulationFrame):      def __setitem__(self,i,value):
405              raise ValueError("assigned value is not a Simulation")          """sets the i-th subsimulation"""
406          for j in range(max(i-len(self.__subsimulations)+1,0)): self.__subsimulations.append(None)          if not isinstance(value,SimulationFrame):
407          self.__subsimulations[i]=value              raise ValueError("assigned value is not a Simulation")
408                for j in range(max(i-len(self.__subsimulations)+1,0)): self.__subsimulations.append(None)
409      def __len__(self):          self.__subsimulations[i]=value
410          """returns the number of subsimulations"""      
411          return len(self.__subsimulations)      def __len__(self):
412            """returns the number of subsimulations"""
413      def toDom(self, document, node):          return len(self.__subsimulations)
414          """ toDom method of Simualtion class """  
415          simulation = document.createElement('Simulation')      def toDom(self, document, node):
416          simulation.setAttribute('type', self.__class__.__name__)          """ toDom method of Simulation class """
417            simulation = document.createElement('Simulation')
418          for rank, sim in enumerate(self.iterSubsimulations()):          simulation.setAttribute('type', self.__class__.__name__)
419              component = document.createElement('Component')  
420              component.setAttribute('rank', str(rank))          for rank, sim in enumerate(self.iterSubsimulations()):
421                component = document.createElement('Component')
422              sim.toDom(document, component)              component.setAttribute('rank', str(rank))
423    
424              simulation.appendChild(component)              sim.toDom(document, component)
425    
426          node.appendChild(simulation)              simulation.appendChild(component)
427    
428      def writeXML(self,ostream=stdout):          node.appendChild(simulation)
429          """writes the object as an XML object into an output stream"""  
430          document, rootnode = esysDoc()      def writeXML(self,ostream=stdout):
431          self.toDom(document, rootnode)          """writes the object as an XML object into an output stream"""
432          ostream.write(document.toprettyxml())          document, rootnode = esysDoc()
433                self.toDom(document, rootnode)
434      def getSafeTimeStepSize(self,dt):          ostream.write(document.toprettyxml())
435          """returns a time step size which can savely be used by all subsimulations"""      
436          out=self.UNDEF_DT      def getSafeTimeStepSize(self,dt):
437          for i in self.iterSubsimulations():          """returns a time step size which can safely be used by all subsimulations"""
438              dt_new=i.getSafeTimeStepSize(dt)          out=self.UNDEF_DT
439              if dt_new!=None: out=min(out,dt_new)          for o in self.iterSubsimulations():
440          return out              dt_new = o.getSafeTimeStepSize(dt)
441                    if dt_new != None:
442      def doInitialization(self,dt):                  out = min(out,dt_new)
443          """initializes all subsimulations """  
444          for i in self.iterSubsimulations(): i.doInitialization(dt)          return out
445            
446      def finalize(self):      def doInitialization(self,dt):
447          """returns True if all subsimulations are finalized"""          """initializes all subsimulations """
448          out=True          for o in self.iterSubsimulations():
449          for i in self.iterSubsimulations(): out = out and i.finalize()              o.doInitialization(dt)
450          return out      
451              def finalize(self):
452      def doFinalization(self):          """returns True if all subsimulations are finalized"""
453          """finalalizes the time stepping for all subsimulations."""          return all([o.finalize() for o in self.iterSubsimulations()])
454          for i in self.iterSubsimulations(): i.doFinalization()        
455            def doFinalization(self):
456      def doIterationInitialization(self,dt):          """finalalizes the time stepping for all subsimulations."""
457          """initializes the iteration at time t for all subsimulations."""          for i in self.iterSubsimulations(): i.doFinalization()
458          self.__iter=0      
459          if self.debug(): print "%s: iteration starts"%self      def doIterationInitialization(self,dt):
460          for i in self.iterSubsimulations(): i.doIterationInitialization(dt)          """initializes the iteration at time t for all subsimulations."""
461                self.__iter=0
462      def terminate(self):          self.trace("%s: iteration starts"%self)
463          """returns True if all iterations for all subsimulations are terminated."""  
464          out=True          for o in self.iterSubsimulations():
465          for i in self.iterSubsimulations(): out = out and i.terminate()              o.doIterationInitialization(dt)
466          return out      
467              def terminate(self):
468      def doIterationFinalization(self,dt):          """returns True if all iterations for all subsimulations are terminated."""
469          """finalalizes the iteration process for each of the subsimulations."""          return all([o.terminate() for o in self.iterSubsimulations()])
470          for i in self.iterSubsimulations(): i.doIterationFinalization(dt)        
471          if self.debug(): print "%s: iteration finalized after %s steps"%(self,self.__iter+1)      def doIterationFinalization(self,dt):
472                """finalalizes the iteration process for each of the subsimulations."""
473      def doIterationStep(self,dt):          for o in self.iterSubsimulations():
474          """executes the iteration step at time step for each subsimulation"""              o.doIterationFinalization(dt)
475          self.__iter+=1          self.trace("%s: iteration finalized after %s steps"%(self,self.__iter+1))
476          if self.debug(): print "%s: iteration step %d"%(self,self.__iter)      
477          for i in self.iterSubsimulations(): i.doIterationStep(dt)      def doIterationStep(self,dt):
478            """executes the iteration step at time step for each subsimulation"""
479            self.__iter+=1
480            self.trace("%s: iteration step %d"%(self,self.__iter))
481            for o in self.iterSubsimulations():
482                o.doIterationStep(dt)
483    
484        def fromDom(cls, doc):
485            """
486            Needs to be filled out.
487            """
488            sims = []
489            for node in doc.childNodes:
490                if isinstance(node, minidom.Text):
491                    continue
492    
493                sims.append(getComponent(node))
494    
495    
496            return cls(sims)
497    
498    
499                
500    
501    
502        fromDom = classmethod(fromDom)
503    
504  class ExplicitSimulation(Simulation):  class ExplicitSimulation(Simulation):
505      """This is a modified form of the Simulation class. In fact it overwrites      """This is a modified form of the Simulation class. In fact it overwrites
506      the doStep method by executing the doStep method of all subsimulations      the doStep method by executing the doStep method of all subsimulations
507      rather then iterating over all subsimulations."""      rather then iterating over all subsimulations."""
508    
509      def doStep(self,dt):      def doStep(self,dt):
510          """executes the time step for all subsimulation"""          """executes the time step for all subsimulation"""
511          for i in self.iterSubsimulations(): i.doStep(dt)          for i in self.iterSubsimulations():
512                i.doStep(dt)
513    
514  class _ParameterIterator:  class _ParameterIterator:
515      def __init__(self,parameterset):      def __init__(self,parameterset):
516          self.__set=parameterset  
517          self.__iter=iter(parameterset.getParameterList())          self.__set=parameterset
518      def next(self):          self.__iter=iter(parameterset.parameters)
519          o=self.__iter.next()  
520          return (o,self.__set.getAttributeObject(o))      def next(self):
521      def __iter__(self):          o=self.__iter.next()
522          return self          return (o,self.__set.getAttributeObject(o))
523    
524        def __iter__(self):
525            return self
526    
527  class ParameterSet(LinkableObject):  class ParameterSet(LinkableObject):
528      """a class which allows to emphazise attributes to be written and read to XML      """a class which allows to emphazise attributes to be written and read to XML
529                
530         Leaves of  an ESySParameters objects can be         Leaves of  an ESySParameters objects can be
531            
532              a real number              a real number
533              a integer number              a integer number
534              a string              a string
535              a boolean value              a boolean value
536              a ParameterSet object              a ParameterSet object
537              a Simulation object              a Simulation object
538              a Model object              a Model object
539              any other object (not considered by writeESySXML and writeXML)              any other object (not considered by writeESySXML and writeXML)
540            
541             Example how to create an ESySParameters object:             Example how to create an ESySParameters object:
542            
543                   p11=ParameterSet(gamma1=1.,gamma2=2.,gamma3=3.)                   p11=ParameterSet(gamma1=1.,gamma2=2.,gamma3=3.)
544                   p1=ParameterSet(dim=2,tol_v=0.001,output_file="/tmp/u.%3.3d.dx",runFlag=True,parm11=p11)                   p1=ParameterSet(dim=2,tol_v=0.001,output_file="/tmp/u.%3.3d.dx",runFlag=True,parm11=p11)
545                   parm=ParameterSet(parm1=p1,parm2=ParameterSet(alpha=Link(p11,"gamma1")))                   parm=ParameterSet(parm1=p1,parm2=ParameterSet(alpha=Link(p11,"gamma1")))
546            
547             This can be accessed as             This can be accessed as
548            
549                   parm.parm1.gamma=0.                   parm.parm1.gamma=0.
550                   parm.parm1.dim=2                   parm.parm1.dim=2
551                   parm.parm1.tol_v=0.001                   parm.parm1.tol_v=0.001
552                   parm.parm1.output_file="/tmp/u.%3.3d.dx"                   parm.parm1.output_file="/tmp/u.%3.3d.dx"
553                   parm.parm1.runFlag=True                   parm.parm1.runFlag=True
554                   parm.parm1.parm11.gamma1=1.                   parm.parm1.parm11.gamma1=1.
555                   parm.parm1.parm11.gamma2=2.                   parm.parm1.parm11.gamma2=2.
556                   parm.parm1.parm11.gamma3=3.                   parm.parm1.parm11.gamma3=3.
557                   parm.parm2.alpha=1. (value of parm.parm1.parm11.gamma1)                   parm.parm2.alpha=1. (value of parm.parm1.parm11.gamma1)
558                            
559      """      """
560      def __init__(self,parameters=[]):      def __init__(self, parameters=[], **kwargs):
561          """creates a ParameterSet with parameters parameters"""          """creates a ParameterSet with parameters parameters"""
562          LinkableObject.__init__(self,debug=False)          LinkableObject.__init__(self, **kwargs)
563          self.__parameters=[]          self.parameters = set()
564          self.declareParameters(parameters)          self.declareParameters(parameters)
565            
566      def getParameterList(self):      def declareParameter(self,**parameters):
567          """returns the list of parameters"""          """declares a new parameter(s) and its (their) inital value."""
568          return self.__parameters          self.declareParameters(parameters)
569            
570      def isParameter(self,name):      def declareParameters(self,parameters):
571          """returns true if name is a parameter"""          """declares a set of parameters. parameters can be a list, a dictonary or a ParameterSet."""
572          return name in self.getParameterList()          if isinstance(parameters,ListType):
573                    parameters = zip(parameters, itertools.repeat(None))
574      def declareParameter(self,**parameters):          if isinstance(parameters,DictType):
575          """declares a new parameter(s) and its (their) inital value."""              parameters = parameters.iteritems()
576          self.declareParameters(parameters)  
577                for prm, value in parameters:
578      def declareParameters(self,parameters):              setattr(self,prm,value)
579          """declares a set of parameters. parameters can be a list, a dictonary or a ParameterSet."""              self.parameters.add(prm)
580          if isinstance(parameters,ListType):  
581              for prm in parameters:              self.trace("%s: parameter %s has been declared."%(self,prm))
582                  setattr(self,prm,None)  
583                  if not self.isParameter(prm): self.getParameterList().append(prm)      def releaseParameters(self,name):
584                  if self.debug(): print "%s: parameter %s has been declared."%(self,prm)          """removes parameter name from the paramameters"""
585          elif isinstance(parameters,DictType):          if self.isParameter(name):
586              for prm,value in parameters.iteritems():              self.parameters.remove(name)
587                  setattr(self,prm,value)              self.debug("%s: parameter %s has been removed."%(self, name))
588                  if not self.isParameter(prm): self.getParameterList().append(prm)        
589                  if self.debug(): print "%s: parameter %s has been declared."%(self,prm)      def __iter__(self):
590          else:          """creates an iterator over the parameter and their values"""
591              for prm,value in parameters:          return _ParameterIterator(self)
592                  setattr(self,prm,value)      
593                  if not self.isParameter(prm): self.getParameterList().append(prm)        def showParameters(self):
594                  if self.debug(): print "%s: parameter %s has been declared."%(self,prm)          """returns a descrition of the parameters"""        
595                out="{"
596      def releaseParameters(self,name):          notfirst=False
597          """removes parameter name from the paramameters"""          for i,v in self:
598          if self.isParameter(name):              if notfirst: out=out+","
599              self.getParameterList().remove(name)              notfirst=True
600              if self.debug(): print "%s: parameter %s has been removed."%(self,prm)              if isinstance(v,ParameterSet):
601                        out="%s\"%s\" : %s"%(out,i,v.showParameters())
602      def __iter__(self):              else:
603          """creates an iterator over the parameter and their values"""                  out="%s\"%s\" : %s"%(out,i,v)
604          return _ParameterIterator(self)          return out+"}"
605            
606      def showParameters(self):      def __delattr__(self,name):
607          """returns a descrition of the parameters"""                  """removes the attribute name."""
608          out="{"          LinkableObject.__delattr__(self,name)
609          notfirst=False          try:
610          for i,v in self:              self.releaseParameter(name)
611              if notfirst: out=out+","          except:
612              notfirst=True              pass
613              if isinstance(v,ParameterSet):  
614                  out="%s\"%s\" : %s"%(out,i,v.showParameters())      def toDom(self, document, node):
615              else:          """ toDom method of ParameterSet class """
616                  out="%s\"%s\" : %s"%(out,i,v)          pset = document.createElement('ParameterSet')
617          return out+"}"          node.appendChild(pset)
618                self._parametersToDom(document, pset)
619      def __delattr__(self,name):  
620          """removes the attribute name."""      def _parametersToDom(self, document, node):
621          LinkableObject.__delattr__(self,name)          node.setAttribute ('id', str(self.id))
622          try:          for name,value in self:
623              self.releaseParameter(name)              param = document.createElement('Parameter')
624          except:              param.setAttribute('type', value.__class__.__name__)
625              pass  
626                param.appendChild(dataNode(document, 'Name', name))
627      def toDom(self, document, node):  
628          """ toDom method of ParameterSet class """              val = document.createElement('Value')
629          pset = document.createElement('ParameterSet')  
630          node.appendChild(pset)              if isinstance(value,ParameterSet):
631          self._parametersToDom(document, pset)                  value.toDom(document, val)
632                    param.appendChild(val)
633      def _parametersToDom(self, document, node):              elif isinstance(value, Link):
634          for name,value in self:                  value.toDom(document, val)
635              param = document.createElement('Parameter')                  param.appendChild(val)
636              param.setAttribute('type', value.__class__.__name__)              elif isinstance(value,StringType):
637                    param.appendChild(dataNode(document, 'Value', value))
638              param.appendChild(dataNode(document, 'Name', name))              else:
639                    param.appendChild(dataNode(document, 'Value', str(value)))
640              val = document.createElement('Value')  
641                node.appendChild(param)
642              if isinstance(value,ParameterSet):  
643                  value.toDom(document, val)    
644              elif isinstance(value, Link):      def fromDom(cls, doc):
645                  value.toDom(document, val)  
646                  param.appendChild(val)          # Define a host of helper functions to assist us.
647              elif isinstance(value,StringType):          def _children(node):
648                  param.appendChild(dataNode(document, 'Value', value))              """
649              else:              Remove the empty nodes from the children of this node
650                  param.appendChild(dataNode(document, 'Value', str(value)))              """
651                return [x for x in node.childNodes
652              node.appendChild(param)                      if not isinstance(x, minidom.Text) or x.nodeValue.strip()]
653        
654      def writeXML(self,ostream=stdout):          def _floatfromValue(doc):
655          """writes the object as an XML object into an output stream"""              return float(doc.nodeValue.strip())
656          # ParameterSet(d) with d[Name]=Value  
657          document, node = esysDoc()          def _stringfromValue(doc):
658          self.toDom(document, node)              return str(doc.nodeValue.strip())
659          ostream.write(document.toprettyxml())        
660            def _intfromValue(doc):
661                return int(doc.nodeValue.strip())
662    
663            def _boolfromValue(doc):
664                return bool(doc.nodeValue.strip())
665          
666            # Mapping from text types in the xml to methods used to process trees of that type
667            ptypemap = {"Simulation": Simulation.fromDom,
668                        "Model":Model.fromDom,
669                        "ParameterSet":ParameterSet.fromDom,
670                        "Link":Link.fromDom,
671                        "float":_floatfromValue,
672                        "int":_intfromValue,
673                        "str":_stringfromValue,
674                        "bool":_boolfromValue
675                        }
676    
677            parameters = {}
678            for node in _children(doc):
679                ptype = node.getAttribute("type")
680    
681                pname = pvalue = None
682                for childnode in _children(node):
683    
684                    if childnode.tagName == "Name":
685                        pname = childnode.firstChild.nodeValue.strip()
686    
687                    if childnode.tagName == "Value":
688                        nodes = _children(childnode)
689                        pvalue = ptypemap[ptype](nodes[0])
690    
691                parameters[pname] = pvalue
692    
693            # Create the instance of ParameterSet
694            o = cls()
695            o.declareParameters(parameters)
696            registerLinkableObject(doc.getAttribute("id"), o)
697            return o
698        
699        fromDom = classmethod(fromDom)
700        
701        def writeXML(self,ostream=stdout):
702            """writes the object as an XML object into an output stream"""
703            # ParameterSet(d) with d[Name]=Value
704            document, node = esysDoc()
705            self.toDom(document, node)
706            ostream.write(document.toprettyxml())
707    
708  class Model(ParameterSet,SimulationFrame):  class Model(ParameterSet,SimulationFrame):
709      """a Model is a SimulationFrame which is also a ParameterSet."""      """a Model is a SimulationFrame which is also a ParameterSet."""
710    
711        def __init__(self,parameters=[],**kwargs):
712            """creates a model"""
713            ParameterSet.__init__(self, parameters=parameters)
714            SimulationFrame.__init__(self,**kwargs)
715    
716        def toDom(self, document, node):
717            """ toDom method of Model class """
718            pset = document.createElement('Model')
719            pset.setAttribute('type', self.__class__.__name__)
720            node.appendChild(pset)
721            self._parametersToDom(document, pset)
722    
     def __init__(self,debug=False):  
         """creates a model"""  
         ParameterSet.__init__(self)  
         SimulationFrame.__init__(self,debug=debug)  
   
     def writeXML(self,ostream=stdout,esysheader=True,modelheader=True):  
         """writes the object as an XML object into an output stream"""  
         # x=type() with x.Name=Value where type is a Model from the model library.  
         if esysheader: ostream.write(_HEADER_START)      
         ostream.write("<Model type=\"%s\" module=\"%s\">\n"%(self.__class__.__name__,self.__class__.__module__))  
         self. writeParameterXML(ostream)        
         ostream.write("</Model>\n")  
         if esysheader: ostream.write(_HEADER_END)  
723    
724  class IterationDivergenceError(Exception):  class IterationDivergenceError(Exception):
725      """excpetion which should be thrown if there is no convergence of the iteration process at a time step but there is a chance taht a smaller step could help      """excpetion which should be thrown if there is no convergence of the iteration process at a time step but there is a chance taht a smaller step could help
726         to reach convergence."""         to reach convergence."""
727      pass      pass
728    
729  class IterationBreakDownError(Exception):  class IterationBreakDownError(Exception):
730      """excpetion which should be thrown if there is no conevregence and there is no chance that a time step reduction would help"""      """excpetion which should be thrown if there is no conevregence and there is no chance that a time step reduction would help"""
731      pass      pass
732    
733  class FailedTimeStepError(Exception):  class FailedTimeStepError(Exception):
734      """excpetion which should be thrown if the time step fails because of a step size that have been choosen to be to large"""      """excpetion which should be thrown if the time step fails because of a step size that have been choosen to be to large"""
735      pass      pass
736    
737  class NonPositiveStepSizeError(Exception):  class NonPositiveStepSizeError(Exception):
738      """excpetion which is thrown if the step size is not positive"""      """excpetion which is thrown if the step size is not positive"""
739      pass      pass
740    
741  #  #
742  #   ignore this text:  #   ignore this text:
743  #  #
744  """ the Model class provides a framework to run a time-dependent simulation. A  """
745    the Model class provides a framework to run a time-dependent simulation. A
746  Model has a set of parameter which may be fixed or altered by the Model itself  Model has a set of parameter which may be fixed or altered by the Model itself
747  or other Models over time.    or other Models over time.  
748    
749         The parameters of a models are declared at instantion, e.g.         The parameters of a models are declared at instantion, e.g.
750    
751             m=Model({"message" : "none" })             m=Model({"message" : "none" })
752    
753         creates a Model with parameters p1 and p2 with inital values 1 and 2.         creates a Model with parameters p1 and p2 with inital values 1 and 2.
754         Typically a particular model is defined as a subclass of Model:         Typically a particular model is defined as a subclass of Model:
755    
756          class Messenger(Model):          class Messenger(Model):
757              def __init__(self):              def __init__(self):
758                 Model.__init__(self,parameters={"message" : "none" })                 Model.__init__(self,parameters={"message" : "none" })
759    
760          m=MyModel()          m=MyModel()
761    
762         There are various ways how model parameters can be changed:         There are various ways how model parameters can be changed:
763    
764         1) use object attributes:         1) use object attributes:
765    
766            m.message="Hello World!"            m.message="Hello World!"
767    
768         2) use setParamter method         2) use setParamter method
769                    
770                        
771            m.setParameters(message="Hello World!")            m.setParameters(message="Hello World!")
772    
773         3) or dictonaries         3) or dictonaries
774        
775             d={ message : "Hello World!" }             d={ message : "Hello World!" }
776             m.setParameters(**d)             m.setParameters(**d)
777    
778    
779         A model executed buy staring the run method of the model:         A model executed buy staring the run method of the model:
780    
781            m=Messenger()            m=Messenger()
782            m.run()            m.run()
783    
784         The run methods marches through time. It first calls the         The run methods marches through time. It first calls the
785         doInitialization() method of the Model to set up the process. In each         doInitialization() method of the Model to set up the process. In each
786         time step the doStep() method is called to get from the current to the         time step the doStep() method is called to get from the current to the
787         next time step. The step size is defined by calling the         next time step. The step size is defined by calling the
788         getSafeTimeStepSize() method.  the time integration process is         getSafeTimeStepSize() method.  the time integration process is
789         terminated when the finalize() methods return true. Final the         terminated when the finalize() methods return true. Final the
790         doFinalization() method is called to finalize the process. To implement         doFinalization() method is called to finalize the process. To implement
791         a particular model a subclass of the Model class is defined. The         a particular model a subclass of the Model class is defined. The
792         subclass overwrites the default methods of Model.         subclass overwrites the default methods of Model.
793    
794         The following class defines a messenger printing in the doStep method         The following class defines a messenger printing in the doStep method
795         what ever the current value of its parameter message is:         what ever the current value of its parameter message is:
796    
797         class Messenger(Model):         class Messenger(Model):
798              def __init__(self):              def __init__(self):
799                 Model.__init__(self,parameters={"message" : "none" })                 Model.__init__(self,parameters={"message" : "none" })
800    
801              def doInitialization(self):              def doInitialization(self):
802                 print "I start talking now!"                 print "I start talking now!"
803    
804              def doStep(self,t):              def doStep(self,t):
805                 print "Message (time %e) : %s "%(t,self.message)                 print "Message (time %e) : %s "%(t,self.message)
806    
807              def doFinalization(self):              def doFinalization(self):
808                 print "I have no more to say!"                 print "I have no more to say!"
809                
810         If a instance of the Messenger class is run, it will print the         If a instance of the Messenger class is run, it will print the
811         initialization and finalization message only.  This is because the         initialization and finalization message only.  This is because the
812         default method for finalize() does always returns True and therefore the         default method for finalize() does always returns True and therefore the
813         transition is terminated startcht away.         transition is terminated startcht away.
814                    
815         Following example for solving the ODE using a forward euler scheme:         Following example for solving the ODE using a forward euler scheme:
816    
817                  u(t=0)=u0                  u(t=0)=u0
818                  u_t=a*u**2       for all 0<t<=ten                  u_t=a*u**2       for all 0<t<=ten
819    
820        exact solution is given by u(t)=1/(1/u0-a*t)        exact solution is given by u(t)=1/(1/u0-a*t)
821    
822        class  Ode1(Model):        class  Ode1(Model):
823           def __init__(self,**args):           def __init__(self,**args):
824              Model.__init__(self,parameters={"tend" : 1., "dt" : 0.0001 ,"a" : 0.1 ,"u" : 1. },name="test",debug=True)              Model.__init__(self,parameters={"tend" : 1., "dt" : 0.0001 ,"a" : 0.1 ,"u" : 1. },name="test",debug=True)
825    
826           def doInitialization(self):           def doInitialization(self):
827               self._tn=0               self._tn=0
828    
829           def doStep(self,t):           def doStep(self,t):
830               self.u=self.u+(t-self._tn)*self.a*self.u**2               self.u=self.u+(t-self._tn)*self.a*self.u**2
831               self._tn=t               self._tn=t
832    
833           def doFinalization(self):           def doFinalization(self):
834               print "all done final error = ",abs(self.u-1./(1./3.-self.a*self._tn))               print "all done final error = ",abs(self.u-1./(1./3.-self.a*self._tn))
835    
836           def getSafeTimeStepSize(self):           def getSafeTimeStepSize(self):
837               return self.dt               return self.dt
838    
839           def finalize(self):           def finalize(self):
840               return self._tn>=self.tend               return self._tn>=self.tend
841    
842         In some cases at a given time step an iteration process has to be         In some cases at a given time step an iteration process has to be
843         performed to get the state of the Model for the next time step. ` In         performed to get the state of the Model for the next time step. ` In
844         this case the doStep() method is replaced by a sequance of methods which         this case the doStep() method is replaced by a sequance of methods which
845         implements this iterative process.  The method then will control the         implements this iterative process.  The method then will control the
846         iteration process by initializing the iteration through calling the         iteration process by initializing the iteration through calling the
847         doIterationInitialization() method. The iteration is preformed by         doIterationInitialization() method. The iteration is preformed by
848         calling the doIterationStep() method until the terminate() method         calling the doIterationStep() method until the terminate() method
849         returns True. The doIterationFinalization() method is called to end the         returns True. The doIterationFinalization() method is called to end the
850         iteration.         iteration.
851         For a particular model these methods have to overwritten by a suitable         For a particular model these methods have to overwritten by a suitable
852         subclass without touching the doStep() method.         subclass without touching the doStep() method.
853    
854         following example is a modification of the example above. Here an         following example is a modification of the example above. Here an
855         implicit euler scheme is used. in each time step the problem         implicit euler scheme is used. in each time step the problem
856                        
857             0= u_{n+1}-u_{n}+a*dt*u_{n+1}**2             0= u_{n+1}-u_{n}+a*dt*u_{n+1}**2
858    
859         has to be solved for u_{n+1}. The Newton scheme is used to solve this non-linear problem.         has to be solved for u_{n+1}. The Newton scheme is used to solve this non-linear problem.
860    
861    
862        class  Ode2(Model):        class  Ode2(Model):
863    
864         def __init__(self,**args):         def __init__(self,**args):
865             Model.__init__(self,{"tend" : 1., "dt" : 0.1 ,"a" : 10. ,"u" : 1. , "tol " : 1.e-8},"test","bla",None,True)             Model.__init__(self,{"tend" : 1., "dt" : 0.1 ,"a" : 10. ,"u" : 1. , "tol " : 1.e-8},"test","bla",None,True)
866    
867         def doInitialization(self):         def doInitialization(self):
868             self.__tn=0             self.__tn=0
869    
870         def doIterationInitialization(self,t):         def doIterationInitialization(self,t):
871              self.__iter=0              self.__iter=0
872              self.u_last=self.u                          self.u_last=self.u            
873              self.current_dt=t-self.tn              self.current_dt=t-self.tn
874              self.__tn=t              self.__tn=t
875    
876         def doIterationStep(self):         def doIterationStep(self):
877            self.__iter+=1            self.__iter+=1
878            self.u_old=self.u            self.u_old=self.u
879            self.u=(self.current_dt*self.a*self.u**2-self.u_last)/(2*self.current_dt*self.a*self.u-1.)            self.u=(self.current_dt*self.a*self.u**2-self.u_last)/(2*self.current_dt*self.a*self.u-1.)
880    
881         def terminate(self):         def terminate(self):
882             return abs(self.u_old-self.u)<self.tol*abs(self.u)             return abs(self.u_old-self.u)<self.tol*abs(self.u)
883    
884         def doIterationFinalization(self)         def doIterationFinalization(self)
885             print "all done"             print "all done"
886    
887         def getSafeTimeStepSize(self):         def getSafeTimeStepSize(self):
888             return self.dt             return self.dt
889    
890         def finalize(self):         def finalize(self):
891              return self.__tn>self.tend              return self.__tn>self.tend
892    
893         A model can be composed from subsimulations. Subsimulations are treated         A model can be composed from subsimulations. Subsimulations are treated
894         as model parameters. If a model parameter is set or a value of a model         as model parameters. If a model parameter is set or a value of a model
895         parameter is requested, the model will search for this parameter its         parameter is requested, the model will search for this parameter its
896         subsimulations in the case the model does not have this parameter         subsimulations in the case the model does not have this parameter
897         itself. The order in which the subsimulations are searched is critical.         itself. The order in which the subsimulations are searched is critical.
898         By default a Model initializes all its subsimulations, is finalized when         By default a Model initializes all its subsimulations, is finalized when
899         all its subsimulations are finalized and finalizes all its         all its subsimulations are finalized and finalizes all its
900         subsimulations. In the case an iterative process is applied on a         subsimulations. In the case an iterative process is applied on a
901         particular time step the iteration is initialized for all         particular time step the iteration is initialized for all
902         subsimulations, then the iteration step is performed for each         subsimulations, then the iteration step is performed for each
903         subsimulation until all subsimulations indicate termination. Then the         subsimulation until all subsimulations indicate termination. Then the
904         iteration is finalized for all subsimulations. Finally teh doStop()         iteration is finalized for all subsimulations. Finally teh doStop()
905         method for all submethods is called.         method for all submethods is called.
906    
907         Here we are creating a model which groups ab instantiation of the Ode2 and the Messenger Model         Here we are creating a model which groups ab instantiation of the Ode2 and the Messenger Model
908    
909         o=Ode2()         o=Ode2()
910         m=Messenger()         m=Messenger()
911         om=Model(subsimulations=[o,m],debug=True)         om=Model(subsimulations=[o,m],debug=True)
912         om.dt=0.01         om.dt=0.01
913         om.u=1.         om.u=1.
914         m.message="it's me!"         m.message="it's me!"
915         om.run()         om.run()
916    
917         Notice that dt and u are parameters of class Ode2 and message is a         Notice that dt and u are parameters of class Ode2 and message is a
918         parameter of the Messenger class. The Model formed from these models         parameter of the Messenger class. The Model formed from these models
919         automatically hand the assignment of new values down to the         automatically hand the assignment of new values down to the
920         subsimulation. om.run() starts this combined model where now the         subsimulation. om.run() starts this combined model where now the
921         soStep() method of the Messenger object printing the value of its         soStep() method of the Messenger object printing the value of its
922         parameter message together with a time stamp is executed in each time         parameter message together with a time stamp is executed in each time
923         step introduced by the Ode2 model.         step introduced by the Ode2 model.
924    
925         A parameter of a Model can be linked to an attribute of onother object,         A parameter of a Model can be linked to an attribute of onother object,
926         typically an parameter of another Model object.             typically an parameter of another Model object.    
927                
928                
929         which is comprised by a set of subsimulations.         which is comprised by a set of subsimulations.
930         The simulation is run through its run method which in the simplest case has the form:         The simulation is run through its run method which in the simplest case has the form:
931    
932            s=Model()            s=Model()
933            s.run()            s.run()
934    
935         The run has an initializion and finalization phase. The latter is called         The run has an initializion and finalization phase. The latter is called
936         if all subsimulations are to be finalized. The simulation is processing         if all subsimulations are to be finalized. The simulation is processing
937         in time through calling the stepForward methods which updates the         in time through calling the stepForward methods which updates the
938         observables of each subsimulation.  A time steps size which is save for         observables of each subsimulation.  A time steps size which is save for
939         all subsimulation is choosen.         all subsimulation is choosen.
940    
941         At given time step an iterative process may be performed to make sure         At given time step an iterative process may be performed to make sure
942         that all observables are consistent across all subsimulations.  In this         that all observables are consistent across all subsimulations.  In this
943         case, similar the time dependence, an initialization and finalization of         case, similar the time dependence, an initialization and finalization of
944         the iteration is performed.         the iteration is performed.
945    
946         A Model has input and output parameters where each input parameter can         A Model has input and output parameters where each input parameter can
947         be constant, time dependent or may depend on an output parameter of         be constant, time dependent or may depend on an output parameter of
948         another model or the model itself. To create a parameter name of a model         another model or the model itself. To create a parameter name of a model
949         and to assign a value to it one can use the statement         and to assign a value to it one can use the statement
950    
951             model.name=object             model.name=object
952    
953    
954         At any time the current value of the parameter name can be obtained by         At any time the current value of the parameter name can be obtained by
955    
956                 value=model.name                 value=model.name
957    
958         If the object that has been assigned to the paramter/attribute name has         If the object that has been assigned to the paramter/attribute name has
959         the attribute/parameter name isself the current value of this attribute         the attribute/parameter name isself the current value of this attribute
960         of the object is returned (e.g. for model.name=object where object has         of the object is returned (e.g. for model.name=object where object has
961         an attribute name, the statement value=model.name whould assign the         an attribute name, the statement value=model.name whould assign the
962         value object.name to value.). If the name of the parameters of a model         value object.name to value.). If the name of the parameters of a model
963         and an object don't match the setParameter method of model can be used.         and an object don't match the setParameter method of model can be used.
964         So         So
965    
966             model.setParameter(name,object,name_for_object)             model.setParameter(name,object,name_for_object)
967    
968         links the parameter name of model with the parameter name_for_object of         links the parameter name of model with the parameter name_for_object of
969         object.         object.
970    
971         The run method initiates checkpointing (it is not clear how to do this         The run method initiates checkpointing (it is not clear how to do this
972         yet)         yet)
973     =====     =====
974                            
975     """     """
976    
977    
978    
979  if __name__=="__main__":  if __name__=="__main__":
980      import math      import math
981      #      #
982      # test for parameter set      # test for parameter set
983      #      #
984      p11=ParameterSet()      p11=ParameterSet()
985      p11.declareParameter(gamma1=1.,gamma2=2.,gamma3=3.)      p11.declareParameter(gamma1=1.,gamma2=2.,gamma3=3.)
986      p1=ParameterSet()      p1=ParameterSet()
987      p1.declareParameter(dim=2,tol_v=0.001,output_file="/tmp/u.%3.3d.dx",runFlag=True,parm11=p11)      p1.declareParameter(dim=2,tol_v=0.001,output_file="/tmp/u.%3.3d.dx",runFlag=True,parm11=p11)
988      parm=ParameterSet({ "parm1" : p1 , "parm2" : ParameterSet(["alpha"])})      parm=ParameterSet({ "parm1" : p1 , "parm2" : ParameterSet(["alpha"])})
989      parm.parm2.alpha=Link(p11,"gamma1")      parm.parm2.alpha=Link(p11,"gamma1")
990      parm.x="that should not be here!"      parm.x="that should not be here!"
991      print parm.showParameters()      print parm.showParameters()
992      # should be something like: {"parm2" : {"alpha" : reference to attribute      # should be something like: {"parm2" : {"alpha" : reference to attribute
993      # gamma1 of <__main__.ParameterSet instance at 0xf6db51cc>},"parm1" : {"dim"      # gamma1 of <__main__.ParameterSet instance at 0xf6db51cc>},"parm1" : {"dim"
994      # : 2,"runFlag" : True,"tol_v": 0.001,"parm11" : {"gamma3" : 3.0,"gamma2" :      # : 2,"runFlag" : True,"tol_v": 0.001,"parm11" : {"gamma3" : 3.0,"gamma2" :
995      # 2.0,"gamma1" : 1.0},"output_file" : /tmp/u.%3.3d.dx}}      # 2.0,"gamma1" : 1.0},"output_file" : /tmp/u.%3.3d.dx}}
996      assert parm.parm2.alpha==1.      assert parm.parm2.alpha==1.
997      parm.writeXML()      parm.writeXML()
998            
999      #=======================      #=======================
1000      class Messenger(Model):      class Messenger(Model):
1001          def __init__(self):          def __init__(self):
1002              Model.__init__(self)              Model.__init__(self)
1003              self.declareParameter(message="none")              self.declareParameter(message="none")
1004    
1005          def doInitialization(self,t):          def doInitialization(self,t):
1006              self.__t=t              self.__t=t
1007              print "I start talking now!"              print "I start talking now!"
1008    
1009          def doStep(self,dt):          def doStep(self,dt):
1010              self.__t+=dt              self.__t+=dt
1011              print "Message (time %e) : %s "%(self.__t,self.message)              print "Message (time %e) : %s "%(self.__t,self.message)
1012    
1013          def doFinalization(self):          def doFinalization(self):
1014              print "I have no more to say!"              print "I have no more to say!"
1015        
1016      class ODETEST(Model):      class ODETEST(Model):
1017          """ implements a solver for the ODE          """ implements a solver for the ODE
1018    
1019                du/dt=a*u+f(t)                du/dt=a*u+f(t)
1020    
1021             we use a implicit euler scheme :             we use a implicit euler scheme :
1022    
1023                u_n-u_{n-1}= dt*a u_n + st*f(t_n)                u_n-u_{n-1}= dt*a u_n + st*f(t_n)
1024    
1025             to get u_n we run an iterative process             to get u_n we run an iterative process
1026    
1027                 u_{n.k}=u_{n-1}+dt*(a u_{n.i-1} + f(t_n))                 u_{n.k}=u_{n-1}+dt*(a u_{n.i-1} + f(t_n))
1028    
1029    
1030             input for this model are step size dt, end time tend and a value for             input for this model are step size dt, end time tend and a value for
1031             a, f and initial value for u. we need also a tolerance tol for a             a, f and initial value for u. we need also a tolerance tol for a
1032             stopping criterion.             stopping criterion.
1033    
1034          """          """
1035    
1036          def __init__(self):          def __init__(self):
1037              Model.__init__(self,debug=True)              Model.__init__(self,debug=True)
1038              self.declareParameter(tend=1.,dt=0.1,a=0.9,u=10.,f=0.,message="",tol=1.e-8)              self.declareParameter(tend=1.,dt=0.1,a=0.9,u=10.,f=0.,message="",tol=1.e-8)
1039    
1040          def doInitialization(self,t):          def doInitialization(self,t):
1041              self.__tn=t              self.__tn=t
1042              self.__iter=0              self.__iter=0
1043                        
1044          def doIterationInitialization(self,dt):          def doIterationInitialization(self,dt):
1045              self.__iter=0              self.__iter=0
1046              self.__u_last=self.u                          self.__u_last=self.u            
1047    
1048          def doIterationStep(self,dt):          def doIterationStep(self,dt):
1049              self.__iter+=1              self.__iter+=1
1050              self.__u_old=self.u              self.__u_old=self.u
1051              self.u=self.__u_last+dt*(self.a*self.__u_old+self.f)              self.u=self.__u_last+dt*(self.a*self.__u_old+self.f)
1052                    
1053          def terminate(self):          def terminate(self):
1054              if self.__iter<1:              if self.__iter<1:
1055                  return False                  return False
1056              else:              else:
1057                  return abs(self.__u_old-self.u)<self.tol*abs(self.u)                  return abs(self.__u_old-self.u)<self.tol*abs(self.u)
1058    
1059          def doIterationFinalization(self,dt):          def doIterationFinalization(self,dt):
1060              self.__tn+=dt              self.__tn+=dt
1061              self.message="current error = %e"%abs(self.u-10.*math.exp((self.a-1.)*self.__tn))              self.message="current error = %e"%abs(self.u-10.*math.exp((self.a-1.)*self.__tn))
1062    
1063          def getSafeTimeStepSize(self,dt):          def getSafeTimeStepSize(self,dt):
1064              return min(self.dt,1./(abs(self.a)+1.))              return min(self.dt,1./(abs(self.a)+1.))
1065    
1066          def finalize(self):          def finalize(self):
1067              return self.__tn>=self.tend              return self.__tn>=self.tend
1068    
1069      #      #
1070      #   s solves the coupled ODE:      #   s solves the coupled ODE:
1071      #      #
1072      #     du/dt=a*u+  v      #     du/dt=a*u+  v
1073      #     dv/dt=  u+a*v      #     dv/dt=  u+a*v
1074      #        #  
1075      #   each equation is treated through the ODETEST class. The equations are      #   each equation is treated through the ODETEST class. The equations are
1076      #   linked and iteration over each time step is performed.  the current      #   linked and iteration over each time step is performed.  the current
1077      #   error of v is reported by the Messenger class.      #   error of v is reported by the Messenger class.
1078      #      #
1079      o1=ODETEST()      o1=ODETEST()
1080      o1.u=10      o1.u=10
1081      o2=ODETEST()      o2=ODETEST()
1082      o2.u=-10.      o2.u=-10.
1083      o1.f=Link(o2,"u")      o1.f=Link(o2,"u")
1084      o2.f=Link(o1,"u")      o2.f=Link(o1,"u")
1085      m=Messenger()      m=Messenger()
1086      o1.dt=0.01      o1.dt=0.01
1087      m.message=Link(o1)      m.message=Link(o1)
1088      s=ExplicitSimulation([Simulation([o1,o2],debug=True),m],debug=True)      s=ExplicitSimulation([Simulation([o1,o2],debug=True),m],debug=True)
1089      s.run()      s.run()
1090      s.writeXML()      s.writeXML()

Legend:
Removed from v.122  
changed lines
  Added in v.126

  ViewVC Help
Powered by ViewVC 1.1.26