Changeset 2448

Timestamp:

Feb 24, 2006, 6:22:02 PM (18 years ago)

Author:

ole

Message:

Created stochastic version of okushiri with multiple realisations

Location:

development/stochastic_study

Files:

• bathymetry.pts (added)
• bathymetry.txt (deleted)
• project.py (modified) (1 diff)
• run_model.py (modified) (4 diffs)

development/stochastic_study/project.py

@@ -4 +4 @@
 # Inputs
 boundary_filename = 'input_wave.tms'
-bathymetry_filename = 'bathymetry.txt'
+bathymetry_filename = 'bathymetry.pts'
 
 # Model name
 basename = 'simulation'
 mesh_filename = basename + '.msh'
+
+# Stats (Suresh ?)
+number_of_samples = 300
+std_dev = 0.01
+mean = 0.0
+blocksize = 100 #How many realisations to fit at a time
+
+
+
 
 

development/stochastic_study/run_model.py

@@ -22 +22 @@
 import os
 import time
+import cPickle
 
 # Related major packages
@@ -33 +34 @@
 # Application specific imports
 import project                 # Definition of file names and polygons
+
+
 
 
@@ -47 +50 @@
 print 'The extent is ', domain.get_extent()
 
-domain.set_name(project.basename)
 domain.set_datadir('.')
 domain.set_quantities_to_be_stored(['stage', 'xmomentum', 'ymomentum'])
 
-domain.check_integrity()
-
-
-#------------------------------------------------------------------------------
-# Setup initial conditions
-#------------------------------------------------------------------------------
-domain.set_quantity('elevation',
-                    filename = project.bathymetry_filename[:-4] + '.pts',
-                    alpha = 0.001,
-                    verbose = True,
-                    use_cache = True)
-
-domain.set_quantity('friction', 0.0)
-domain.set_quantity('stage', 0.0)
+#domain.check_integrity()
 
 
@@ -80 +69 @@
 
 #------------------------------------------------------------------------------
-# Evolve system through time
+# Setup initial conditions
 #------------------------------------------------------------------------------
-t0 = time.time()
-for t in domain.evolve(yieldstep = 0.05, finaltime = 22.5):
-    domain.write_time()
-    domain.write_boundary_statistics(tags = 'wave')         
+domain.set_quantity('friction', 0.0)
+domain.set_quantity('stage', 0.0)
 
-print 'That took %.2f seconds' %(time.time()-t0)
+# Get prefitted realisations
+realisation = 0
+for filename in os.listdir('.'):
+    if filename.startswith(project.basename) and filename.endswith('.pck'):
+        print 'Reading %s' %filename
+        fid = open(filename)
+        V = cPickle.load(fid)
+        fid.close()
+
+        # For each column (each realisation)
+        for i in range(V.shape[1]):
+            name = project.basename + '_realisation_%d' %realisation
+            domain.set_name(name)                     #Output name
+            domain.set_quantity('elevation', V[:,i])  #Assign bathymetry
+
+            #---------------------------------------------------
+            # Evolve system through time
+            #---------------------------------------------------
+            print 'Running realisation %d of %d in block %s'\
+                  %(i, V.shape[1], filename)
+            t0 = time.time()
+            for t in domain.evolve(yieldstep = 0.05, finaltime = 22.5):
+                domain.write_time()
+                #domain.write_boundary_statistics(tags = 'wave')         
+        
+            print 'Realisation %d took %.2f seconds'\
+                  %(realisation, time.time()-t0)
+
+            realisation += 1