/[escript]/trunk/modellib/py_src/input.py
ViewVC logotype

Annotation of /trunk/modellib/py_src/input.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1098 - (hide annotations)
Mon Apr 16 23:15:23 2007 UTC (14 years, 7 months ago) by gross
File MIME type: text/x-python
File size: 20010 byte(s)
add a few #pragma ivdep which should speed up MV but cannot confiirm this on Pentium




1 jgs 127 # $Id$
2 elspeth 628
3     __copyright__=""" Copyright (c) 2006 by ACcESS MNRF
4     http://www.access.edu.au
5     Primary Business: Queensland, Australia"""
6     __license__="""Licensed under the Open Software License version 3.0
7     http://www.opensource.org/licenses/osl-3.0.php"""
8    
9 jgs 149 from esys.escript import *
10     from esys.escript.modelframe import Model,ParameterSet
11 gross 927 from esys.escript.linearPDEs import LinearPDE
12 jgs 147 from math import log
13 jgs 127
14 jgs 147 class Sequencer(Model):
15 jgs 148 """
16 gross 814 Runs through time until t_end is reached.
17 gross 720
18 gross 814 @ivar t_end: model is terminated when t_end is passed, default 1 (in).
19     @type t_end: C{float}
20     @ivar dt_max: maximum time step size, default L{Model.UNDEF_DT} (in)
21     @type dt_max: C{float}
22     @ivar t: current time stamp (in/out). By default it is initialized with zero.
23     @type t: C{float}
24    
25 jgs 148 """
26 gross 918 def __init__(self,**kwargs):
27 jgs 148 """
28 jgs 149 """
29 gross 918 super(Sequencer,self).__init__(**kwargs)
30 gross 423 self.declareParameter(t=0.,
31     t_end=1.,
32     dt_max=Model.UNDEF_DT)
33 jgs 127
34 jgs 147 def doInitialization(self):
35 jgs 148 """
36 jgs 149 initialize time integration
37 jgs 148 """
38     self.__t_old = self.t
39 jgs 147
40 jgs 148 def doStepPreprocessing(self, dt):
41     self.t = self.__t_old+dt
42 jgs 147
43 jgs 148 def doStepPostprocessing(self, dt):
44     self.__t_old = self.t
45 jgs 147
46     def finalize(self):
47 jgs 148 """
48 gross 814 returns true when L{t} has reached L{t_end}
49 jgs 148 """
50     return self.t >= self.t_end
51 jgs 147
52 jgs 148 def getSafeTimeStepSize(self, dt):
53     """
54 gross 814 returns L{dt_max}
55 jgs 148 """
56     return self.dt_max
57 jgs 147
58 jgs 148 class GaussianProfile(ParameterSet):
59     """
60 jgs 149 Generates a Gaussian profile at center x_c, width width and height A
61 jgs 148 over a domain
62 jgs 127
63 gross 814 @ivar domain: domain
64     @ivar x_c: center of the Gaussian profile (default [0.,0.,0.])
65     @ivar A: (in) height of the profile. A maybe a vector. (default 1.)
66     @ivar width: (in) width of the profile (default 0.1)
67     @ivar r: (in) radius of the circle (default = 0)
68 jgs 127
69 jgs 148 In the case that the spatial dimension is two, The third component of
70 gross 814 x_c is dropped.
71 jgs 127 """
72 gross 918 def __init__(self,**kwargs):
73     super(GaussianProfile, self).__init__(**kwargs)
74 jgs 147 self.declareParameter(domain=None,
75     x_c=numarray.zeros([3]),
76     A=1.,
77     width=0.1,
78     r=0)
79 jgs 127
80     def out(self):
81 jgs 149 """
82     Generate the Gaussian profile
83 gross 814
84     Link against this method to get the output of this model.
85 jgs 149 """
86 jgs 148 x = self.domain.getX()
87     dim = self.domain.getDim()
88     l = length(x-self.x_c[:dim])
89 gross 323 m = whereNegative(l-self.r)
90 jgs 148
91 jgs 127 return (m+(1.-m)*exp(-log(2.)*(l/self.width)**2))*self.A
92    
93 jgs 148 class InterpolateOverBox(ParameterSet):
94 jgs 147 """
95 jgs 148 Returns values at each time. The values are defined through given values
96 gross 814 at time node. For two dimensional domains back values are ignored.
97 jgs 147
98 gross 814 @ivar domain: domain
99     @ivar value_left_bottom_front: (in) value at left,bottom,front corner
100     @ivar value_right_bottom_front: (in) value at right, bottom, front corner
101     @ivar value_left_top_front: (in) value at left,top,front corner
102     @ivar value_right_top_front: (in) value at right,top,front corner
103     @ivar value_left_bottom_back: (in) value at left,bottom,back corner
104     @ivar value_right_bottom_back: (in) value at right,bottom,back corner
105     @ivar value_left_top_back: (in) value at left,top,back corner
106     @ivar value_right_top_back: (in) value at right,top,back corner
107 jgs 147 """
108    
109 gross 918 def __init__(self, **kwargs):
110 gross 911 super(InterpolateOverBox, self).__init__(self)
111 jgs 147 self.declareParameter(domain=None,
112     value_left_bottom_front=0.,
113     value_right_bottom_front=0.,
114     value_left_top_front=0.,
115     value_right_top_front=0.,
116     value_left_bottom_back=0.,
117     value_right_bottom_back=0.,
118     value_left_top_back=0.,
119     value_right_top_back=0.)
120    
121    
122     def out(self):
123 gross 814 """
124     values at domain locations by bilinear interpolation of the given values.
125    
126     Link against this method to get the output of this model.
127     """
128 jgs 148 x = self.domain.getX()
129     if self.domain.getDim() == 2:
130 gross 819 x0,x1=x[0],x[1]
131     left_bottom_front0,right_top_back0=inf(x0),sup(x0)
132 gross 820 left_bottom_front1,right_top_back1=inf(x1),sup(x1)
133     f_right = (x0 - left_bottom_front0)/(right_top_back0 -left_bottom_front0)
134 jgs 148 f_left = 1. - f_right
135 gross 820 f_top = (x1 - left_bottom_front1)/(right_top_back1 - left_bottom_front1)
136 jgs 148 f_bottom = 1. - f_top
137 gross 816 out = f_left * f_bottom * self.value_left_bottom_front \
138     + f_right * f_bottom * self.value_right_bottom_front \
139     + f_left * f_top * self.value_left_top_front \
140     + f_right * f_top * self.value_right_top_front
141 jgs 147 else:
142 gross 819 x0,x1,x2=x[0],x[1],x[2]
143     left_bottom_front0,right_top_back0=inf(x0),sup(x0)
144 gross 820 left_bottom_front1,right_top_back1=inf(x1),sup(x1)
145     left_bottom_front2,right_top_back2=inf(x2),sup(x2)
146     f_right = (x0 - left_bottom_front0)/(right_top_back0 - left_bottom_front0)
147 jgs 148 f_left = 1. - f_right
148 gross 820 f_top = (x1 - left_bottom_front1)/(right_top_back1 - left_bottom_front1)
149 jgs 148 f_bottom = 1. - f_top
150 gross 820 f_back = (x2 - left_bottom_front1)/(right_top_back2 - left_bottom_front2)
151 jgs 148 f_front = 1. - f_back
152 gross 816 out = f_left * f_bottom * f_front * self.value_left_bottom_front\
153     + f_right * f_bottom * f_front * self.value_right_bottom_front\
154     + f_left * f_top * f_front * self.value_left_top_front\
155     + f_right * f_top * f_front * self.value_right_top_front\
156     + f_left * f_bottom * f_back * self.value_left_bottom_back\
157     + f_right * f_bottom * f_back * self.value_right_bottom_back\
158     + f_left * f_top * f_back * self.value_left_top_back\
159     + f_right * f_top * f_back * self.value_right_top_back
160 jgs 147 return out
161    
162    
163 jgs 148 class InterpolatedTimeProfile(ParameterSet):
164     """
165    
166     Returns values at each time. The values are defined through given
167     values at time node.
168 jgs 147
169 jgs 148 value[i] defines the value at time nodes[i]. Between nodes linear
170     interpolation is used.
171    
172     For time t<nodes[0], value[0] is used and for t>nodes[l], values[l]
173     is used where l=len(nodes)-1.
174 jgs 147
175 gross 814 @ivar t: (in) current time
176     @ivar node: (in) list of time nodes
177     @ivar values: (in) list of values at time nodes
178 jgs 148 """
179 jgs 127
180 gross 918 def __init__(self,**kwargs):
181     super( InterpolatedTimeProfile, self).__init__(**kwargs)
182 jgs 147 self.declareParameter(t=0., \
183     nodes=[0.,1.],\
184     values=[1.,1.])
185     def out(self):
186 gross 814 """
187     current value
188    
189     Link against this method to get the output of this model.
190     """
191 jgs 148 l = len(self.nodes) - 1
192     t = self.t
193     if t <= self.nodes[0]:
194     return self.values[0]
195 jgs 147 else:
196     for i in range(1,l):
197 jgs 148 if t < self.nodes[i]:
198     m = (self.values[i-1] - self.values[i])/\
199     (self.nodes[i-1] - self.nodes[i])
200     return m*(t-self.nodes[i-1]) + self.values[i-1]
201 jgs 147 return self.values[l]
202 jgs 127
203 gross 903 class ScalarDistributionFromTags(ParameterSet):
204 jgs 148 """
205 gross 953 creates a scalar distribution on a domain from tags, If tag_map is given
206     the tags can be given a names and tag_map is used to map it into domain tags.
207 gross 903
208     @ivar domain: domain
209     @type domain: L{esys.escript.Domain}
210     @ivar default: default value
211     @ivar tag0: tag 0
212     @type tag0: C{int}
213     @ivar value0: value for tag 0
214     @type value0: C{float}
215     @ivar tag1: tag 1
216     @type tag1: C{int}
217     @ivar value1: value for tag 1
218     @type value1: C{float}
219     @ivar tag2: tag 2
220     @type tag2: C{int}
221     @ivar value2: value for tag 2
222     @type value2: C{float}
223     @ivar tag3: tag 3
224     @type tag3: C{int}
225     @ivar value3: value for tag 3
226     @type value3: C{float}
227     @ivar tag4: tag 4
228     @type tag4: C{int}
229     @ivar value4: value for tag 4
230     @type value4: C{float}
231     @ivar tag5: tag 5
232     @type tag5: C{int}
233     @ivar value5: value for tag 5
234     @type value5: C{float}
235     @ivar tag6: tag 6
236     @type tag6: C{int}
237     @ivar value6: value for tag 6
238     @type value6: C{float}
239     @ivar tag7: tag 7
240     @type tag7: C{int}
241     @ivar value7: value for tag 7
242     @type value7: C{float}
243     @ivar tag8: tag 8
244     @type tag8: C{int}
245     @ivar value8: value for tag 8
246     @type value8: C{float}
247     @ivar tag9: tag 9
248     @type tag9: C{int}
249     @ivar value9: value for tag 9
250     @type value9: C{float}
251     """
252 gross 918 def __init__(self,**kwargs):
253     super(ScalarDistributionFromTags, self).__init__(**kwargs)
254 gross 903 self.declareParameter(domain=None,
255     default=0.,
256     tag0=None,
257     value0=0.,
258     tag1=None,
259     value1=0.,
260     tag2=None,
261     value2=0.,
262     tag3=None,
263     value3=0.,
264     tag4=None,
265     value4=0.,
266     tag5=None,
267     value5=0.,
268     tag6=None,
269     value6=0.,
270     tag7=None,
271     value7=0.,
272     tag8=None,
273     value8=0.,
274     tag9=None,
275     value9=0.)
276    
277    
278     def out(self):
279     """
280     returns a L{esys.escript.Data} object
281     Link against this method to get the output of this model.
282     """
283     d=Scalar(self.default,Function(self.domain))
284 gross 1098 if not self.tag0 == None: d.setTaggedValue(self.tag0,self.value0)
285     if not self.tag1 == None: d.setTaggedValue(self.tag1,self.value1)
286     if not self.tag2 == None: d.setTaggedValue(self.tag2,self.value2)
287     if not self.tag3 == None: d.setTaggedValue(self.tag3,self.value3)
288     if not self.tag4 == None: d.setTaggedValue(self.tag4,self.value4)
289     if not self.tag5 == None: d.setTaggedValue(self.tag5,self.value5)
290     if not self.tag6 == None: d.setTaggedValue(self.tag6,self.value6)
291     if not self.tag7 == None: d.setTaggedValue(self.tag7,self.value7)
292     if not self.tag8 == None: d.setTaggedValue(self.tag8,self.value8)
293     if not self.tag9 == None: d.setTaggedValue(self.tag9,self.value9)
294 gross 903 return d
295    
296 gross 927 class SmoothScalarDistributionFromTags(ParameterSet):
297     """
298     creates a smooth scalar distribution on a domain from region tags
299    
300     @ivar domain: domain
301     @type domain: L{esys.escript.Domain}
302     @ivar default: default value
303     @ivar tag0: tag 0
304     @type tag0: C{int}
305     @ivar value0: value for tag 0
306     @type value0: C{float}
307     @ivar tag1: tag 1
308     @type tag1: C{int}
309     @ivar value1: value for tag 1
310     @type value1: C{float}
311     @ivar tag2: tag 2
312     @type tag2: C{int}
313     @ivar value2: value for tag 2
314     @type value2: C{float}
315     @ivar tag3: tag 3
316     @type tag3: C{int}
317     @ivar value3: value for tag 3
318     @type value3: C{float}
319     @ivar tag4: tag 4
320     @type tag4: C{int}
321     @ivar value4: value for tag 4
322     @type value4: C{float}
323     @ivar tag5: tag 5
324     @type tag5: C{int}
325     @ivar value5: value for tag 5
326     @type value5: C{float}
327     @ivar tag6: tag 6
328     @type tag6: C{int}
329     @ivar value6: value for tag 6
330     @type value6: C{float}
331     @ivar tag7: tag 7
332     @type tag7: C{int}
333     @ivar value7: value for tag 7
334     @type value7: C{float}
335     @ivar tag8: tag 8
336     @type tag8: C{int}
337     @ivar value8: value for tag 8
338     @type value8: C{float}
339     @ivar tag9: tag 9
340     @type tag9: C{int}
341     @ivar value9: value for tag 9
342     @type value9: C{float}
343     """
344     def __init__(self,**kwargs):
345     super(SmoothScalarDistributionFromTags, self).__init__(**kwargs)
346     self.declareParameter(domain=None,
347     default=0.,
348     tag0=None,
349     value0=0.,
350     tag1=None,
351     value1=0.,
352     tag2=None,
353     value2=0.,
354     tag3=None,
355     value3=0.,
356     tag4=None,
357     value4=0.,
358     tag5=None,
359     value5=0.,
360     tag6=None,
361     value6=0.,
362     tag7=None,
363     value7=0.,
364     tag8=None,
365     value8=0.,
366     tag9=None,
367     value9=0.)
368    
369    
370     def __update(self,tag,tag_value,value):
371     if self.__pde==None:
372     self.__pde=LinearPDE(self.domain,numSolutions=1)
373     mask=Scalar(0.,Function(self.domain))
374 gross 1098 mask.setTaggedValue(tag,1.)
375 gross 927 self.__pde.setValue(Y=mask)
376     mask=wherePositive(abs(self.__pde.getRightHandSide()))
377     value*=(1.-mask)
378     value+=tag_value*mask
379     return value
380    
381     def out(self):
382     """
383     returns a L{esys.escript.Data} object
384     Link against this method to get the output of this model.
385     """
386     d=Scalar(self.default,Solution(self.domain))
387     self.__pde=None
388     if not self.tag0 == None: d=self.__update(self.tag0,self.value0,d)
389     if not self.tag1 == None: d=self.__update(self.tag1,self.value1,d)
390     if not self.tag2 == None: d=self.__update(self.tag2,self.value2,d)
391     if not self.tag3 == None: d=self.__update(self.tag3,self.value3,d)
392     if not self.tag4 == None: d=self.__update(self.tag4,self.value4,d)
393     if not self.tag5 == None: d=self.__update(self.tag5,self.value5,d)
394     if not self.tag6 == None: d=self.__update(self.tag6,self.value6,d)
395     if not self.tag7 == None: d=self.__update(self.tag7,self.value7,d)
396     if not self.tag8 == None: d=self.__update(self.tag8,self.value8,d)
397     if not self.tag9 == None: d=self.__update(self.tag9,self.value9,d)
398     return d
399    
400 gross 903 class LinearCombination(ParameterSet):
401     """
402 jgs 148 Returns a linear combination of the f0*v0+f1*v1+f2*v2+f3*v3+f4*v4
403 jgs 147
404 gross 819 @ivar f0: numerical object or None, default=None (in)
405     @ivar v0: numerical object or None, default=None (in)
406     @ivar f1: numerical object or None, default=None (in)
407     @ivar v1: numerical object or None, default=None (in)
408     @ivar f2: numerical object or None, default=None (in)
409     @ivar v2: numerical object or None, default=None (in)
410     @ivar f3: numerical object or None, default=None (in)
411     @ivar v3: numerical object or None, default=None (in)
412     @ivar f4: numerical object or None, default=None (in)
413     @ivar v4: numerical object or None, default=None (in)
414 jgs 148 """
415 gross 918 def __init__(self,**kwargs):
416     super(LinearCombination, self).__init__(**kwargs)
417 jgs 148 self.declareParameter(f0=None, \
418     v0=None, \
419     f1=None, \
420     v1=None, \
421     f2=None, \
422     v2=None, \
423     f3=None, \
424     v3=None, \
425     f4=None, \
426     v4=None)
427 jgs 147
428 jgs 148 def out(self):
429 gross 814 """
430     returns f0*v0+f1*v1+f2*v2+f3*v3+f4*v4.
431     Link against this method to get the output of this model.
432     """
433 jgs 148 if not self.f0 == None and not self.v0 == None:
434     fv0 = self.f0*self.v0
435     else:
436     fv0 = None
437 jgs 147
438 jgs 148 if not self.f1 == None and not self.v1 == None:
439     fv1 = self.f1*self.v1
440     else:
441     fv1 = None
442 jgs 147
443 jgs 148 if not self.f2 == None and not self.v2 == None:
444     fv2 = f2*v2
445     else:
446     fv2 = None
447 jgs 147
448 jgs 148 if not self.f3 == None and not self.v3 == None:
449     fv3 = self.f3*self.v3
450     else:
451     fv3 = None
452 jgs 147
453 jgs 148 if not self.f4 == None and not self.v4 == None:
454     fv4 = self.f4*self.v4
455     else:
456     fv4 = None
457 jgs 147
458 jgs 148 if fv0 == None:
459     out = 0.
460     else:
461     out = fv0
462     if not fv1 == None:
463     out += fv1
464     if not fv2 == None:
465     out += fv2
466     if not fv3 == None:
467     out += fv3
468     return out
469 jgs 147
470 gross 904 class MergeConstraints(ParameterSet):
471     """
472     Returns a linear combination of the f0*v0+f1*v1+f2*v2+f3*v3+f4*v4
473     """
474 gross 918 def __init__(self,**kwargs):
475     super(MergeConstraints, self).__init__(**kwargs)
476 gross 904 self.declareParameter(location_of_constraint0=None, \
477     value_of_constraint0=None, \
478     location_of_constraint1=None, \
479     value_of_constraint1=None, \
480     location_of_constraint2=None, \
481     value_of_constraint2=None, \
482     location_of_constraint3=None, \
483     value_of_constraint3=None, \
484     location_of_constraint4=None, \
485     value_of_constraint4=None)
486     def location_of_constraint(self):
487     """
488     return the values used to constrain a solution
489    
490     @return: the mask marking the locations of the constraints
491 gross 906 @rtype: L{escript.Scalar}
492 gross 904 """
493     out_loc=0
494 gross 906 if not self.location_of_constraint0 == None:
495 gross 904 out_loc=wherePositive(out_loc+wherePositive(self.location_of_constraint0))
496 gross 906 if not self.location_of_constraint1 == None:
497 gross 904 out_loc=wherePositive(out_loc+wherePositive(self.location_of_constraint1))
498 gross 906 if not self.location_of_constraint2 == None:
499 gross 904 out_loc=wherePositive(out_loc+wherePositive(self.location_of_constraint2))
500 gross 906 if not self.location_of_constraint3 == None:
501 gross 904 out_loc=wherePositive(out_loc+wherePositive(self.location_of_constraint3))
502     return out_loc
503    
504     def value_of_constraint(self):
505     """
506     return the values used to constrain a solution
507    
508     @return: values to be used at the locations of the constraints. If
509     L{value} is not given C{None} is rerturned.
510 gross 906 @rtype: L{escript.Scalar}
511 gross 904 """
512     out_loc=0
513     out=0
514 gross 906 if not self.location_of_constraint0 == None:
515 gross 904 tmp=wherePositive(self.location_of_constraint0)
516 gross 906 out=out*(1.-tmp)+self.value_of_constraint0*tmp
517 gross 904 out_loc=wherePositive(out_loc+tmp)
518 gross 906 if not self.location_of_constraint1 == None:
519 gross 904 tmp=wherePositive(self.location_of_constraint1)
520 gross 906 out=out*(1.-tmp)+self.value_of_constraint1*tmp
521 gross 904 out_loc=wherePositive(out_loc+tmp)
522 gross 906 if not self.location_of_constraint2 == None:
523 gross 904 tmp=wherePositive(self.location_of_constraint2)
524 gross 906 out=out*(1.-tmp)+self.value_of_constraint2*tmp
525 gross 904 out_loc=wherePositive(out_loc+tmp)
526 gross 906 if not self.location_of_constraint3 == None:
527 gross 904 tmp=wherePositive(self.location_of_constraint3)
528 gross 906 out=out*(1.-tmp)+self.value_of_constraint3*tmp
529 gross 904 out_loc=wherePositive(out_loc+tmp)
530     return out
531 jgs 148 # vim: expandtab shiftwidth=4:

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26