/[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 6939 - (hide annotations)
Mon Jan 20 03:37:18 2020 UTC (9 months, 1 week ago) by uqaeller
File MIME type: text/x-python
File size: 21145 byte(s)
Updated the copyright header.


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