# Diff of /trunk/dudley/test/python/FCT_benchmark.py

revision 3891 by jfenwick, Wed Jan 18 06:29:34 2012 UTC revision 3892 by jfenwick, Tue Apr 10 08:57:23 2012 UTC
33  #  #
34  #  the solution is given as   u(x,t)=U0*s^dim/(s**2+E*t)^{dim/2} * exp ( - |x-x_0(t)|^2/(4*(s**2+E*t)) )  #  the solution is given as   u(x,t)=U0*s^dim/(s**2+E*t)^{dim/2} * exp ( - |x-x_0(t)|^2/(4*(s**2+E*t)) )
35  #  #
36  #                 with x_0(t) = X0 + (v+w)*t  #                             with x_0(t) = X0 + (v+w)*t
37  #  #
38  #  #
39  #    the region |x-x_0(t)|^2/(4*(s**2+E*t)) < - log(TAU)  is within the domain for all time  #    the region |x-x_0(t)|^2/(4*(s**2+E*t)) < - log(TAU)  is within the domain for all time
41  #  #
42  #    this holds if  #    this holds if
43  #  #
44  #       |x_i-X0_i-(v_i+w_i)*tend | < sqrt(- log(TAU) * 4*(s**2+E*tend))=b0  and  #       |x_i-X0_i-(v_i+w_i)*tend | < sqrt(- log(TAU) * 4*(s**2+E*tend))=b0  and
45  #       |x_i-X0_i | < sqrt(- log(TAU)) * 2*s = b1 implies 0<=x_i<=l_i  #       |x_i-X0_i | < sqrt(- log(TAU)) * 2*s = b1 implies 0<=x_i<=l_i
46  #  #
47  from math import pi, ceil  from math import pi, ceil
48  from time import time as clock  from time import time as clock
# Line 149  def XXX(dim,tend,dt, s, h,b,c,d,c_dir="x Line 149  def XXX(dim,tend,dt, s, h,b,c,d,c_dir="x
149      if dim == 2:      if dim == 2:
150          if VERBOSITY: print("%d x %d grid over %e  x %e with element size %e."%(NE_0,NE_1,l_0,l_1,h))          if VERBOSITY: print("%d x %d grid over %e  x %e with element size %e."%(NE_0,NE_1,l_0,l_1,h))
151          if NE_0*NE_1 > NE_MAX:          if NE_0*NE_1 > NE_MAX:
152         raise ValueError("too many elements %s."%(NE_0*NE_1,))             raise ValueError("too many elements %s."%(NE_0*NE_1,))
153          dom=Rectangle(n0=NE_0,n1=NE_1,l0=l_0,l1=l_1)          dom=Rectangle(n0=NE_0,n1=NE_1,l0=l_0,l1=l_1)
154          x0=[X0_0, X0_1]          x0=[X0_0, X0_1]
155      else:      else:
# Line 158  def XXX(dim,tend,dt, s, h,b,c,d,c_dir="x Line 158  def XXX(dim,tend,dt, s, h,b,c,d,c_dir="x
158         NE_2=max(int(l_2/h+0.5),1)         NE_2=max(int(l_2/h+0.5),1)
159         if VERBOSITY: print("%d x %d x %d grid over %e  x %e x %e with element size %e."%(NE_0,NE_1,NE_2,l_0,l_1,l_2,h))         if VERBOSITY: print("%d x %d x %d grid over %e  x %e x %e with element size %e."%(NE_0,NE_1,NE_2,l_0,l_1,l_2,h))
160         if NE_0*NE_1*NE_2 > NE_MAX:         if NE_0*NE_1*NE_2 > NE_MAX:
161        raise ValueError("too many elements %s."%(NE_0*NE_1*NE_2,))            raise ValueError("too many elements %s."%(NE_0*NE_1*NE_2,))
162         dom=Brick(n0=NE_0,n1=NE_1, ne2=NE_2, l0=l_0,l1=l_1, l2=l_2)         dom=Brick(n0=NE_0,n1=NE_1, ne2=NE_2, l0=l_0,l1=l_1, l2=l_2)
163         x0=[X0_0, X0_1, X0_2]         x0=[X0_0, X0_1, X0_2]
164      if VERBOSITY:      if VERBOSITY:

Legend:
 Removed from v.3891 changed lines Added in v.3892