Changeset 7581


Ignore:
Timestamp:
Dec 7, 2009, 9:27:44 PM (14 years ago)
Author:
ole
Message:

Further simplifications

File:
1 edited

Legend:

Unmodified
Added
Removed
  • anuga_work/development/classroom/ripcurrent.py

    r7580 r7581  
    1 """Simple water flow example using ANUGA
     1"""ANUGA simulation of simple rip current.
     2
     3Source: Reference to paper in here
    24"""
    35
     
    1820# Parameters
    1921#------------------------------------------------------------------------------
     22
    2023filename = "WORKING-RIP-LAB_Expt-Geometry_Triangular_Mesh"
     24
    2125location_of_shore = 140 # The position along the y axis of the shorefront
    2226sandbar = 1.2           # Height of sandbar
     
    6064                                    verbose = True)
    6165                                   
    62 domain.set_name(filename)                  # Output name
     66domain.set_name(filename) # Output name
    6367print domain.statistics()
    6468
     
    7175
    7276    # General slope and buildings
    73     z=0.05*(y-(location_of_shore))
     77    z=0.05*(y-location_of_shore)
    7478   
    7579    N = len(x)
    7680    for i in range(N):
    7781        if y[i] < 25:
    78             z[i] = (0.2*(y[i]-25)) + 0.05*(y[i]-(location_of_shore))
     82            z[i] = 0.2*(y[i]-25) + 0.05*(y[i]-location_of_shore)
    7983    for i in range(N):
    8084        if y[i]>150:
    81             z[i] = (0.1*(y[i]-150)) + 0.05*(y[i]-(location_of_shore))
     85            z[i] = 0.1*(y[i]-150) + 0.05*(y[i]-location_of_shore)
    8286
    8387    return z
     
    9195    # It would be great with a comment about what this does
    9296    for i in range(N):
    93         ymin = -1*(bank_slope)*x[i] + 112
    94         ymax = -1*(bank_slope)*x[i] + 124
     97        ymin = -bank_slope*x[i] + 112
     98        ymax = -bank_slope*x[i] + 124
    9599        xmin = 0
    96         xmax = (length/2)-halfchannelwidth
     100        xmax = length/2-halfchannelwidth
    97101        if ymin < y[i] < ymax and xmin < x[i]< xmax:
    98             z[i] += sandbar*((cos((y[i]-118)/steepness)))
     102            z[i] += sandbar*cos((y[i]-118)/steepness)
    99103   
    100104    # It would be great with a comment about what this does       
    101105    for i in range(N):
    102         ymin = -1*(bank_slope)*(x[i]-(length/2)) + (-1*(bank_slope)*(length/2)+112)
    103         ymax = -1*(bank_slope)*(x[i]-(length/2)) + (-1*(bank_slope)*(length/2)+124)
    104         xmin = (length/2)+halfchannelwidth
     106        ymin = -bank_slope*(x[i]-length/2) - bank_slope*length/2 + 112
     107        ymax = -bank_slope*(x[i]-length/2) - bank_slope*length/2 + 124
     108        xmin = length/2+halfchannelwidth
    105109        xmax = 183
    106110        if ymin < y[i] < ymax and xmin < x[i] < xmax:
    107             z[i] += sandbar*(cos((y[i]-118)/steepness))
     111            z[i] += sandbar*cos((y[i]-118)/steepness)
    108112           
    109113    return z
     
    151155gauges = numpy.array(G)
    152156number_of_gauges = len(gauges)
    153 print gauges                               
    154157
    155158# Allocate space for velocity values
     
    168171    u += uh/depth
    169172    v += vh/depth
    170 print 'Evolution took %.2f seconds' % (time.time()-t0)
    171173
    172174   
     
    215217v_average = v/n_time_intervals
    216218
    217 #print "there were", n_time_intervals, "time steps"
     219print 'There were %i time steps' % n_time_intervals
    218220
    219221#print "sum y velocity", v
    220 #print "average y velocity", v_average
     222print "average y velocity", v_average
    221223#print "sum x velocity", u
    222224#print "average x velocity", u_average
    223 
    224 x_output = file('x_velocity.txt', 'w')
    225 y_output = file('y_velocity.txt', 'w')
    226 
    227 #print >> x_output, " "
    228 #print >> y_output, " "
    229 
    230 #print >> x_output, u_average
    231 #print >> y_output, v_average
    232225
    233226
     
    235228Y = gauges[:,1]
    236229                             
    237 U = u_average.tolist()
    238 V = v_average.tolist()
     230U = u_average #.tolist()
     231V = v_average #.tolist()
    239232
    240233#print "U = ", U
    241234#print "U has type", type(U)
     235
     236print 'Computation took %.2f seconds' % (time.time()-t0)
    242237
    243238
Note: See TracChangeset for help on using the changeset viewer.