| 1 | """View results of runup.py |
|---|
| 2 | """ |
|---|
| 3 | #--------------- |
|---|
| 4 | # Import Modules |
|---|
| 5 | #--------------- |
|---|
| 6 | import anuga |
|---|
| 7 | #import struct |
|---|
| 8 | import numpy |
|---|
| 9 | #import scipy |
|---|
| 10 | from matplotlib import pyplot as pyplot |
|---|
| 11 | import util |
|---|
| 12 | |
|---|
| 13 | p=util.get_output('parabola_v2.sww', 0.01) |
|---|
| 14 | p2=util.get_centroids(p, velocity_extrapolation=True) |
|---|
| 15 | |
|---|
| 16 | # Define some 'lines' along which to plot |
|---|
| 17 | v=(p2.y==p2.y[0]) |
|---|
| 18 | |
|---|
| 19 | #-------------------- |
|---|
| 20 | # Make plot animation |
|---|
| 21 | #-------------------- |
|---|
| 22 | pyplot.close() #If the plot is open, there will be problems |
|---|
| 23 | pyplot.ion() |
|---|
| 24 | |
|---|
| 25 | if False: |
|---|
| 26 | line, = pyplot.plot( (p2.x[v].min(),p2.x[v].max()) ,(p2.stage[:,v].min(),p2.stage[:,v].max() ) ) |
|---|
| 27 | for i in range(700,p2.xmom.shape[0]): |
|---|
| 28 | line.set_xdata(p2.x[v]) |
|---|
| 29 | line.set_ydata(p2.stage[i,v]) |
|---|
| 30 | pyplot.draw() |
|---|
| 31 | #pyplot.plot( (0,1),(0,0), 'r' ) |
|---|
| 32 | pyplot.plot(p2.x[v],p2.elev[v]) |
|---|
| 33 | pyplot.title(str(i)+'/200') # : velocity does not converge to zero' ) |
|---|
| 34 | pyplot.xlabel('x') |
|---|
| 35 | pyplot.ylabel('stage (m/s)') |
|---|
| 36 | |
|---|
| 37 | if False: |
|---|
| 38 | pyplot.clf() |
|---|
| 39 | |
|---|
| 40 | line, = pyplot.plot( (p2.x[v].min(),p2.x[v].max()) ,(p2.xvel[:,v].min(),p2.xvel[:,v].max() ), 'ro' ) |
|---|
| 41 | for i in range(700,p2.xmom.shape[0]): |
|---|
| 42 | line.set_xdata(p2.x[v]) |
|---|
| 43 | line.set_ydata(p2.xvel[i,v]) |
|---|
| 44 | pyplot.draw() |
|---|
| 45 | pyplot.plot( (p2.x[v].min(),p2.x[v].max()),(0,0), 'ro' ) |
|---|
| 46 | #pyplot.plot(x[v],elev[v]) |
|---|
| 47 | pyplot.title(str(i)+'/200') # : velocity does not converge to zero' ) |
|---|
| 48 | pyplot.xlabel('x') |
|---|
| 49 | pyplot.ylabel('xvel (m/s)') |
|---|
| 50 | |
|---|
| 51 | D0=4.0 |
|---|
| 52 | L=10.0 |
|---|
| 53 | A=2.0 |
|---|
| 54 | g=9.8 |
|---|
| 55 | #t=time[30] |
|---|
| 56 | |
|---|
| 57 | omega=numpy.sqrt(2*D0*g)/L |
|---|
| 58 | T= 2*numpy.pi/omega |
|---|
| 59 | ppp= (abs(p2.x-4.*L/2.)).argmin() |
|---|
| 60 | ppp2= (abs(p2.x-4.*L/4.)).argmin() |
|---|
| 61 | |
|---|
| 62 | w = D0 + 2*A*D0/(L**2)*numpy.cos(omega*p2.time)*( (p2.x[ppp]-4.*L/2.) -A/2.*numpy.cos(omega*p2.time)) |
|---|
| 63 | w2 = D0 + 2*A*D0/(L**2)*numpy.cos(omega*p2.time)*( (p2.x[ppp2]-4.*L/2.) -A/2.*numpy.cos(omega*p2.time)) |
|---|
| 64 | w2 = w2*(w2>p2.elev[ppp2])+p2.elev[ppp2]*(w2<=p2.elev[ppp2]) |
|---|
| 65 | |
|---|
| 66 | pyplot.clf() |
|---|
| 67 | pyplot.plot(p2.time,w, color='blue') |
|---|
| 68 | pyplot.plot(p2.time,p2.stage[:,ppp], color='green') |
|---|
| 69 | pyplot.savefig('Stage_centre_v2.png') |
|---|
| 70 | # pyplot.savefig('runup_x_velocities.png') |
|---|
| 71 | pyplot.clf() |
|---|
| 72 | pyplot.plot(p2.time,w2, color='blue') |
|---|
| 73 | pyplot.plot(p2.time,p2.stage[:,ppp2], color='green') |
|---|
| 74 | pyplot.savefig('Stage_centre_v3.png') |
|---|
| 75 | |
|---|
| 76 | |
|---|
| 77 | pltind=numpy.argmin(abs(p2.time-T*3)) |
|---|
| 78 | |
|---|
| 79 | pyplot.clf() |
|---|
| 80 | pyplot.plot(p2.x[v],p2.xvel[pltind,v]) |
|---|
| 81 | pyplot.title('Velocity near to time=3T') |
|---|
| 82 | pyplot.savefig('Vel_3T_v2.png') |
|---|
| 83 | |
|---|
| 84 | pltind=numpy.argmin(abs(p2.time-T*3.25)) |
|---|
| 85 | |
|---|
| 86 | pyplot.clf() |
|---|
| 87 | pyplot.plot(p2.x[v],p2.xvel[pltind,v]) |
|---|
| 88 | pyplot.title('Velocity near to time=3.25T') |
|---|
| 89 | pyplot.savefig('Vel_3_5T_v2.png') |
|---|
| 90 | #pyplot.clf() |
|---|
| 91 | #pyplot.close() |
|---|
