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

Contents of /trunk/escript/py_src/levelset.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3981 - (show annotations)
Fri Sep 21 02:47:54 2012 UTC (7 years ago) by jfenwick
File MIME type: text/x-python
File size: 21909 byte(s)
First pass of updating copyright notices
1
2 ##############################################################################
3 #
4 # Copyright (c) 2003-2012 by University of Queensland
5 # http://www.uq.edu.au
6 #
7 # Primary Business: Queensland, Australia
8 # Licensed under the Open Software License version 3.0
9 # http://www.opensource.org/licenses/osl-3.0.php
10 #
11 # Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 # Development since 2012 by School of Earth Sciences
13 #
14 ##############################################################################
15
16 __copyright__="""Copyright (c) 2003-2012 by University of Queensland
17 http://www.uq.edu.au
18 Primary Business: Queensland, Australia"""
19 __license__="""Licensed under the Open Software License version 3.0
20 http://www.opensource.org/licenses/osl-3.0.php"""
21 __url__="https://launchpad.net/escript-finley"
22
23 #from esys.escript import *
24 from esys.escript.linearPDEs import LinearPDE, SingleTransportPDE
25 from esys.escript.pdetools import Projector
26
27 import math
28
29 USE_OLD_VERSION=False
30 USE_OLD_VERSION_REINIT=False
31
32
33 class LevelSet:
34 """
35 The level set method tracking an interface defined by the zero contour of the
36 level set function phi which defines the signed distance of a point x from the
37 interface. The contour phi(x)=0 defines the interface.
38
39 It is assumed that phi(x)<0 defines the volume of interest,
40 phi(x)>0 the outside world.
41 """
42 def __init__(self,phi,reinit_max=10,reinitialize_after=20,smooth=2., useReducedOrder=False):
43 """
44 Sets up the level set method.
45
46 :param phi: the initial level set function
47 :param reinit_max: maximum number of reinitialization steps
48 :param reinitialize_after: ``phi`` is reinitialized every ``reinit_after`` step
49 :param smooth: smoothing width
50 """
51 self.__domain = phi.getDomain()
52 self.__phi = phi
53
54 if USE_OLD_VERSION:
55 self.__PDE = LinearPDE(self.__domain)
56 if useReducedOrder: self.__PDE.setReducedOrderOn()
57 self.__PDE.setSolverMethod(solver=LinearPDE.PCG)
58 self.__PDE.setValue(D=1.0)
59 else:
60 self.__PDE=SingleTransportPDE(self.__domain)
61 if useReducedOrder: self.__PDE.setReducedOrderOn()
62 self.__PDE.setValue(M=1.0)
63
64 if USE_OLD_VERSION_REINIT:
65 self.__reinitPDE = LinearPDE(self.__domain, numEquations=1)
66 self.__reinitPDE.getSolverOption().setSolverMethod(solver=LinearPDE.LUMPING)
67 if useReducedOrder: self.__reinitPDE.setReducedOrderOn()
68 self.__reinitPDE.setValue(D=1.0)
69 else:
70 self.__reinitPDE=SingleTransportPDE(self.__domain)
71 if useReducedOrder: self.__reinitPDE.setReducedOrderOn()
72 self.__reinitPDE.setValue(M=1.0)
73
74 # revise:
75 self.__reinit_max = reinit_max
76 self.__reinit_after = reinitialize_after
77 self.__h = inf(self.__domain.getSize())
78 self.__smooth = smooth
79 self.__n_step=0
80
81 def getH(self):
82 """
83 Returns the mesh size.
84 """
85 return self.__h
86
87 def getDomain(self):
88 """
89 Returns the domain.
90 """
91 return self.__domain
92
93 def getAdvectionSolverOptions(self):
94 """
95 Returns the solver options for the interface advective.
96 """
97 return self.__PDE.getSolverOptions()
98
99 def getReinitializationSolverOptions(self):
100 """
101 Returns the options of the solver for the reinitialization
102 """
103 return self.__reinitPDE.getSolverOption()
104
105 def getLevelSetFunction(self):
106 """
107 Returns the level set function.
108 """
109 return self.__phi
110
111 def getTimeStepSize(self,velocity):
112 """
113 Returns a new ``dt`` for a given ``velocity`` using the Courant condition.
114
115 :param velocity: velocity field
116 """
117 if USE_OLD_VERSION:
118 self.__velocity=velocity
119 dt=0.5*self.__h/sup(length(velocity))
120 else:
121 self.__PDE.setValue(C=-velocity)
122 dt=self.__PDE.getSafeTimeStepSize()
123
124 return dt
125
126 def update(self,dt):
127 """
128 Sets a new velocity and updates the level set function.
129
130 :param dt: time step forward
131 """
132 self.__phi=self.__advect(dt)
133 self.__n_step+=1
134 if self.__n_step%self.__reinit_after ==0: self.__phi = self.__reinitialise()
135 return self.__phi
136
137 def update_phi(self, velocity, dt):
138 """
139 Updates ``phi`` under the presence of a velocity field.
140
141 If dt is small this call is equivalent to::
142
143 dt=LevelSet.getTimeStepSize(velocity)
144 phi=LevelSet.update(dt)
145
146 otherwise substepping is used.
147
148 :param velocity: velocity field
149 :param dt: time step forward
150 """
151 dt2=self.getTimeStepSize(velocity)
152 n=int(math.ceil(dt/dt2))
153 dt_new=dt/n
154 for i in range(n):
155 phi=self.update(dt_new)
156 return phi
157
158 def __advect(self, dt):
159 """
160 Advects the level set function in the presence of a velocity field.
161
162 :param dt: time increment
163 :return: the advected level set function
164 """
165 if USE_OLD_VERSION:
166 velocity=self.__velocity
167 Y = self.__phi-dt/2.*inner(velocity,grad(self.__phi))
168 self.__PDE.setValue(Y=Y)
169 phi_half = self.__PDE.getSolution()
170 Y = self.__phi-dt*inner(velocity,grad(phi_half))
171 self.__PDE.setValue(Y=Y)
172 phi = self.__PDE.getSolution()
173 else:
174 phi = self.__PDE.getSolution(dt, self.__phi)
175 print("LevelSet: Advection step done")
176 return phi
177
178 #==============================================================================================
179 def __reinitialise(self):
180 """
181 Reinitializes the level set.
182
183 It solves the PDE...
184
185 :return: reinitialized level set
186 """
187 phi=self.__phi
188 s = sign(phi.interpolate(Function(self.__domain)))
189 g=grad(phi)
190 w = s*g/(length(g)+1e-10)
191 dtau = 0.5*self.__h
192 iter =0
193 mask = whereNegative(abs(phi)-1.*self.__h)
194 self.__reinitPDE.setValue(q=mask, r=phi)
195 g=grad(phi)
196 while (iter<=self.__reinit_max):
197 Y = phi+(dtau/2.)*(s-inner(w,g))
198 self.__reinitPDE.setValue(Y = Y)
199 phi_half = self.__reinitPDE.getSolution()
200 Y = phi+dtau*(s-inner(w,grad(phi_half)))
201 self.__reinitPDE.setValue(Y = Y)
202 phi = self.__reinitPDE.getSolution()
203 g=grad(phi)
204 error = Lsup(length(g)*whereNegative(abs(phi.interpolate(g.getFunctionSpace()))-3.0*self.__h))
205 print(("LevelSet:reinitialization: iteration :", iter, " error:", error))
206 iter +=1
207 return phi
208
209
210 def getVolume(self):
211 """
212 Returns the volume of the *phi(x)<0* region.
213 """
214 return integrate(whereNegative(self.__phi.interpolate(Function(self.__domain))))
215
216
217 def update_parameter_sharp(self, param_neg=-1, param_pos=1, phi=None):
218 """
219 Creates a function with ``param_neg`` where ``phi<0`` and ``param_pos``
220 where ``phi>0`` (no smoothing).
221
222 :param param_neg: value of parameter on the negative side (phi<0)
223 :param param_pos: value of parameter on the positive side (phi>0)
224 :param phi: level set function to be used. If not present the current
225 level set is used.
226 """
227 mask_neg = whereNegative(self.__phi)
228 mask_pos = whereNonNegative(self.__phi)
229 param = param_pos*mask_pos + param_neg*mask_neg
230 return param
231
232 def update_parameter(self, param_neg=-1, param_pos=1, phi=None, smoothing_width=None):
233 """
234 Creates a smoothed function with ``param_neg`` where ``phi<0`` and
235 ``param_pos`` where ``phi>0`` which is smoothed over a length
236 ``smoothing_width`` across the interface.
237
238 :param smoothing_width: width of the smoothing zone relative to mesh size.
239 If not present the initial value of ``smooth`` is
240 used.
241 """
242 if smoothing_width==None: smoothing_width = self.__smooth
243 if phi==None: phi=self.__phi
244 s=self.__makeInterface(phi,smoothing_width)
245 return ((param_pos-param_neg)*s+param_pos+param_neg)/2
246
247 def __makeInterface(self,phi,smoothing_width):
248 """
249 Creates a smooth interface from -1 to 1 over the length
250 *2*h*smoothing_width* where -1 is used where the level set is negative
251 and 1 where the level set is 1.
252 """
253 s=smoothing_width*self.__h
254 phi_on_h=interpolate(phi,Function(self.__domain))
255 mask_neg = whereNonNegative(-s-phi_on_h)
256 mask_pos = whereNonNegative(phi_on_h-s)
257 mask_interface = 1.-mask_neg-mask_pos
258 interface=phi_on_h/s
259 return - mask_neg + mask_pos + mask_interface * interface
260
261 def makeCharacteristicFunction(self, contour=0, phi=None, positiveSide=True, smoothing_width=None):
262 """
263 Makes a smooth characteristic function of the region ``phi(x)>contour`` if
264 ``positiveSide`` and ``phi(x)<contour`` otherwise.
265
266 :param phi: level set function to be used. If not present the current
267 level set is used.
268 :param smoothing_width: width of the smoothing zone relative to mesh size.
269 If not present the initial value of ``smooth`` is
270 used.
271 """
272 if phi==None: phi=self.__phi
273 if smoothing_width == None: smoothing_width=self.__smooth
274 s=self.__makeInterface(phi=phi-contour,smoothing_width=smoothing_width)
275 if positiveSide:
276 return (1+s)/2
277 else:
278 return (1-s)/2
279
280
281 class LevelSet2(object):
282 def __init__(self,phi,reinit_max=10,reinit_after=2,smooth=2.):
283 """
284 Initializes the model.
285 """
286 self.__domain = phi.getDomain()
287 x=self.__domain.getX()
288 diam=0
289 for i in range(self.__domain.getDim()):
290 xi=x[i]
291 diam+=(inf(xi)-sup(xi))**2
292 self.__diam=sqrt(diam)
293 self.__h = sup(Function(self.__domain).getSize())
294 self.__reinit_after=reinit_after
295 self.__reinit_max = reinit_max
296 self.__smooth = smooth
297 self.__phi = phi
298 self.__update_count=0
299 self.velocity = None
300
301 #==================================================
302 self.__FC=True
303 if self.__FC:
304 self.__fc=TransportPDE(self.__domain,num_equations=1)
305 # self.__fc.setReducedOrderOn()
306 self.__fc.setValue(M=Scalar(1.,Function(self.__domain)))
307 self.__fc.setInitialSolution(phi+self.__diam)
308 else:
309 self.__fcpde = LinearPDE(self.__domain,numEquations=1, numSolutions=1)
310 self.__fcpde.setSolverMethod(solver=LinearPDE.LUMPING)
311 self.__fcpde.setReducedOrderOn()
312 self.__fcpde.setSymmetryOn()
313 self.__fcpde.setValue(D=1.)
314
315 #=======================================
316 self.__reinitFC=False
317 if self.__reinitFC:
318 self.__reinitfc=TransportPDE(self.__domain,num_equations=1)
319 self.__reinitfc.setValue(M=Scalar(1.,Function(self.__domain)))
320 self.__reinitfc.setTolerance(1.e-5)
321 else:
322 self.__reinitPde = LinearPDE(self.__domain,numEquations=1, numSolutions=1)
323 # self.__reinitPde.setSolverMethod(solver=LinearPDE.LUMPING)
324 self.__reinitPde.setReducedOrderOn()
325 self.__reinitPde.setSymmetryOn()
326 self.__reinitPde.setValue(D=1.)
327 self.__reinitPde.setTolerance(1.e-2)
328 # self.__reinitPde.setSolverMethod(preconditioner=LinearPDE.ILU0)
329 #=======================================
330 self.__updateInterface()
331 print(("phi range:",inf(phi), sup(phi)))
332
333 def setFixedRegion(self,mask,contour=0):
334 q=whereNonPositive(abs(self.__phi-contour)-self.__h)*mask
335 # self.__fc.setValue(q=q)
336 self.__reinitPde.setValue(q=q)
337
338
339 def __updateInterface(self):
340 self.__smoothed_char=self.__makeInterface(self.__phi)
341
342
343 def update(self,dt):
344 """
345 Sets a new velocity and updates the level set function.
346
347 :param dt: time step forward
348 """
349 if self.__FC:
350 self.__phi=self.__fc.solve(dt, verbose=False)-self.__diam
351 else:
352 self.__fcpde.setValue(Y = self.__phi-(dt/2.)*inner(self.velocity,grad(self.__phi)))
353 phi_half = self.__fcpde.getSolution()
354 self.__fcpde.setValue(Y = self.__phi-dt*inner(self.velocity,grad(phi_half)))
355 self.__phi= self.__fcpde.getSolution()
356 self.__update_count += 1
357 if self.__update_count%self.__reinit_after == 0:
358 self.__phi=self.__reinitialise(self.__phi)
359 if self.__FC:
360 self.__fc.setInitialSolution(self.__phi+self.__diam)
361 self.__updateInterface()
362
363 def setTolerance(self,tolerance=1e-3):
364 if self.__FC:
365 self.__fc.setTolerance(tolerance)
366
367 def __reinitialise(self,phi):
368 print("reintialization started:")
369 s=self.__makeInterface(phi,1.)
370 print(("phi range:",inf(phi), sup(phi)))
371 g=grad(phi)
372 w=s*g/(length(g)+EPSILON)
373 #=======================================================
374 # # positive part:
375 # phi_p=wherePositive(phi)*phi
376 # self.__reinitfc.setInitialSolution(phi_p)
377 # self.__reinitfc.setValue(C=-wherePositive(s)*w,Y=wherePositive(s)*s,q=whereNonPositive(phi))
378 # dtau=self.__h
379 # print "step size: dt (pos)= ",dtau
380 # print "phi_p range:",inf(phi_p), sup(phi_p)
381 # iter=0
382 # while (iter<=self.__reinit_max):
383 # phi_p=self.__reinitfc.solve(dtau)
384 # print "phi_p range:",inf(phi_p), sup(phi_p)
385 # iter+=1
386 # # negative part:
387 # phi_n=-whereNegative(phi)*phi
388 # self.__reinitfc.setInitialSolution(phi_n)
389 # self.__reinitfc.setValue(C=-whereNegative(s)*w,Y=-whereNegative(s)*s,q=whereNonNegative(phi))
390 # # dtau=self.__reinitfc.getSafeTimeStepSize()
391 # dtau=self.__h
392 # print "step size: dt (neg)= ",dtau
393 # print "phi_n range:",inf(phi_n), sup(phi_n)
394 # iter=0
395 # while (iter<=self.__reinit_max):
396 # phi_n=self.__reinitfc.solve(dtau)
397 # print "phi_n range:",inf(phi_n), sup(phi_n)
398 # iter+=1
399 # phi=phi_p-phi_n
400 # print "phi range:",inf(phi), sup(phi)
401 # print "reintialization completed."
402 # return phi
403
404 #=======================================================
405 if self.__reinitFC:
406 self.__reinitfc.setValue(C=-w,Y=s)
407 self.__reinitfc.setInitialSolution(phi+self.__diam)
408 dtau=self.__reinitfc.getSafeTimeStepSize()
409 # dtau = 0.5*inf(Function(self.__domain).getSize())
410 else:
411 dtau = 0.5*inf(Function(self.__domain).getSize())
412 print(("step size: dt = ",dtau,inf(abs(phi.interpolate(Function(self.__domain)))/abs(s-inner(w,grad(phi))))))
413 iter =0
414 # self.__reinitPde.setValue(q=whereNegative(abs(phi)-2*self.__h), r=phi)
415 # self.__reinitPde.setValue(r=phi)
416 while (iter<=self.__reinit_max):
417 phi_old=phi
418 if self.__reinitFC:
419 phi = self.__reinitfc.solve(dtau)-self.__diam
420 else:
421 self.__reinitPde.setValue(Y = phi+(dtau/2.)*(s-inner(w,grad(phi))))
422 phi_half = self.__reinitPde.getSolution()
423 self.__reinitPde.setValue(Y = phi+dtau*(s-inner(w,grad(phi_half))))
424 # g=grad(phi)
425 # S=inner(w,grad(phi))
426 # self.__reinitPde.setValue(Y = phi+dtau*(s-S),X=dtau**2/2*w*(s-S))
427 phi = self.__reinitPde.getSolution()
428 change = Lsup(phi-phi_old)/self.__diam
429 print(("phi range:",inf(phi), sup(phi)))
430 print(("iteration :", iter, " change:", change))
431 iter +=1
432 print("reintialization completed.")
433 return phi
434
435 def createParameter(self,value_negative=-1.,value_positive=1):
436 out = (value_negative+value_positive)/2. + self.__smoothed_char * ((value_positive-value_negative)/2.)
437 return out
438
439
440 #=========== things from here onwards are not used nor tested: ===========
441
442
443 def getCharacteristicFunction(self):
444 return self.__smoothed_char
445
446
447 def RK2(self,L):
448 k0=L(phi)
449 phi_1=phi+dt*k0
450 k1=L(phi_1)
451 phi=phi+dt/2*(k0+k1)
452
453 def RK3(self,L):
454 k0=L(phi)
455 phi_1=phi+dt*L(phi)
456 k1=L(phi_1)
457 phi_2=phi+dt/4*(k0+k1)
458 k2=L(phi_2)
459 phi=phi+dt/6*(k0+4*k2+K1)
460
461 def TG(self,L):
462 k0=L(phi)
463 phi_1=phi+dt/2*k0
464 k1=L(phi_1)
465 phi=phi+dt*k1
466
467
468 def __reinitialise_old(self):
469 #=============================================
470 f=0.1
471 f2=0.3
472 h=self.__h
473 fs=h.getFunctionSpace()
474 vol=self.getVolumeOfNegativeDomain()
475 self.__reinitPde.setValue(D=1.0)
476 #============
477 grad_phi=grad(self.__phi,fs)
478 len_grad_phi=length(grad_phi)
479 w=grad_phi/len_grad_phi
480
481 # s=wherePositive(self. __makeInterface(2.))-whereNegative(self. __makeInterface(2.))
482 s=self. __makeInterface(2.)
483 # s2=whereNegative(len_grad_phi-f2)*(1-(1-len_grad_phi/f2)**2)+(1-whereNegative(len_grad_phi-f2))
484 s2=1.
485 s=s*s2
486 dtau=f*inf(h/abs(s))
487
488 #========================
489 diff=1.
490 d=1.
491 c =0
492 #==============================================
493 TVD=integrate(length(grad_phi))
494 print(("initial range ",inf(self.__phi),sup(self.__phi),"error:",Lsup(1.-len_grad_phi),"volume =",vol,TVD))
495 # saveVTK("test.%s.vtu"%c,l=length(grad(self.__phi,fs))-1,s=s,phi=self.__phi)
496
497 dtau=f*inf(h/abs(s))
498 while c < self.__reinit_max: # and abs(diff) >= 0.01:
499 #
500 grad_phi=grad(self.__phi,fs)
501 len_grad_phi=length(grad_phi)
502 #
503 # s=self.__makeInterface(1.)
504 # s2=whereNegative(len_grad_phi-f2)*(1-(1-len_grad_phi/f2)**2)+(1-whereNegative(len_grad_phi-f2))
505 # s=s*s2
506 # phi_on_h=interpolate(self.__phi,fs)
507 # self.__reinitPde.setValue(Y =self.__phi+dtau/2*s*(1.-len_grad_phi))
508 # phi_half=self.__reinitPde.getSolution()
509 # self.__reinitPde.setValue(Y =self.__phi+dtau*s*(1.-length(grad(phi_half,fs))))
510 # self.__reinitPde.setValue(Y =self.__phi, Y_reduced=dtau*s*(1.-inner(w,grad(phi_half,ReducedFunction(self.__domain)))))
511 # self.__reinitPde.setValue(Y =self.__phi+dtau*(s-inner(w,grad_phi)),X=dtau/2*h*(s-inner(w,grad_phi))*w)
512 # self.__reinitPde.setValue(Y =self.__phi+dtau*s*(1.-len_grad_phi),X=-dtau/2*h*s**2*grad_phi)
513 self.__reinitPde.setValue(Y=self.__phi+dtau*s*(1.-len_grad_phi),X=f*dtau/2*h*abs(s)*(1.-len_grad_phi)*grad_phi/len_grad_phi)
514
515 # self.__reinitPde.setValue(Y=self.__phi+dtau*s*(1.-len_grad_phi),X=f*dtau/2*h*abs(s)*(1.-len_grad_phi)*grad_phi/len_grad_phi)
516 self.__phi, previous = self.__reinitPde.getSolution(), self.__phi
517 # self.__updateInterface()
518
519 vol,vol_old=self.getVolumeOfNegativeDomain(),vol
520 diff=(vol-vol_old)
521 r=Lsup(length(grad(self.__phi))-1.)
522 TVD=integrate(length(grad(self.__phi,fs)))
523 print(("iteration :", c, "range ",inf(self.__phi),sup(self.__phi),"error :",r,"volume change:",diff,TVD))
524 # saveVTK("test.%s.vtu"%(c+1),l=length(grad(self.__phi,fs)),s=s,phi=self.__phi,v=grad(self.__phi,fs))
525 c += 1
526 return
527 #==============================================
528 f2=0.7
529 h=self.__h
530 fs=h.getFunctionSpace()
531 vol=self.getVolumeOfNegativeDomain()
532 s=abs(self. __makeInterface(2.))
533 grad_phi=grad(self.__phi,fs)
534 len_grad_phi=length(grad_phi)
535
536 #----
537 # aphi_on_h=abs(interpolate(self.__phi,self.__h.getFunctionSpace()))
538 # q=2*self.__h
539 # mask_neg = whereNegative(aphi_on_h-q)
540 # mask_pos = wherePositive(aphi_on_h-q-self.__h)
541 # mask_interface = 1.-mask_neg-mask_pos
542 # s2=mask_neg + mask_interface * (q+self.__h-aphi_on_h)/self.__h
543 #----
544 m=whereNonPositive(len_grad_phi-f2)
545 s2=m*len_grad_phi/f2+(1.-m)
546 #----
547 # s2=1.
548 #----
549 self.__reinitPde.setValue(D=(1-s)/h**2,Y=(1-s)/h**2*self.__phi)
550 self.__reinitPde.setValue(A=(s2*len_grad_phi+(1.-s2))*kronecker(3),X=grad_phi)
551 # self.__reinitPde.setValue(A=(len_grad_phi-1)/len_grad_phi*kronecker(3))
552 # self.__reinitPde.setValue(A=kronecker(3), X=grad_phi/len_grad_phi)
553 self.__phi = self.__reinitPde.getSolution()
554 r=Lsup(length(grad(self.__phi))-1.)
555 vol,vol_old=self.getVolumeOfNegativeDomain(),vol
556 diff=(vol-vol_old)/vol
557 print(("iteration :", inf(self.__phi),sup(self.__phi),r,diff))
558 # saveVTK("test.%s.vtu"%0,l=length(grad(self.__phi,fs)),s=s,phi=self.__phi,v=grad(self.__phi,fs),s2=s2)
559 return
560 #=============================================
561
562 #============
563
564
565 def getVolumeOfNegativeDomain(self):
566 """
567 Returns the current volume of domain with phi<0.
568 """
569 return integrate((1.-self.__makeInterface(1.))/2.)
570
571
572 def getBoundingBoxOfNegativeDomain(self):
573 """
574 Returns the height of the region with phi<0.
575 """
576 fs=self.__h.getFunctionSpace()
577 mask_phi1=wherePositive(interpolate(self.__phi,fs))
578 mask_phi2=wherePositive(self.__phi)
579 x1=fs.getX()
580 x2=self.__domain.getX()
581 out=[]
582 for i in range(fs.getDim()):
583 x1_i=x1[i]
584 x2_i=x2[i]
585 d=2*(sup(x2_i)-inf(x2_i))
586 offset1=d*mask_phi1
587 offset2=d*mask_phi2
588 out.append((min(inf(x1_i+offset1),inf(x2_i+offset2)),max(sup(x1_i-offset1),sup(x2_i-offset2))))
589 return tuple(out)
590

  ViewVC Help
Powered by ViewVC 1.1.26