Updated documentation.

Added OpenMP, performance tips, how to use the variable cache.

doc/source/basic_usage.rst

@@ -1925,6 +1925,181 @@ Here is an example of the loop-and-fill technique:
        f.close()
 
       
+The *cache* Argument for :meth:`wrf.getvar`
+*********************************************     
+
+If you have read through the documentation, you may have noticed that the 
+:meth:`wrf.getvar` routine contains a *cache* argument. What is this for?
+
+Internally, if metadata is turned on, a variable is extracted from the NetCDF 
+file and its metadata is copied to form the result's metadata. Often this 
+variable is one of the computation's function arguments, so rather than 
+spend time extracting the variable again for the computation, it is placed 
+in a cache (dictionary) and passed on to the computational function.
+
+What isn't widely known is that this cache argument can also be supplied by 
+end users wishing to speed up their application. This can be useful in 
+situations where numerous calculations are being performed on the same 
+data set. For many algorithms, the cost to extract the arrays from the 
+NetCDF file is on par with the time to perform the calculation. If you are 
+computing numerous diagnostics, extracting the variables up front allows you 
+to only pay this extraction penalty once, rather than inside of each call 
+to :meth:`wrf.getvar`.
+
+The cache is nothing more than a dictionary where each key is the variable 
+name (e.g. "P") and the value is the :class:`xarray.DataArray` or 
+:class:`numpy.ndarray` variable. Creating the cache dictionary is easy, 
+since the :meth:`wrf.extract_vars` routine returns a dictionary for a 
+sequence of variables. 
 
+.. note:: 
+
+   The *timeidx* parameter supplied to :meth:`extract_vars` 
+   must be the same *timeidx* parameter that you plan to use for 
+   :meth:`wrf.getvar`. Otherwise, it will crash with dimension mismatch errors.
+
+Some common variables that you can use to create an effective cache are: P, PB, 
+PH, PHB, T, QVAPOR, HGT, PSFC, U, V, W.
+
+Below is an example showing the same computations done with and without the 
+cache. The execution time is printed. The hardware used is a 2.8 GHz Intel Core 
+i7, which contains 4 CPU cores with 2 hyper threads (8 total threads). This 
+will be interpreted as 8 CPUs for OpenMP.
+
+.. code:: python
+
+   from __future__ import print_function
+
+   import time
+   from netCDF4 import Dataset
+   from wrf import getvar, ALL_TIMES, extract_vars
+
+   # The first two files contain four times, the last file contains only one.
+   wrf_filenames = ["/path/to/wrfout_d02_2005-08-28_00:00:00",
+                    "/path/to/wrfout_d02_2005-08-28_12:00:00", 
+                    "/path/to/wrfout_d02_2005-08-29_00:00:00"]
+
+   wrfin = [Dataset(x) for x in wrf_filenames]
+
+   start = time.time()
+   my_cache = extract_vars(wrfin, ALL_TIMES, ("P", "PSFC", "PB", "PH", "PHB", 
+                                              "T", "QVAPOR", "HGT", "U", "V", 
+                                              "W"))
+   end = time.time()
+   print ("Time taken to build cache: ", (end-start), "s")
+    
+   vars = ("avo", "eth", "cape_2d", "cape_3d", "ctt", "dbz", "mdbz", 
+           "geopt", "helicity", "lat", "lon", "omg", "p", "pressure", 
+           "pvo", "pw", "rh2", "rh", "slp", "ter", "td2", "td", "tc", 
+           "theta", "tk", "tv", "twb", "updraft_helicity", "ua", "va", 
+           "wa", "uvmet10", "uvmet", "z", "cfrac", "zstag", "geopt_stag")
+   
+   # No cache
+   start = time.time()
+   for var in vars:
+       v = getvar(wrfin, var, ALL_TIMES)
+   end = time.time()
+   no_cache_time = (end-start)
+
+   print ("Time taken without variable cache: ", no_cache_time, "s")
+
+   # With a cache
+   start = time.time()
+   for var in vars:
+       v = getvar(wrfin, var, ALL_TIMES, cache=my_cache)
+   end = time.time()
+   cache_time = (end-start)
+
+   print ("Time taken with variable cache: ", cache_time, "s")
+
+   improvement = ((no_cache_time-cache_time)/no_cache_time) * 100 
+   print ("The cache decreased computation time by: ", improvement, "%")
+
+
+Result:
+
+.. code:: none
+
+   Time taken to build cache:  0.28154706955 s
+   Time taken without variable cache:  11.0905270576 s
+   Time taken with variable cache:  8.25931215286 s
+   The cache decreased computation time by:  25.5282268378 %
+   
+By removing the repeated extraction of common variables in the getvar routine,
+for the single threaded case, the computation time has been reduced by 
+25.5% in the particular example.
+
+Things get more interesting when OpenMP is turned on, and set to use the 
+maximum number of processors (in this case 8 threads are used).  
+
+.. code:: python
+
+   from __future__ import print_function
+
+   import time
+   from netCDF4 import Dataset
+   from wrf import (getvar, ALL_TIMES, extract_vars, 
+                    omp_set_num_threads, omp_get_num_procs)
+   
+   # The first two files contain four times, the last file contains only one.
+   wrf_filenames = ["/path/to/wrfout_d02_2005-08-28_00:00:00",
+                    "/path/to/wrfout_d02_2005-08-28_12:00:00", 
+                    "/path/to/wrfout_d02_2005-08-29_00:00:00"]
+
+   wrfin = [Dataset(x) for x in wrf_filenames]
+
+   start = time.time()
+   my_cache = extract_vars(wrfin, ALL_TIMES, ("P", "PSFC", "PB", "PH", "PHB", 
+                                              "T", "QVAPOR", "HGT", "U", "V", 
+                                              "W"))
+   end = time.time()
+   print ("Time taken to build cache: ", (end-start), "s")
+
+   omp_set_num_threads(omp_get_num_procs())
+   
+   vars = ("avo", "eth", "cape_2d", "cape_3d", "ctt", "dbz", "mdbz", 
+           "geopt", "helicity", "lat", "lon", "omg", "p", "pressure", 
+           "pvo", "pw", "rh2", "rh", "slp", "ter", "td2", "td", "tc", 
+           "theta", "tk", "tv", "twb", "updraft_helicity", "ua", "va", 
+           "wa", "uvmet10", "uvmet", "z", "cfrac", "zstag", "geopt_stag")
+   
+   # No cache
+   start = time.time()
+   for var in vars:
+       v = getvar(wrfin, var, ALL_TIMES)
+   end = time.time()
+   no_cache_time = (end-start)
+
+   print ("Time taken without variable cache: ", no_cache_time, "s")
+
+   # With a cache
+   start = time.time()
+   for var in vars:
+       v = getvar(wrfin, var, ALL_TIMES, cache=my_cache)
+   end = time.time()
+   cache_time = (end-start)
+
+   print ("Time taken with variable cache: ", cache_time, "s")
+
+   improvement = ((no_cache_time-cache_time)/no_cache_time) * 100 
+   print ("The cache decreased computation time by: ", improvement, "%")
+
+Result:
+
+.. code:: none
+
+   Time taken to build cache:  0.2700548172 s
+   Time taken without variable cache:  6.02652812004 s
+   Time taken with variable cache:  3.27777099609 s
+   The cache decreased computation time by:  45.6109565772 %
    
+In this example, 4 CPU cores (8 total threads) are used. When the cache is 
+used, the computation time drops by 45%, so almost half the time was spent 
+simply extracting variables from the NetCDF file. When compared to the 
+11.09 s needed to compute the single threaded case with no variable cache, the
+computation time drops by roughly 70% (compared to 45% with 8 threads but 
+no cache). 
 
+In summary, if you are computing a lot of diagnostic variables, consider using
+the *cache* argument to improve performance, particularly if you want to 
+maximize your multithreaded performance with OpenMP.

src/wrf/extension.py

@@ -976,7 +976,7 @@ def omp_get_num_threads():
 
 
 def omp_get_max_threads():
-    """Return the maximum number of threads that can be used in a parallel
+    """Return the maximum number of threads that can be used in a parallel \
     region.
     
     The omp_get_max_threads routine returns an upper bound on the number of 
@@ -997,7 +997,7 @@ def omp_get_max_threads():
 
 
 def omp_get_thread_num():
-    """Return the thread number, within the current team, of the 
+    """Return the thread number, within the current team, of the \
     calling thread.
     
     The omp_get_thread_num routine returns the thread number of the calling 
@@ -1041,7 +1041,7 @@ def omp_get_num_procs():
 
 
 def omp_in_parallel():
-    """Return 1 if the active-levels-var ICV is greater than zero; 
+    """Return 1 if the active-levels-var ICV is greater than zero; \
     otherwise, return 0.
     
     The effect of the omp_in_parallel routine is to return 1 if the current 
@@ -1063,8 +1063,8 @@ def omp_in_parallel():
 
 
 def omp_set_dynamic(dynamic_threads):
-    """Enable or disable dynamic adjustment of the number of threads 
+    """Enable or disable dynamic adjustment of the number of threads \
-    available for the execution of subsequent parallel regions by setting the 
+    available for the execution of subsequent parallel regions by setting the \
     value of the dyn-var ICV.
     
     For implementations that support dynamic adjustment of the number of 
@@ -1092,7 +1092,7 @@ def omp_set_dynamic(dynamic_threads):
 
 
 def omp_get_dynamic():
-    """Return the value of the dyn-var ICV, which determines whether
+    """Return the value of the dyn-var ICV, which determines whether \
     dynamic adjustment of the number of threads is enabled or disabled.
     
     This routine returns 1 if dynamic adjustment of the number of threads 
@@ -1140,7 +1140,7 @@ def omp_set_nested(nested):
 
 
 def omp_get_nested():
-    """Return the value of the nest-var ICV, which determines if nested 
+    """Return the value of the nest-var ICV, which determines if nested \
     parallelism is enabled or disabled
     
     This routine returns 1 if nested parallelism is enabled for the current 
@@ -1160,7 +1160,7 @@ def omp_get_nested():
 
 
 def omp_set_schedule(kind, modifier=0):
-    """Set the schedule that is applied when *runtime* is used as 
+    """Set the schedule that is applied when *runtime* is used as \
     schedule kind, by setting the value of the run-sched-var ICV.
     
     The effect of this routine is to set the value of the run-sched-var ICV 
@@ -1221,8 +1221,8 @@ def omp_get_schedule():
 
 
 def omp_get_thread_limit(): 
-    """Return the maximum number of OpenMP threads available to participate in 
+    """Return the maximum number of OpenMP threads available to participate \
-    the current contention group.
+    in the current contention group.
     
     The omp_get_thread_limit routine returns the value of the thread-limit-var 
     ICV.
@@ -1241,7 +1241,7 @@ def omp_get_thread_limit():
 
 
 def omp_set_max_active_levels(max_levels):
-    """Limit the number of nested active parallel regions on the device, 
+    """Limit the number of nested active parallel regions on the device, \
     by setting the max-active-levels-var ICV.
     
     The effect of this routine is to set the value of the max-active-levels-var 
@@ -1271,7 +1271,7 @@ def omp_set_max_active_levels(max_levels):
 
 
 def omp_get_max_active_levels():
-    """Return the value of the max-active-levels-var ICV, which determines 
+    """Return the value of the max-active-levels-var ICV, which determines \
     the maximum number of nested active parallel regions on the device
     
     The omp_get_max_active_levels routine returns the value of the 
@@ -1311,7 +1311,7 @@ def omp_get_level():
    
                           
 def omp_get_ancestor_thread_num(level):
-    """Return, for a given nested level of the current thread, the thread 
+    """Return, for a given nested level of the current thread, the thread \
     number of the ancestor of the current thread.
     
     The omp_get_ancestor_thread_num routine returns the thread number of the 
@@ -1338,7 +1338,7 @@ def omp_get_ancestor_thread_num(level):
 
 
 def omp_get_team_size(level):
-    """Return, for a given nested level of the current thread, the size 
+    """Return, for a given nested level of the current thread, the size \
     of the thread team to which the ancestor or the current thread belongs
     
     The omp_get_team_size routine returns the size of the thread team to which 
@@ -1385,7 +1385,7 @@ def omp_get_active_level():
 
 
 def omp_in_final():
-    """Return 1 (True) if the routine is executed in a final task region; 
+    """Return 1 (True) if the routine is executed in a final task region; \
     otherwise, it returns 0 (False).
     
     Returns:

src/wrf/util.py

@@ -2233,9 +2233,6 @@ def extract_vars(wrfin, timeidx, varnames, method="cat", squeeze=True,
             :data:`wrf.ALL_TIMES` (an alias for None) to return 
             all times in the file or sequence. The default is 0.
             
-        is_moving (:obj:`bool`): A boolean type that indicates if the 
-            sequence is a moving nest.
-            
         method (:obj:`str`, optional): The aggregation method to use for 
             sequences.  Must be either 'cat' or 'join'.  
             'cat' combines the data along the Time dimension.  

src/wrf/version.py

@@ -1,2 +1,2 @@
-__version__ = "1.0.5"
+__version__ = "1.1.0"
 

test/ipynb/Doc_Examples.ipynb

@@ -2,10 +2,7 @@
  "cells": [
   {
    "cell_type": "markdown",
-   "metadata": {
+   "metadata": {},
-    "deletable": true,
-    "editable": true
-   },
    "source": [
     "# Cartopy Examples"
    ]
@@ -13,11 +10,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "%matplotlib inline"
@@ -27,9 +20,6 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "collapsed": false,
-    "deletable": true,
-    "editable": true,
     "scrolled": false
    },
    "outputs": [],
@@ -115,11 +105,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "%matplotlib inline\n",
@@ -203,9 +189,6 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "collapsed": false,
-    "deletable": true,
-    "editable": true,
     "scrolled": false
    },
    "outputs": [],
@@ -343,9 +326,6 @@
    "cell_type": "code",
    "execution_count": null,
    "metadata": {
-    "collapsed": false,
-    "deletable": true,
-    "editable": true,
     "scrolled": false
    },
    "outputs": [],
@@ -411,11 +391,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "%matplotlib inline\n",
@@ -484,11 +460,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "from __future__ import (absolute_import, division, print_function, unicode_literals)\n",
@@ -568,11 +540,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "from __future__ import print_function\n",
@@ -606,11 +574,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "from __future__ import print_function\n",
@@ -631,10 +595,7 @@
   },
   {
    "cell_type": "markdown",
-   "metadata": {
+   "metadata": {},
-    "deletable": true,
-    "editable": true
-   },
    "source": [
     "# Basemap Examples"
    ]
@@ -642,11 +603,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "\n",
@@ -704,11 +661,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "%matplotlib inline\n",
@@ -781,11 +734,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "%matplotlib inline\n",
@@ -921,11 +870,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "from __future__ import print_function\n",
@@ -960,11 +905,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "from __future__ import print_function\n",
@@ -986,11 +927,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "# SLP\n",
@@ -1041,11 +978,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "from __future__ import print_function\n",
@@ -1079,11 +1012,7 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "metadata": {
+   "metadata": {},
-    "collapsed": false,
-    "deletable": true,
-    "editable": true
-   },
    "outputs": [],
    "source": [
     "from __future__ import print_function\n",
@@ -1101,6 +1030,233 @@
     "\n",
     "print (bounds)"
    ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# OpenMP Routines"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from __future__ import print_function\n",
+    "\n",
+    "from wrf import omp_enabled\n",
+    "\n",
+    "print(omp_enabled())"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from __future__ import print_function\n",
+    "\n",
+    "from wrf import omp_get_num_procs\n",
+    "\n",
+    "print(omp_get_num_procs())"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from __future__ import print_function\n",
+    "\n",
+    "from wrf import omp_set_num_threads, omp_get_max_threads\n",
+    "\n",
+    "omp_set_num_threads(4)\n",
+    "\n",
+    "print(omp_get_max_threads())"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from __future__ import print_function\n",
+    "\n",
+    "from wrf import omp_set_schedule, omp_get_schedule, OMP_SCHED_GUIDED\n",
+    "\n",
+    "omp_set_schedule(OMP_SCHED_GUIDED, 0)\n",
+    "\n",
+    "sched, modifier = omp_get_schedule()\n",
+    "\n",
+    "print(sched, modifier)"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Loop and Fill Technique"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from __future__ import print_function, division\n",
+    "\n",
+    "import numpy as np\n",
+    "from netCDF4 import Dataset\n",
+    "from wrf import getvar, ALL_TIMES\n",
+    "\n",
+    "filename_list = [\"/Users/ladwig/Documents/wrf_files/wrf_vortex_multi/wrfout_d02_2005-08-28_00:00:00\",\n",
+    "                 \"/Users/ladwig/Documents/wrf_files/wrf_vortex_multi/wrfout_d02_2005-08-28_12:00:00\", \n",
+    "                 \"/Users/ladwig/Documents/wrf_files/wrf_vortex_multi/wrfout_d02_2005-08-29_00:00:00\"]\n",
+    "\n",
+    "# Result shape (hardcoded for this example, modify as necessary)\n",
+    "result_shape = (9, 29, 96, 96)\n",
+    "\n",
+    "# Only need 4-byte floats\n",
+    "z_final = np.empty(result_shape, np.float32)\n",
+    "\n",
+    "# Modify this number if using more than 1 time per file\n",
+    "times_per_file = 4\n",
+    "\n",
+    "for timeidx in xrange(result_shape[0]):\n",
+    "   # Compute the file index and the time index inside the file\n",
+    "   fileidx = timeidx // times_per_file\n",
+    "   file_timeidx = timeidx % times_per_file\n",
+    "\n",
+    "   f = Dataset(filename_list[fileidx])  \n",
+    "   z = getvar(f, \"z\", file_timeidx)\n",
+    "\n",
+    "   z_final[timeidx,:] = z[:]\n",
+    "   f.close()\n"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Using the cache argument"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from __future__ import print_function\n",
+    "\n",
+    "import time\n",
+    "from netCDF4 import Dataset\n",
+    "from wrf import getvar, ALL_TIMES, extract_vars\n",
+    "\n",
+    "wrf_filenames = [\"/Users/ladwig/Documents/wrf_files/wrf_vortex_multi/wrfout_d02_2005-08-28_00:00:00\",\n",
+    "                 \"/Users/ladwig/Documents/wrf_files/wrf_vortex_multi/wrfout_d02_2005-08-28_12:00:00\", \n",
+    "                 \"/Users/ladwig/Documents/wrf_files/wrf_vortex_multi/wrfout_d02_2005-08-29_00:00:00\"]\n",
+    "\n",
+    "wrfin = [Dataset(x) for x in wrf_filenames]\n",
+    "\n",
+    "start = time.time()\n",
+    "my_cache = extract_vars(wrfin, ALL_TIMES, (\"P\", \"PSFC\", \"PB\", \"PH\", \"PHB\", \"T\", \"QVAPOR\", \n",
+    "                                           \"HGT\", \"U\", \"V\", \"W\"))\n",
+    "end = time.time()\n",
+    "print (\"Time taken to build cache: \", (end-start), \"s\")\n",
+    "\n",
+    "vars = (\"avo\", \"eth\", \"cape_2d\", \"cape_3d\", \"ctt\", \"dbz\", \"mdbz\", \n",
+    "        \"geopt\", \"helicity\", \"lat\", \"lon\", \"omg\", \"p\", \"pressure\", \n",
+    "        \"pvo\", \"pw\", \"rh2\", \"rh\", \"slp\", \"ter\", \"td2\", \"td\", \"tc\", \n",
+    "        \"theta\", \"tk\", \"tv\", \"twb\", \"updraft_helicity\", \"ua\", \"va\", \n",
+    "        \"wa\", \"uvmet10\", \"uvmet\", \"z\", \"cfrac\", \"zstag\", \"geopt_stag\")\n",
+    "\n",
+    "start = time.time()\n",
+    "for var in vars:\n",
+    "    v = getvar(wrfin, var, ALL_TIMES)\n",
+    "end = time.time()\n",
+    "no_cache_time = (end-start)\n",
+    "\n",
+    "print (\"Time taken without variable cache: \", no_cache_time, \"s\")\n",
+    "\n",
+    "start = time.time()\n",
+    "for var in vars:\n",
+    "    v = getvar(wrfin, var, ALL_TIMES, cache=my_cache)\n",
+    "end = time.time()\n",
+    "cache_time = (end-start)\n",
+    "\n",
+    "print (\"Time taken with variable cache: \", cache_time, \"s\")\n",
+    "\n",
+    "improvement = ((no_cache_time-cache_time)/no_cache_time) * 100 \n",
+    "print (\"The cache decreased computation time by: \", improvement, \"%\")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "# Using the cache argument with OpenMP"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from __future__ import print_function\n",
+    "\n",
+    "import time\n",
+    "from netCDF4 import Dataset\n",
+    "from wrf import getvar, ALL_TIMES, extract_vars, omp_set_num_threads, omp_get_num_procs\n",
+    "\n",
+    "wrf_filenames = [\"/Users/ladwig/Documents/wrf_files/wrf_vortex_multi/wrfout_d02_2005-08-28_00:00:00\",\n",
+    "                 \"/Users/ladwig/Documents/wrf_files/wrf_vortex_multi/wrfout_d02_2005-08-28_12:00:00\", \n",
+    "                 \"/Users/ladwig/Documents/wrf_files/wrf_vortex_multi/wrfout_d02_2005-08-29_00:00:00\"]\n",
+    "\n",
+    "wrfin = [Dataset(x) for x in wrf_filenames]\n",
+    "\n",
+    "start = time.time()\n",
+    "my_cache = extract_vars(wrfin, ALL_TIMES, (\"P\", \"PSFC\", \"PB\", \"PH\", \"PHB\", \"T\", \"QVAPOR\", \n",
+    "                                           \"HGT\", \"U\", \"V\", \"W\"))\n",
+    "end = time.time()\n",
+    "print (\"Time taken to build cache: \", (end-start), \"s\")\n",
+    "\n",
+    "omp_set_num_threads(omp_get_num_procs())\n",
+    "\n",
+    "vars = (\"avo\", \"eth\", \"cape_2d\", \"cape_3d\", \"ctt\", \"dbz\", \"mdbz\", \n",
+    "        \"geopt\", \"helicity\", \"lat\", \"lon\", \"omg\", \"p\", \"pressure\", \n",
+    "        \"pvo\", \"pw\", \"rh2\", \"rh\", \"slp\", \"ter\", \"td2\", \"td\", \"tc\", \n",
+    "        \"theta\", \"tk\", \"tv\", \"twb\", \"updraft_helicity\", \"ua\", \"va\", \n",
+    "        \"wa\", \"uvmet10\", \"uvmet\", \"z\", \"cfrac\", \"zstag\", \"geopt_stag\")\n",
+    "\n",
+    "start = time.time()\n",
+    "for var in vars:\n",
+    "    v = getvar(wrfin, var, ALL_TIMES)\n",
+    "end = time.time()\n",
+    "no_cache_time = (end-start)\n",
+    "\n",
+    "print (\"Time taken without variable cache: \", no_cache_time, \"s\")\n",
+    "\n",
+    "start = time.time()\n",
+    "for var in vars:\n",
+    "    v = getvar(wrfin, var, ALL_TIMES, cache=my_cache)\n",
+    "end = time.time()\n",
+    "cache_time = (end-start)\n",
+    "\n",
+    "print (\"Time taken with variable cache: \", cache_time, \"s\")\n",
+    "\n",
+    "improvement = ((no_cache_time-cache_time)/no_cache_time) * 100 \n",
+    "print (\"The cache decreased computation time by: \", improvement, \"%\")\n",
+    "\n",
+    "omp_set_num_threads(1)"
+   ]
   }
  ],
  "metadata": {
@@ -1119,7 +1275,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython2",
-   "version": "2.7.12"
+   "version": "2.7.13"
   }
  },
  "nbformat": 4,

test/varcache.py

@@ -0,0 +1,36 @@
+from __future__ import print_function
+
+import time
+from netCDF4 import Dataset
+from wrf import getvar, ALL_TIMES, extract_vars
+
+wrf_filenames = ["/Users/ladwig/Documents/wrf_files/wrf_vortex_multi/wrfout_d02_2005-08-28_00:00:00",
+                 "/Users/ladwig/Documents/wrf_files/wrf_vortex_multi/wrfout_d02_2005-08-28_12:00:00", 
+                 "/Users/ladwig/Documents/wrf_files/wrf_vortex_multi/wrfout_d02_2005-08-29_00:00:00"]
+
+wrfin = [Dataset(x) for x in wrf_filenames]
+
+my_cache = extract_vars(wrfin, ALL_TIMES, ("P", "PB", "PH", "PHB", "T", "QVAPOR", "HGT", "U", "V", "W", "PSFC"))
+
+start = time.time()
+for var in ("avo", "eth", "cape_2d", "cape_3d", "ctt", "dbz", "mdbz", 
+            "geopt", "helicity", "lat", "lon", "omg", "p", "pressure", 
+            "pvo", "pw", "rh2", "rh", "slp", "ter", "td2", "td", "tc", 
+            "theta", "tk", "tv", "twb", "updraft_helicity", "ua", "va", 
+            "wa", "uvmet10", "uvmet", "z", "cfrac", "zstag", "geopt_stag"):
+    v = getvar(wrfin, var, ALL_TIMES)
+end = time.time()
+
+print ("Time taken without variable cache: ", (end-start))
+
+start = time.time()
+for var in ("avo", "eth",  "cape_2d", "cape_3d", "ctt", "dbz", "mdbz", 
+            "geopt", "helicity", "lat", "lon", "omg", "p", "pressure", 
+            "pvo", "pw", "rh2", "rh", "slp", "ter", "td2", "td", "tc", 
+            "theta", "tk", "tv", "twb", "updraft_helicity", "ua", "va", 
+            "wa", "uvmet10", "uvmet", "z", "cfrac", "zstag", "geopt_stag"):
+    v = getvar(wrfin, var, ALL_TIMES, cache=my_cache)
+end = time.time()
+
+print ("Time taken with variable cache: ", (end-start))
+