Merge pull request #31 from danielru/add/images_boussinesq

Add/images boussinesq

Author: pancetta

examples/boussinesq_2d_imex/HookClass.py

@@ -20,8 +20,8 @@ def __init__(self):
 
         # add figure object for further use
         #        self.fig = plt.figure(figsize=(18,6))
-        self.fig = plt.figure(figsize=(9,9))
-
+        self.fig = plt.figure(figsize=(15,5))
+        self.counter = 0
 
     def dump_step(self,status):
         """
@@ -40,23 +40,34 @@ def dump_step(self,status):
         #plt.ylim([-1.1, 1.1])
         #plt.show(block=False)
         #plt.pause(0.00001)        
-          
-        if True:
+        
+        fs = 18
+        
+        if False:
           yplot = self.level.uend.values
           xx    = self.level.prob.xx
           zz    = self.level.prob.zz
-          self.fig.clear()
+          #self.fig.clear()
           levels = [-5e-3, -4e-3, -3e-3, -2e-3, -1e3, 1e-3, 2e-3, 3e-3, 4e-3, 5e-3]
           CS = plt.contourf(xx, zz, yplot[2,:,:], rstride=1, cstride=1, cmap=cm.coolwarm, linewidth=0, antialiased=False)
 #          CS = plt.contour(xx, zz, yplot[2,:,:], rstride=1, cstride=1, linewidth=0, antialiased=False)
-          cbar = plt.colorbar(CS)
+#          cbar = plt.colorbar(CS)
           plt.axes().set_xlim(xmin = self.level.prob.x_bounds[0], xmax = self.level.prob.x_bounds[1])
           plt.axes().set_ylim(ymin = self.level.prob.z_bounds[0], ymax = self.level.prob.z_bounds[1])
-          #plt.axes().set_aspect('equal')
-          plt.xlabel('x')
-          plt.ylabel('z')
+          plt.axes().set_aspect('auto')
+          plt.xlabel('x in km', fontsize=fs)
+          plt.ylabel('z in km', fontsize=fs)
+          plt.tick_params(axis='both', which='major', labelsize=fs)
           #plt.tight_layout()
-          plt.show(block=False)
-          plt.pause(0.00001)
-
+          #plt.show(block=False)
+          #plt.show()
+          plt.draw()
+          plt.gcf().savefig('images/boussinesq-'+"%04d" % self.counter +'.png', bbox_inches='tight')
+          self.counter = self.counter + 1
+          plt.pause(0.001)
+        
+          # NOTE: Can use ffmpeg to collate images into movie.
+          # HELPFUL WEBSITE: http://hamelot.co.uk/visualization/using-ffmpeg-to-convert-a-set-of-images-into-a-video/
+          # USEFUL COMMAND: ffmpeg -r 25 -i boussinesq-%04d.jpeg -vcodec libx264 -crf 25 test1800.avi
+          # WEBSITE TO CONVERT TO MP4: http://video.online-convert.com/convert-to-mp4
         return None

examples/boussinesq_2d_imex/playground.py

@@ -28,34 +28,36 @@
 
     # This comes as read-in for the level class
     lparams = {}
-    lparams['restol'] = 1E-12
+    lparams['restol'] = 1E-14
 
     swparams = {}
     swparams['collocation_class'] = collclass.CollGaussLobatto
     swparams['num_nodes'] = 3
     swparams['do_LU'] = True
 
     sparams = {}
-    sparams['maxiter'] = 15
+    sparams['maxiter'] = 4
 
     # setup parameters "in time"
     t0     = 0
-    Tend   = 3
-    Nsteps =  1
+    Tend   = 3000
+    Nsteps =  500
+    Tend   = 30
+    Nsteps =  5
     dt = Tend/float(Nsteps)
 
     # This comes as read-in for the problem class
     pparams = {}
-    pparams['nvars']    = [(4,300,10)]
+    pparams['nvars']    = [(4,450,30)]
     pparams['u_adv']    = 0.02
     pparams['c_s']      = 0.3
     pparams['Nfreq']    = 0.01
     pparams['x_bounds'] = [(-150.0, 150.0)]
     pparams['z_bounds'] = [(   0.0,  10.0)]
-    pparams['order']    = [0, 0] # [fine_level, coarse_level]
-    pparams['gmres_maxiter'] = [50, 5]
+    pparams['order']    = [0] # [fine_level, coarse_level]
+    pparams['gmres_maxiter'] = [50]
     pparams['gmres_restart'] = 20
-    pparams['gmres_tol']     = 1e-14
+    pparams['gmres_tol']     = 1e-6
 
     # This comes as read-in for the transfer operations
     tparams = {}
@@ -71,15 +73,22 @@
     description['sweeper_class']     = imex_1st_order
     description['level_params']      = lparams
     description['hook_class']        = plot_solution
-    description['transfer_class']    = mesh_to_mesh_2d
-    description['transfer_params']   = tparams
+    #description['transfer_class']    = mesh_to_mesh_2d
+    #description['transfer_params']   = tparams
 
     # quickly generate block of steps
     MS = mp.generate_steps(num_procs,sparams,description)
 
     # get initial values on finest level
     P = MS[0].levels[0].prob
     uinit = P.u_exact(t0)
+    
+    cfl_advection = pparams['u_adv']*dt/P.h[0]
+    cfl_acoustic_hor = pparams['c_s']*dt/P.h[0]
+    cfl_acoustic_ver = pparams['c_s']*dt/P.h[1]
+    print ("CFL number of advection: %4.2f" % cfl_advection)
+    print ("CFL number of acoustics (horizontal): %4.2f" % cfl_acoustic_hor)
+    print ("CFL number of acoustics (vertical):   %4.2f" % cfl_acoustic_ver)
 
     # call main function to get things done...
     uend,stats = mp.run_pfasst(MS,u0=uinit,t0=t0,dt=dt,Tend=Tend)
@@ -92,8 +101,8 @@
 
     P.report_log()
 
-    plt.show()
+    #plt.show()
 
-    # extract_stats = grep_stats(stats,iter=-1,type='residual')
-    # sortedlist_stats = sort_stats(extract_stats,sortby='step')
-    # print(extract_stats,sortedlist_stats)
+    #extract_stats = grep_stats(stats,iter=-1,type='residual')
+    #sortedlist_stats = sort_stats(extract_stats,sortby='step')
+    #print(extract_stats,sortedlist_stats)