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


1
2 ##############################################################################
3 #
4 # Copyright (c) 2003-2020 by The University of Queensland
5 # http://www.uq.edu.au
6 #
7 # Primary Business: Queensland, Australia
8 # Licensed under the Apache License, version 2.0
9 # http://www.apache.org/licenses/LICENSE-2.0
10 #
11 # Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 # Development 2012-2013 by School of Earth Sciences
13 # Development from 2014 by Centre for Geoscience Computing (GeoComp)
14 # Development from 2019 by School of Earth and Environmental Sciences
15 #
16 ##############################################################################
17
18 from __future__ import division, print_function
19
20 __copyright__="""Copyright (c) 2003-2020 by The University of Queensland
21 http://www.uq.edu.au
22 Primary Business: Queensland, Australia"""
23 __license__="""Licensed under the Apache License, version 2.0
24 http://www.apache.org/licenses/LICENSE-2.0"""
25 __url__="https://launchpad.net/escript-finley"
26
27 from esys.escript import length, wherePositive, whereNegative, exp, inf, sup
28 from esys.escript.modelframe import Model,ParameterSet
29 from esys.escript.linearPDEs import LinearPDE
30 from math import log
31 import numpy
32
33 class Sequencer(Model):
34 """
35 Runs through time until t_end is reached.
36
37 :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
44 """
45 def __init__(self,**kwargs):
46 """
47 """
48 super(Sequencer,self).__init__(**kwargs)
49 self.declareParameter(t=0.,
50 t_end=1.,
51 dt_max=Model.UNDEF_DT)
52
53 def doInitialization(self):
54 """
55 initialize time integration
56 """
57 self.__t_old = self.t
58
59 def doStepPreprocessing(self, dt):
60 self.t = self.__t_old+dt
61
62 def doStepPostprocessing(self, dt):
63 self.__t_old = self.t
64
65 def finalize(self):
66 """
67 returns true when `t` has reached `t_end`
68 """
69 return self.t >= self.t_end
70
71 def getSafeTimeStepSize(self, dt):
72 """
73 returns `dt_max`
74 """
75 return self.dt_max
76
77 class GaussianProfile(ParameterSet):
78 """
79 Generates a Gaussian profile at center x_c, width width and height A
80 over a domain
81
82 :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
88 In the case that the spatial dimension is two, The third component of
89 x_c is dropped.
90 """
91 def __init__(self,**kwargs):
92 super(GaussianProfile, self).__init__(**kwargs)
93 self.declareParameter(domain=None,
94 x_c=numpy.zeros([3]),
95 A=1.,
96 width=0.1,
97 r=0)
98
99 def out(self):
100 """
101 Generate the Gaussian profile
102
103 Link against this method to get the output of this model.
104 """
105 x = self.domain.getX()
106 dim = self.domain.getDim()
107 l = length(x-self.x_c[:dim])
108 m = whereNegative(l-self.r)
109
110 return (m+(1.-m)*exp(-log(2.)*(l/self.width)**2))*self.A
111
112 class InterpolateOverBox(ParameterSet):
113 """
114 Returns values at each time. The values are defined through given values
115 at time node. For two dimensional domains back values are ignored.
116
117 :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 """
127
128 def __init__(self, **kwargs):
129 super(InterpolateOverBox, self).__init__(self)
130 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 """
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 x = self.domain.getX()
148 if self.domain.getDim() == 2:
149 x0,x1=x[0],x[1]
150 left_bottom_front0,right_top_back0=inf(x0),sup(x0)
151 left_bottom_front1,right_top_back1=inf(x1),sup(x1)
152 f_right = (x0 - left_bottom_front0)/(right_top_back0 -left_bottom_front0)
153 f_left = 1. - f_right
154 f_top = (x1 - left_bottom_front1)/(right_top_back1 - left_bottom_front1)
155 f_bottom = 1. - f_top
156 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 else:
161 x0,x1,x2=x[0],x[1],x[2]
162 left_bottom_front0,right_top_back0=inf(x0),sup(x0)
163 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 f_left = 1. - f_right
167 f_top = (x1 - left_bottom_front1)/(right_top_back1 - left_bottom_front1)
168 f_bottom = 1. - f_top
169 f_back = (x2 - left_bottom_front1)/(right_top_back2 - left_bottom_front2)
170 f_front = 1. - f_back
171 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 return out
180
181
182 class InterpolatedTimeProfile(ParameterSet):
183 """
184
185 Returns values at each time. The values are defined through given
186 values at time node.
187
188 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
194 :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 """
198
199 def __init__(self,**kwargs):
200 super( InterpolatedTimeProfile, self).__init__(**kwargs)
201 self.declareParameter(t=0., \
202 nodes=[0.,1.],\
203 values=[1.,1.])
204 def out(self):
205 """
206 current value
207
208 Link against this method to get the output of this model.
209 """
210 l = len(self.nodes) - 1
211 t = self.t
212 if t <= self.nodes[0]:
213 return self.values[0]
214 else:
215 for i in range(1,l):
216 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 return self.values[l]
221
222 class ScalarDistributionFromTags(ParameterSet):
223 """
224 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
227 :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 """
271 def __init__(self,**kwargs):
272 super(ScalarDistributionFromTags, self).__init__(**kwargs)
273 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 returns a `esys.escript.Data` object
300 Link against this method to get the output of this model.
301 """
302 d=Scalar(self.default,Function(self.domain))
303 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 return d
314
315 class SmoothScalarDistributionFromTags(ParameterSet):
316 """
317 creates a smooth scalar distribution on a domain from region tags
318
319 :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 """
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 mask.setTaggedValue(tag,1.)
394 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 returns a `esys.escript.Data` object
403 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 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 return d
418
419 class LinearCombination(ParameterSet):
420 """
421 Returns a linear combination of the f0*v0+f1*v1+f2*v2+f3*v3+f4*v4
422
423 :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 """
434 def __init__(self,**kwargs):
435 super(LinearCombination, self).__init__(**kwargs)
436 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
447 def out(self):
448 """
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 if not self.f0 is None and not self.v0 is None:
453 fv0 = self.f0*self.v0
454 else:
455 fv0 = None
456
457 if not self.f1 is None and not self.v1 is None:
458 fv1 = self.f1*self.v1
459 else:
460 fv1 = None
461
462 if not self.f2 is None and not self.v2 is None:
463 fv2 = f2*v2
464 else:
465 fv2 = None
466
467 if not self.f3 is None and not self.v3 is None:
468 fv3 = self.f3*self.v3
469 else:
470 fv3 = None
471
472 if not self.f4 is None and not self.v4 is None:
473 fv4 = self.f4*self.v4
474 else:
475 fv4 = None
476
477 if fv0 is None:
478 out = 0.
479 else:
480 out = fv0
481 if not fv1 is None:
482 out += fv1
483 if not fv2 is None:
484 out += fv2
485 if not fv3 is None:
486 out += fv3
487 return out
488
489 class MergeConstraints(ParameterSet):
490 """
491 Returns a linear combination of the f0*v0+f1*v1+f2*v2+f3*v3+f4*v4
492 """
493 def __init__(self,**kwargs):
494 super(MergeConstraints, self).__init__(**kwargs)
495 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 :return: the mask marking the locations of the constraints
510 :rtype: `escript.Scalar`
511 """
512 out_loc=0
513 if not self.location_of_constraint0 is None:
514 out_loc=wherePositive(out_loc+wherePositive(self.location_of_constraint0))
515 if not self.location_of_constraint1 is None:
516 out_loc=wherePositive(out_loc+wherePositive(self.location_of_constraint1))
517 if not self.location_of_constraint2 is None:
518 out_loc=wherePositive(out_loc+wherePositive(self.location_of_constraint2))
519 if not self.location_of_constraint3 is None:
520 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 :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 """
531 out_loc=0
532 out=0
533 if not self.location_of_constraint0 is None:
534 tmp=wherePositive(self.location_of_constraint0)
535 out=out*(1.-tmp)+self.value_of_constraint0*tmp
536 out_loc=wherePositive(out_loc+tmp)
537 if not self.location_of_constraint1 is None:
538 tmp=wherePositive(self.location_of_constraint1)
539 out=out*(1.-tmp)+self.value_of_constraint1*tmp
540 out_loc=wherePositive(out_loc+tmp)
541 if not self.location_of_constraint2 is None:
542 tmp=wherePositive(self.location_of_constraint2)
543 out=out*(1.-tmp)+self.value_of_constraint2*tmp
544 out_loc=wherePositive(out_loc+tmp)
545 if not self.location_of_constraint3 is None:
546 tmp=wherePositive(self.location_of_constraint3)
547 out=out*(1.-tmp)+self.value_of_constraint3*tmp
548 out_loc=wherePositive(out_loc+tmp)
549 return out
550 # 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