# Changeset 7842 for anuga_work/development/2010-projects/anuga_1d/channel/test_padarn.py

Ignore:
Timestamp:
Jun 15, 2010, 5:34:24 PM (12 years ago)
Message:

Adding in a few new files

File:
1 edited

Unmodified
Added
Removed
• ## anuga_work/development/2010-projects/anuga_1d/channel/test_padarn.py

 r7839 import random from math import sqrt, pow, pi ,sin, cos from channel_domain import * from Numeric import allclose, array, zeros, ones, Float, take, sqrt from config import g, epsilon import numpy from anuga_1d.channel.channel_domain import * from anuga_1d.config import g, epsilon omega = sqrt(2*g*z_infty)/L_x ## angular frequency of osccilation t=0.0 y = zeros(len(x),Float) y = numpy.zeros(len(x),numpy.float) for i in range(len(x)): #y[i] = z_infty+2*A0*z_infty/L_x*cos(omega*t)*(x[i]/L_x-0.5*A0/(L_x)*cos(omega*t)) N = len(x) z_infty = 10.0 z = zeros(N,Float) z = numpy.zeros(N,numpy.float) L_x = 2500.0 A0 = 0.5*L_x def initial_area(x): y = zeros(len(x),Float) y = numpy.zeros(len(x),numpy.float) for i in range(len(x)): y[i]=initial_stage([x[i]])-bed([x[i]]) print "Evaluating domain with %d cells" %N cell_len = 4*L/N # Origin = 0.0 points = zeros(N+1,Float) points = numpy.zeros(N+1,numpy.float) for i in range(N+1): points[i] = -2*L +i*cell_len N = float(N) HeightC = domain.quantities['height'].centroid_values HeightC    = domain.quantities['height'].centroid_values DischargeC = domain.quantities['discharge'].centroid_values C = domain.centroids C          = domain.centroids print 'That took %.2f seconds' %(time.time()-t0) X = domain.vertices HeightQ = domain.quantities['area'].vertex_values VelocityQ = domain.quantities['velocity'].vertex_values x = X.flat z = domain.quantities['elevation'].vertex_values.flat stage=HeightQ.flat+z X          = domain.vertices HeightQ    = domain.quantities['area'].vertex_values VelocityQ  = domain.quantities['velocity'].vertex_values Z          = domain.quantities['elevation'].vertex_values Stage      = HeightQ + Z plot1 = subplot(211) plot(x,z,x,stage) plot(X.flat,Z.flat, X.flat,Stage.flat) plot1.set_ylim([-1,35]) 'upper right', shadow=True) plot2 = subplot(212) plot(x,VelocityQ.flat) plot(X.flat,VelocityQ.flat) plot2.set_ylim([-10,10])
Note: See TracChangeset for help on using the changeset viewer.