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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 142 - (hide annotations)
Mon Jul 25 05:28:20 2005 UTC (13 years, 8 months ago) by jgs
Original Path: trunk/esys2/escript/py_src/modelframe.py
File MIME type: text/x-python
File size: 38069 byte(s)
Merge of development branch back to main trunk on 2005-07-25

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

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26