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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 953 - (show annotations)
Fri Feb 9 08:55:28 2007 UTC (12 years, 7 months ago) by gross
File MIME type: text/x-python
File size: 21138 byte(s)
TagMap is supported by modellib now
1 # $Id$
2
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 from esys.escript import *
10 from esys.escript.modelframe import Model,ParameterSet
11 from esys.escript.linearPDEs import LinearPDE
12 from math import log
13
14 class Sequencer(Model):
15 """
16 Runs through time until t_end is reached.
17
18 @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 """
26 def __init__(self,**kwargs):
27 """
28 """
29 super(Sequencer,self).__init__(**kwargs)
30 self.declareParameter(t=0.,
31 t_end=1.,
32 dt_max=Model.UNDEF_DT)
33
34 def doInitialization(self):
35 """
36 initialize time integration
37 """
38 self.__t_old = self.t
39
40 def doStepPreprocessing(self, dt):
41 self.t = self.__t_old+dt
42
43 def doStepPostprocessing(self, dt):
44 self.__t_old = self.t
45
46 def finalize(self):
47 """
48 returns true when L{t} has reached L{t_end}
49 """
50 return self.t >= self.t_end
51
52 def getSafeTimeStepSize(self, dt):
53 """
54 returns L{dt_max}
55 """
56 return self.dt_max
57
58 class GaussianProfile(ParameterSet):
59 """
60 Generates a Gaussian profile at center x_c, width width and height A
61 over a domain
62
63 @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
69 In the case that the spatial dimension is two, The third component of
70 x_c is dropped.
71 """
72 def __init__(self,**kwargs):
73 super(GaussianProfile, self).__init__(**kwargs)
74 self.declareParameter(domain=None,
75 x_c=numarray.zeros([3]),
76 A=1.,
77 width=0.1,
78 r=0)
79
80 def out(self):
81 """
82 Generate the Gaussian profile
83
84 Link against this method to get the output of this model.
85 """
86 x = self.domain.getX()
87 dim = self.domain.getDim()
88 l = length(x-self.x_c[:dim])
89 m = whereNegative(l-self.r)
90
91 return (m+(1.-m)*exp(-log(2.)*(l/self.width)**2))*self.A
92
93 class InterpolateOverBox(ParameterSet):
94 """
95 Returns values at each time. The values are defined through given values
96 at time node. For two dimensional domains back values are ignored.
97
98 @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 """
108
109 def __init__(self, **kwargs):
110 super(InterpolateOverBox, self).__init__(self)
111 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 """
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 x = self.domain.getX()
129 if self.domain.getDim() == 2:
130 x0,x1=x[0],x[1]
131 left_bottom_front0,right_top_back0=inf(x0),sup(x0)
132 left_bottom_front1,right_top_back1=inf(x1),sup(x1)
133 f_right = (x0 - left_bottom_front0)/(right_top_back0 -left_bottom_front0)
134 f_left = 1. - f_right
135 f_top = (x1 - left_bottom_front1)/(right_top_back1 - left_bottom_front1)
136 f_bottom = 1. - f_top
137 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 else:
142 x0,x1,x2=x[0],x[1],x[2]
143 left_bottom_front0,right_top_back0=inf(x0),sup(x0)
144 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 f_left = 1. - f_right
148 f_top = (x1 - left_bottom_front1)/(right_top_back1 - left_bottom_front1)
149 f_bottom = 1. - f_top
150 f_back = (x2 - left_bottom_front1)/(right_top_back2 - left_bottom_front2)
151 f_front = 1. - f_back
152 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 return out
161
162
163 class InterpolatedTimeProfile(ParameterSet):
164 """
165
166 Returns values at each time. The values are defined through given
167 values at time node.
168
169 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
175 @ivar t: (in) current time
176 @ivar node: (in) list of time nodes
177 @ivar values: (in) list of values at time nodes
178 """
179
180 def __init__(self,**kwargs):
181 super( InterpolatedTimeProfile, self).__init__(**kwargs)
182 self.declareParameter(t=0., \
183 nodes=[0.,1.],\
184 values=[1.,1.])
185 def out(self):
186 """
187 current value
188
189 Link against this method to get the output of this model.
190 """
191 l = len(self.nodes) - 1
192 t = self.t
193 if t <= self.nodes[0]:
194 return self.values[0]
195 else:
196 for i in range(1,l):
197 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 return self.values[l]
202
203 class ScalarDistributionFromTags(ParameterSet):
204 """
205 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
208 @ivar domain: domain
209 @type domain: L{esys.escript.Domain}
210 @ivar tag_map: maping from names to tags
211 @type tag_map: L{esys.pycad.TagMap}
212 @ivar default: default value
213 @ivar tag0: tag 0
214 @type tag0: C{int}
215 @ivar value0: value for tag 0
216 @type value0: C{float}
217 @ivar tag1: tag 1
218 @type tag1: C{int}
219 @ivar value1: value for tag 1
220 @type value1: C{float}
221 @ivar tag2: tag 2
222 @type tag2: C{int}
223 @ivar value2: value for tag 2
224 @type value2: C{float}
225 @ivar tag3: tag 3
226 @type tag3: C{int}
227 @ivar value3: value for tag 3
228 @type value3: C{float}
229 @ivar tag4: tag 4
230 @type tag4: C{int}
231 @ivar value4: value for tag 4
232 @type value4: C{float}
233 @ivar tag5: tag 5
234 @type tag5: C{int}
235 @ivar value5: value for tag 5
236 @type value5: C{float}
237 @ivar tag6: tag 6
238 @type tag6: C{int}
239 @ivar value6: value for tag 6
240 @type value6: C{float}
241 @ivar tag7: tag 7
242 @type tag7: C{int}
243 @ivar value7: value for tag 7
244 @type value7: C{float}
245 @ivar tag8: tag 8
246 @type tag8: C{int}
247 @ivar value8: value for tag 8
248 @type value8: C{float}
249 @ivar tag9: tag 9
250 @type tag9: C{int}
251 @ivar value9: value for tag 9
252 @type value9: C{float}
253 """
254 def __init__(self,**kwargs):
255 super(ScalarDistributionFromTags, self).__init__(**kwargs)
256 self.declareParameter(domain=None,
257 tag_map=None,
258 default=0.,
259 tag0=None,
260 value0=0.,
261 tag1=None,
262 value1=0.,
263 tag2=None,
264 value2=0.,
265 tag3=None,
266 value3=0.,
267 tag4=None,
268 value4=0.,
269 tag5=None,
270 value5=0.,
271 tag6=None,
272 value6=0.,
273 tag7=None,
274 value7=0.,
275 tag8=None,
276 value8=0.,
277 tag9=None,
278 value9=0.)
279
280
281 def out(self):
282 """
283 returns a L{esys.escript.Data} object
284 Link against this method to get the output of this model.
285 """
286 d=Scalar(self.default,Function(self.domain))
287 if self.tag_map == None:
288 if not self.tag0 == None: d.setTaggedValue(self.tag0,self.value0)
289 if not self.tag1 == None: d.setTaggedValue(self.tag1,self.value1)
290 if not self.tag2 == None: d.setTaggedValue(self.tag2,self.value2)
291 if not self.tag3 == None: d.setTaggedValue(self.tag3,self.value3)
292 if not self.tag4 == None: d.setTaggedValue(self.tag4,self.value4)
293 if not self.tag5 == None: d.setTaggedValue(self.tag5,self.value5)
294 if not self.tag6 == None: d.setTaggedValue(self.tag6,self.value6)
295 if not self.tag7 == None: d.setTaggedValue(self.tag7,self.value7)
296 if not self.tag8 == None: d.setTaggedValue(self.tag8,self.value8)
297 if not self.tag9 == None: d.setTaggedValue(self.tag9,self.value9)
298 else:
299 args={}
300 if not self.tag0 == None: args[self.tag0]=self.value0
301 if not self.tag1 == None: args[self.tag1]=self.value1
302 if not self.tag2 == None: args[self.tag2]=self.value2
303 if not self.tag3 == None: args[self.tag3]=self.value3
304 if not self.tag4 == None: args[self.tag4]=self.value4
305 if not self.tag5 == None: args[self.tag5]=self.value5
306 if not self.tag6 == None: args[self.tag6]=self.value6
307 if not self.tag7 == None: args[self.tag7]=self.value7
308 if not self.tag8 == None: args[self.tag8]=self.value8
309 if not self.tag9 == None: args[self.tag9]=self.value9
310 self.tag_map.insert(d,**args)
311 return d
312
313 class SmoothScalarDistributionFromTags(ParameterSet):
314 """
315 creates a smooth scalar distribution on a domain from region tags
316
317 @ivar domain: domain
318 @type domain: L{esys.escript.Domain}
319 @ivar tag_map: maping from names to tags
320 @type tag_map: L{esys.pycad.TagMap}
321 @ivar default: default value
322 @ivar tag0: tag 0
323 @type tag0: C{int}
324 @ivar value0: value for tag 0
325 @type value0: C{float}
326 @ivar tag1: tag 1
327 @type tag1: C{int}
328 @ivar value1: value for tag 1
329 @type value1: C{float}
330 @ivar tag2: tag 2
331 @type tag2: C{int}
332 @ivar value2: value for tag 2
333 @type value2: C{float}
334 @ivar tag3: tag 3
335 @type tag3: C{int}
336 @ivar value3: value for tag 3
337 @type value3: C{float}
338 @ivar tag4: tag 4
339 @type tag4: C{int}
340 @ivar value4: value for tag 4
341 @type value4: C{float}
342 @ivar tag5: tag 5
343 @type tag5: C{int}
344 @ivar value5: value for tag 5
345 @type value5: C{float}
346 @ivar tag6: tag 6
347 @type tag6: C{int}
348 @ivar value6: value for tag 6
349 @type value6: C{float}
350 @ivar tag7: tag 7
351 @type tag7: C{int}
352 @ivar value7: value for tag 7
353 @type value7: C{float}
354 @ivar tag8: tag 8
355 @type tag8: C{int}
356 @ivar value8: value for tag 8
357 @type value8: C{float}
358 @ivar tag9: tag 9
359 @type tag9: C{int}
360 @ivar value9: value for tag 9
361 @type value9: C{float}
362 """
363 def __init__(self,**kwargs):
364 super(SmoothScalarDistributionFromTags, self).__init__(**kwargs)
365 self.declareParameter(domain=None,
366 default=0.,
367 tag0=None,
368 value0=0.,
369 tag1=None,
370 value1=0.,
371 tag2=None,
372 value2=0.,
373 tag3=None,
374 value3=0.,
375 tag4=None,
376 value4=0.,
377 tag5=None,
378 value5=0.,
379 tag6=None,
380 value6=0.,
381 tag7=None,
382 value7=0.,
383 tag8=None,
384 value8=0.,
385 tag9=None,
386 value9=0.)
387
388
389 def __update(self,tag,tag_value,value):
390 if self.__pde==None:
391 self.__pde=LinearPDE(self.domain,numSolutions=1)
392 mask=Scalar(0.,Function(self.domain))
393 if self.tag_map == None:
394 mask.setTaggedValue(tag,1.)
395 else:
396 self.tag_map.insert(mask,**{tag:tag_value})
397 self.__pde.setValue(Y=mask)
398 mask=wherePositive(abs(self.__pde.getRightHandSide()))
399 value*=(1.-mask)
400 value+=tag_value*mask
401 return value
402
403 def out(self):
404 """
405 returns a L{esys.escript.Data} object
406 Link against this method to get the output of this model.
407 """
408 d=Scalar(self.default,Solution(self.domain))
409 self.__pde=None
410 if not self.tag0 == None: d=self.__update(self.tag0,self.value0,d)
411 if not self.tag1 == None: d=self.__update(self.tag1,self.value1,d)
412 if not self.tag2 == None: d=self.__update(self.tag2,self.value2,d)
413 if not self.tag3 == None: d=self.__update(self.tag3,self.value3,d)
414 if not self.tag4 == None: d=self.__update(self.tag4,self.value4,d)
415 if not self.tag5 == None: d=self.__update(self.tag5,self.value5,d)
416 if not self.tag6 == None: d=self.__update(self.tag6,self.value6,d)
417 if not self.tag7 == None: d=self.__update(self.tag7,self.value7,d)
418 if not self.tag8 == None: d=self.__update(self.tag8,self.value8,d)
419 if not self.tag9 == None: d=self.__update(self.tag9,self.value9,d)
420 return d
421
422 class LinearCombination(ParameterSet):
423 """
424 Returns a linear combination of the f0*v0+f1*v1+f2*v2+f3*v3+f4*v4
425
426 @ivar f0: numerical object or None, default=None (in)
427 @ivar v0: numerical object or None, default=None (in)
428 @ivar f1: numerical object or None, default=None (in)
429 @ivar v1: numerical object or None, default=None (in)
430 @ivar f2: numerical object or None, default=None (in)
431 @ivar v2: numerical object or None, default=None (in)
432 @ivar f3: numerical object or None, default=None (in)
433 @ivar v3: numerical object or None, default=None (in)
434 @ivar f4: numerical object or None, default=None (in)
435 @ivar v4: numerical object or None, default=None (in)
436 """
437 def __init__(self,**kwargs):
438 super(LinearCombination, self).__init__(**kwargs)
439 self.declareParameter(f0=None, \
440 v0=None, \
441 f1=None, \
442 v1=None, \
443 f2=None, \
444 v2=None, \
445 f3=None, \
446 v3=None, \
447 f4=None, \
448 v4=None)
449
450 def out(self):
451 """
452 returns f0*v0+f1*v1+f2*v2+f3*v3+f4*v4.
453 Link against this method to get the output of this model.
454 """
455 if not self.f0 == None and not self.v0 == None:
456 fv0 = self.f0*self.v0
457 else:
458 fv0 = None
459
460 if not self.f1 == None and not self.v1 == None:
461 fv1 = self.f1*self.v1
462 else:
463 fv1 = None
464
465 if not self.f2 == None and not self.v2 == None:
466 fv2 = f2*v2
467 else:
468 fv2 = None
469
470 if not self.f3 == None and not self.v3 == None:
471 fv3 = self.f3*self.v3
472 else:
473 fv3 = None
474
475 if not self.f4 == None and not self.v4 == None:
476 fv4 = self.f4*self.v4
477 else:
478 fv4 = None
479
480 if fv0 == None:
481 out = 0.
482 else:
483 out = fv0
484 if not fv1 == None:
485 out += fv1
486 if not fv2 == None:
487 out += fv2
488 if not fv3 == None:
489 out += fv3
490 return out
491
492 class MergeConstraints(ParameterSet):
493 """
494 Returns a linear combination of the f0*v0+f1*v1+f2*v2+f3*v3+f4*v4
495 """
496 def __init__(self,**kwargs):
497 super(MergeConstraints, self).__init__(**kwargs)
498 self.declareParameter(location_of_constraint0=None, \
499 value_of_constraint0=None, \
500 location_of_constraint1=None, \
501 value_of_constraint1=None, \
502 location_of_constraint2=None, \
503 value_of_constraint2=None, \
504 location_of_constraint3=None, \
505 value_of_constraint3=None, \
506 location_of_constraint4=None, \
507 value_of_constraint4=None)
508 def location_of_constraint(self):
509 """
510 return the values used to constrain a solution
511
512 @return: the mask marking the locations of the constraints
513 @rtype: L{escript.Scalar}
514 """
515 out_loc=0
516 if not self.location_of_constraint0 == None:
517 out_loc=wherePositive(out_loc+wherePositive(self.location_of_constraint0))
518 if not self.location_of_constraint1 == None:
519 out_loc=wherePositive(out_loc+wherePositive(self.location_of_constraint1))
520 if not self.location_of_constraint2 == None:
521 out_loc=wherePositive(out_loc+wherePositive(self.location_of_constraint2))
522 if not self.location_of_constraint3 == None:
523 out_loc=wherePositive(out_loc+wherePositive(self.location_of_constraint3))
524 return out_loc
525
526 def value_of_constraint(self):
527 """
528 return the values used to constrain a solution
529
530 @return: values to be used at the locations of the constraints. If
531 L{value} is not given C{None} is rerturned.
532 @rtype: L{escript.Scalar}
533 """
534 out_loc=0
535 out=0
536 if not self.location_of_constraint0 == None:
537 tmp=wherePositive(self.location_of_constraint0)
538 out=out*(1.-tmp)+self.value_of_constraint0*tmp
539 out_loc=wherePositive(out_loc+tmp)
540 if not self.location_of_constraint1 == None:
541 tmp=wherePositive(self.location_of_constraint1)
542 out=out*(1.-tmp)+self.value_of_constraint1*tmp
543 out_loc=wherePositive(out_loc+tmp)
544 if not self.location_of_constraint2 == None:
545 tmp=wherePositive(self.location_of_constraint2)
546 out=out*(1.-tmp)+self.value_of_constraint2*tmp
547 out_loc=wherePositive(out_loc+tmp)
548 if not self.location_of_constraint3 == None:
549 tmp=wherePositive(self.location_of_constraint3)
550 out=out*(1.-tmp)+self.value_of_constraint3*tmp
551 out_loc=wherePositive(out_loc+tmp)
552 return out
553 # 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