diff --git a/docs/source/examples/example_2_wake_comparison.rst b/docs/source/examples/example_2_wake_comparison.rst
new file mode 100644
index 0000000000000000000000000000000000000000..dd8b2e2cdf9f4d6055a21458cd1d140c0d5d15c2
--- /dev/null
+++ b/docs/source/examples/example_2_wake_comparison.rst
@@ -0,0 +1,56 @@
+.. Example 2
+
+Example 2: Optimization with Different Wake Models
+===================================================
+
+This example demonstrates the optimization of a 3-turbine farm with
+the N. O. Jensen (NOJ) and G. C. Larsen (GCL) wake models implemented in
+FUSED-Wake.
+
+Specifications
+--------------
+
+- The farm is offshore.
+- The wind resource is defined by arrays of the frequency and Weibull
+ A and k parameters taken from a Horns Rev calculation.
+- There is a boundary beyond which the turbines cannot go.
+- The turbines cannot be closer than 2D.
+
+Results
+-------
+
+The optimization results are visualized in the figure below, and the
+script also outputs the following information on the resultings costs::
+
+ Comparison of cost models vs. optimized locations:
+
+ Cost | GCL_aep NOJ_aep
+ ---------------------------------
+ GCL_loc | -2704.40 -2688.57 (0.59%)
+ NOJ_loc | -2703.90 -2688.59 (0.57%)
+ (0.02%) (-0.00%)
+
+The subscript ``_aep`` indicates the AEP calculated using the specified wake
+model. The subscript ``_loc`` indicates the optimal locations determined when
+optimizing with the specified wake model. For example, the cost corresponding
+to the ``GCL_aep`` column and the ``NOJ_loc`` column is the AEP calculated
+using the GCL wake model for the optimal location results when the N. O.
+Jensen wake model was used in the optimization. The percentages in parentheses
+are the percent differences for the column or row.
+
+Important observations:
+
+- The upper and lower turbines are pushed away from each other, towards
+ the corners.
+- There is little difference in the AEP value when calculated using either of
+ the optimal locations, since the optimal locations themselves are quite
+ close.
+- The N. O. Jensen wake model has an AEP that is approximately 0.6% lower than
+ the AEP with the G. C. Larsen wake model for this farm.
+
+.. figure:: /../../examples/example_2_wake_comparison.png
+
+Code
+----
+
+.. literalinclude:: /../../examples/example_2_wake_comparison.py
diff --git a/docs/source/index.rst b/docs/source/index.rst
index 7757c90e4a5f4f9838493c86e5294421d2aa4080..b6edc02de16d8303d51950fd5c2c120d3408c668 100644
--- a/docs/source/index.rst
+++ b/docs/source/index.rst
@@ -16,3 +16,4 @@ Welcome to TOPFARM, the wind-farm optimizer
:maxdepth: 2
examples/example_1_constrained_layout_optimization
+ examples/example_2_wake_comparison
diff --git a/examples/example_2_wake_comparison.png b/examples/example_2_wake_comparison.png
new file mode 100644
index 0000000000000000000000000000000000000000..bf916cb9f311d27821d60d6b3c2ad49756fd2f2c
Binary files /dev/null and b/examples/example_2_wake_comparison.png differ
diff --git a/examples/example_2_wake_comparison.py b/examples/example_2_wake_comparison.py
new file mode 100644
index 0000000000000000000000000000000000000000..4ce725c0994adf1f70a095a15c82773227850cf9
--- /dev/null
+++ b/examples/example_2_wake_comparison.py
@@ -0,0 +1,130 @@
+"""Example: optimization with different wake models
+
+This example uses a dummy cost function to optimize a simple wind turbine
+layout that is subject to constraints. The optimization pushes the wind turbine
+locations to specified locations in the farm.
+"""
+import os
+import warnings
+
+from matplotlib.patches import Polygon
+import matplotlib.pyplot as plt
+import numpy as np
+
+from topfarm.plotting import PlotComp
+from topfarm import TopFarm
+from topfarm.cost_models.utils.wind_resource import WindResource
+from topfarm.cost_models.cost_model_wrappers import AEPCostModelComponent
+from topfarm.cost_models.fused_wake_wrappers import FusedWakeGCLWakeModel, \
+ FusedWakeNOJWakeModel
+from topfarm.cost_models.utils.aep_calculator import AEPCalculator
+
+
+# ------------------------ INPUTS ------------------------
+
+# paths to input files
+test_files_dir = os.path.join(os.path.dirname(__file__), '..', 'topfarm',
+ 'tests', 'test_files') # file locations
+wf_path = os.path.join(test_files_dir, 'wind_farms',
+ '3tb.yml') # path to wind farm
+f = [3.597152, 3.948682, 5.167395, 7.000154, 8.364547, 6.43485, 8.643194,
+ 11.77051, 15.15757, 14.73792, 10.01205, 5.165975] # horns rev
+a = [9.176929, 9.782334, 9.531809, 9.909545, 10.04269, 9.593921, 9.584007,
+ 10.51499, 11.39895, 11.68746, 11.63732, 10.08803] # horns rev
+k = [2.392578, 2.447266, 2.412109, 2.591797, 2.755859, 2.595703, 2.583984,
+ 2.548828, 2.470703, 2.607422, 2.626953, 2.326172] # horns rev
+rot_diam = 80.0 # rotor diameter [m]
+bnd_size = 2 * rot_diam + 10 # boundary size
+init_pos = np.array([(0, 2 * rot_diam), (0, 0),
+ (0, -2 * rot_diam)]) # initial turbine positions
+boundary = [(-bnd_size, bnd_size), (bnd_size, bnd_size),
+ (bnd_size, -bnd_size), (-bnd_size, -bnd_size),
+ (-bnd_size, bnd_size)] # wind farm boundary
+min_spacing = 2.0 * rot_diam # minimum spacing between turbines [m]
+
+# ------------------------ DEFINE WIND RESOURCE ------------------------
+
+wind_res = WindResource(f, a, k, np.zeros_like(k))
+
+# ------------------------ OPTIMIZATION ------------------------
+
+warnings.filterwarnings('ignore') # temporarily disable fusedwake warnings
+
+# GCL: define the wake model, aep calculator, and optimization problem
+wake_mod_gcl = FusedWakeGCLWakeModel(wf_path)
+aep_calc_gcl = AEPCalculator(wind_res, wake_mod_gcl)
+tf_gcl = TopFarm(init_pos, aep_calc_gcl.get_TopFarm_cost_component(),
+ min_spacing, boundary=boundary)
+
+# NOJ: define the wake model, aep calculator, and optimization problem
+wake_mod_noj = FusedWakeNOJWakeModel(wf_path)
+aep_calc_noj = AEPCalculator(wind_res, wake_mod_noj)
+tf_noj = TopFarm(init_pos, aep_calc_noj.get_TopFarm_cost_component(),
+ min_spacing, boundary=boundary)
+
+# run the optimization
+cost_gcl, state_gcl, recorder_gcl = tf_gcl.optimize()
+cost_noj, state_noj, recorder_noj = tf_noj.optimize()
+
+# ------------------------ POST-PROCESS ------------------------
+
+# get the optimized locations
+opt_gcl = tf_gcl.turbine_positions
+opt_noj = tf_noj.turbine_positions
+
+# create the array of costs for easier printing
+costs = np.diag([cost_gcl, cost_noj])
+costs[0, 1] = TopFarm(opt_gcl, aep_calc_noj.get_TopFarm_cost_component(),
+ min_spacing,
+ boundary=boundary).evaluate()[0] # noj cost, gcl locs
+costs[1, 0] = TopFarm(opt_noj, aep_calc_gcl.get_TopFarm_cost_component(),
+ min_spacing,
+ boundary=boundary).evaluate()[0] # gcl cost, noj locs
+
+warnings.filterwarnings('default') # re-enable warnings warnings
+
+# ------------------------ PRINT STATS ------------------------
+
+aep_diffs = 200 * (costs[:, 0] - costs[:, 1]) / (costs[:, 0] + costs[:, 1])
+loc_diffs = 200 * (costs[0, :] - costs[1, :]) / (costs[0, :] + costs[1, :])
+
+print('\nComparison of cost models vs. optimized locations:')
+print('\nCost | GCL_aep NOJ_aep')
+print('---------------------------------')
+print(f'GCL_loc |{costs[0,0]:11.2f} {costs[0,1]:11.2f}' +
+ f' ({aep_diffs[0]:.2f}%)')
+print(f'NOJ_loc |{costs[1,0]:11.2f} {costs[1,1]:11.2f}' +
+ f' ({aep_diffs[1]:.2f}%)')
+print(f' ({loc_diffs[0]:.2f}%) ({loc_diffs[1]:.2f}%)')
+
+# ------------------------ PLOT (if possible) ------------------------
+
+try:
+
+ # initialize the figure and axes
+ fig = plt.figure(1, figsize=(7, 5))
+ plt.clf()
+ ax = plt.axes()
+
+ # plot the boundary and desired locations
+ ax.add_patch(Polygon(boundary, fill=False,
+ label='Boundary')) # boundary
+ ax.plot(init_pos[:, 0], init_pos[:, 1], 'xk',
+ label='Initial')
+ ax.plot(opt_gcl[:, 0], opt_gcl[:, 1], 'o',
+ label='GCL')
+ ax.plot(opt_noj[:, 0], opt_noj[:, 1], '^',
+ label='NOJ')
+
+ # make a few adjustments to the plot
+ ax.autoscale_view() # autoscale the boundary
+ plt.legend(bbox_to_anchor=(0., 1.02, 1., .102), loc=3,
+ ncol=4, mode='expand', borderaxespad=0.) # add a legend
+ plt.tight_layout() # zoom the plot in
+ plt.axis('off') # remove the axis
+
+ # save the png
+ fig.savefig(os.path.basename(__file__).replace('.py', '.png'))
+
+except RuntimeError:
+ pass