/[escript]/branches/prepy3/doc/examples/cookbook/example09c.py
ViewVC logotype

Diff of /branches/prepy3/doc/examples/cookbook/example09c.py

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 3662 by jfenwick, Tue Nov 8 04:24:06 2011 UTC revision 3663 by jfenwick, Tue Nov 15 01:12:13 2011 UTC
# Line 45  from esys.finley import ReadMesh Line 45  from esys.finley import ReadMesh
45  ########################################################MPI WORLD CHECK  ########################################################MPI WORLD CHECK
46  if getMPISizeWorld() > 1:  if getMPISizeWorld() > 1:
47      import sys      import sys
48      print "This example will not run in an MPI world."      print("This example will not run in an MPI world.")
49      sys.exit(0)      sys.exit(0)
50    
51  #################################################ESTABLISHING VARIABLES  #################################################ESTABLISHING VARIABLES
# Line 68  vel=Scalar(0,Function(domain)) Line 68  vel=Scalar(0,Function(domain))
68  mu=Scalar(0,Function(domain))  mu=Scalar(0,Function(domain))
69  lam=Scalar(0,Function(domain))  lam=Scalar(0,Function(domain))
70    
71  print 0.5*np.sqrt(prho/(plam+2*pmu))*0.5  print(0.5*np.sqrt(prho/(plam+2*pmu))*0.5)
72    
73  for i in range(0,nlayers):  for i in range(0,nlayers):
74      rho.setTaggedValue('lblock%d'%i,prho[i])      rho.setTaggedValue('lblock%d'%i,prho[i])
# Line 89  lam.setTaggedValue('fault',plam[i]) Line 89  lam.setTaggedValue('fault',plam[i])
89  # Time related variables.  # Time related variables.
90  testing=False  testing=False
91  if testing:  if testing:
92      print 'The testing end time is currently selected. This severely limits the number of time iterations.'      print('The testing end time is currently selected. This severely limits the number of time iterations.')
93      print "Try changing testing to False for more iterations."      print("Try changing testing to False for more iterations.")
94      tend=0.1      tend=0.1
95  else:  else:
96      tend=0.1    # end time      tend=0.1    # end time
# Line 100  h=0.00001    # time step Line 100  h=0.00001    # time step
100  rtime=0.0 # first time to record  rtime=0.0 # first time to record
101  rtime_inc=tend/750.0 # time increment to record  rtime_inc=tend/750.0 # time increment to record
102  #Check to make sure number of time steps is not too large.  #Check to make sure number of time steps is not too large.
103  print "Time step size= ",h, "Expected number of outputs= ",tend/h  print("Time step size= ",h, "Expected number of outputs= ",tend/h)
104    
105  U0=0.1 # amplitude of point source  U0=0.1 # amplitude of point source
106  ls=500   # length of the source  ls=500   # length of the source
# Line 115  a = 2.0 * (np.pi * dfeq)**2.0 Line 115  a = 2.0 * (np.pi * dfeq)**2.0
115  t0 = 5.0 / (2.0 * np.pi * dfeq)  t0 = 5.0 / (2.0 * np.pi * dfeq)
116  srclength = 5. * t0  srclength = 5. * t0
117  ls = int(srclength/h)  ls = int(srclength/h)
118  print 'source length',ls  print('source length',ls)
119  source=np.zeros(ls,'float') # source array  source=np.zeros(ls,'float') # source array
120  ampmax=0  ampmax=0
121  for it in range(0,ls):  for it in range(0,ls):
# Line 140  mypde.setValue(D=rho*kmat) #set the gene Line 140  mypde.setValue(D=rho*kmat) #set the gene
140    
141  ############################################FIRST TIME STEPS AND SOURCE  ############################################FIRST TIME STEPS AND SOURCE
142  # define small radius around point xc  # define small radius around point xc
143  src_length = 10; print "src_length = ",src_length  src_length = 10; print("src_length = ",src_length)
144  # set initial values for first two time steps with source terms  # set initial values for first two time steps with source terms
145  xb=FunctionOnBoundary(domain).getX()  xb=FunctionOnBoundary(domain).getX()
146  yx=(cos(length(xb-xc)*3.1415/src_length)+1)*whereNegative(length(xb-xc)-src_length)  yx=(cos(length(xb-xc)*3.1415/src_length)+1)*whereNegative(length(xb-xc)-src_length)
# Line 176  while t<tend: Line 176  while t<tend:
176      t=t+h; n=n+1      t=t+h; n=n+1
177      if (n < ls):      if (n < ls):
178          mypde.setValue(y=source[n]*yx*src_dir*stop) #set the source as a function on the boundary          mypde.setValue(y=source[n]*yx*src_dir*stop) #set the source as a function on the boundary
179      print n,"-th time step t ",t      print(n,"-th time step t ",t)

Legend:
Removed from v.3662  
changed lines
  Added in v.3663

  ViewVC Help
Powered by ViewVC 1.1.26