Changeset 2966


Ignore:
Timestamp:
May 25, 2006, 11:45:56 AM (18 years ago)
Author:
ole
Message:

Work on stochastic study

Location:
development/stochastic_study
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • development/stochastic_study/create_realisations.py

    r2959 r2966  
    77"""
    88
    9 from RandomArray import normal
     9from RandomArray import normal, seed
    1010from Numeric import concatenate, NewAxis
    1111import cPickle
     
    3333elevation = base_bathymetry['attributelist']['elevation']
    3434N = len(elevation)
    35 print 'read', N, 'points'
     35print 'read', N, 'points:' + 'z in [%f, %f]' %(min(elevation),
     36                                               max(elevation))
    3637
    3738# Create different bathymetries based on sampling (Suresh to fill this in)
     
    3940Z = None
    4041block_number = 0
    41 for n in range(project.number_of_samples):
     42seed(13,17)
     43for n in range(project.number_of_realisations):
    4244
    4345    print 'Creating sample #%d' %n
     
    4749                           N)
    4850
    49 
    5051    #Concatenate into array with each realisation a column
    5152    if Z is None:
     
    5455        Z = concatenate((Z, z[:, NewAxis]), axis=1)
    5556
    56     if (n+1)%project.blocksize == 0 or n == project.number_of_samples:
     57    if (n+1)%project.blocksize == 0 or n == project.number_of_realisations-1:
    5758        #Fit all of them to computational mesh
    5859        print 'fit block %d to the mesh' %block_number
     
    6768                        use_cache = False)
    6869
    69 
     70        print V.shape
     71        print V[:2,:]
    7072        Z = None
    7173
  • development/stochastic_study/project.py

    r2957 r2966  
    1515gauge_names = ['ch5', 'ch7', 'ch9']
    1616
     17# Constants
     18number_of_timesteps = 451 # Known from problem description
    1719
    1820# Stats (Suresh ?)
    19 number_of_samples = 300
    20 std_dev = 0.01
     21number_of_realisations = 4
     22#std_dev = 0.0026  #Range is 26.035 cm
     23std_dev = 0.0013  #Range is 26.035 cm
    2124mean = 0.0
    2225blocksize = 100 #How many realisations to fit at a time
  • development/stochastic_study/run_model.py

    r2959 r2966  
    4444print 'Creating domain from', project.mesh_filename
    4545
    46 #domain = pmesh_to_domain_instance(project.mesh_filename, Domain,
    47 #                                  use_cache=True,
    48 #                                  verbose=True)
    49 
    5046domain = Domain(project.mesh_filename,
    51                 use_cache=True,
     47                use_cache=False,
    5248                verbose=True)               
    5349               
    5450
    5551print 'Number of triangles = ', len(domain)
    56 print 'The extent is ', domain.get_extent()
    5752print domain.statistics()
    5853
     
    6055domain.set_datadir('.')
    6156domain.set_quantities_to_be_stored(['stage', 'xmomentum', 'ymomentum'])
    62 
    63 #domain.check_integrity()
    6457
    6558
     
    10093            domain.set_name(project.basename)         #Output name
    10194            domain.set_quantity('elevation', V[:,i])  #Assign bathymetry
    102             domain.starttime = 0.0                    #Reset time
     95            domain.set_time(0.0)                      #Reset time
    10396
    10497            #---------------------------------------------------
     
    107100            print 'Running realisation %d of %d in block %s'\
    108101                  %(i, V.shape[1], filename)
     102            print 'Z Range:', min(V[:,i]), max(V[:,i])
     103
    109104            t0 = time.time()
    110105            for t in domain.evolve(yieldstep = timestep,
     
    113108               
    114109       
    115             print 'Realisation %d took %.2f seconds'\
     110            print 'Simulation of realisation %d took %.2f seconds'\
    116111                  %(realisation, time.time()-t0)
    117112
Note: See TracChangeset for help on using the changeset viewer.