Fixed calc and plot. Added power station to README.

README.md

@@ -50,18 +50,60 @@ The basic steps in evaluating a model are:
 * evaluation of the contract model.
 
 A convenience function `calc_and_plot()` can be used to perform all the steps.
-
 ```python
 from quantdsl.interfaces.calcandplot import calc_and_plot
+```
 
-source_code = """from quantdsl.lib.storage2 import GasStorage
+Here's an evaluation of a gas storage facility.
+
+```pyt hon
+from quantdsl.interfaces.calcandplot import calc_and_plot
 
-GasStorage(Date('2011-6-1'), Date('2011-9-1'), 'GAS', 0, 0, 50000, TimeDelta('1m'), 'monthly')
-"""
 
 calc_and_plot(
     title="Gas Storage",
-    source_code=source_code,
+    source_code="""def GasStorage(start, end, commodity_name, quantity, target, limit, step, period):
+    if ((start < end) and (limit > 0)):
+        if quantity <= 0:
+            return Wait(start, Choice(
+                Continue(start, end, commodity_name, quantity, limit, step, period, target),
+                Inject(start, end, commodity_name, quantity, limit, step, period, target, 1),
+            ))
+        elif quantity >= limit:
+            return Wait(start, Choice(
+                Continue(start, end, commodity_name, quantity, limit, step, period, target),
+                Inject(start, end, commodity_name, quantity, limit, step, period, target, -1),
+            ))
+        else:
+            return Wait(start, Choice(
+                Continue(start, end, commodity_name, quantity, limit, step, period, target),
+                Inject(start, end, commodity_name, quantity, limit, step, period, target, 1),
+                Inject(start, end, commodity_name, quantity, limit, step, period, target, -1),
+            ))
+    else:
+        if target < 0 or target == quantity:
+            return 0
+        else:
+            return BreachOfContract()
+
+
+@inline
+def BreachOfContract():
+    -10000000000000000
+
+@inline
+def Continue(start, end, commodity_name, quantity, limit, step, period, target):
+    GasStorage(start + step, end, commodity_name, quantity, target, limit, step, period)
+
+
+@inline
+def Inject(start, end, commodity_name, quantity, limit, step, period, target, vol):
+    Continue(start, end, commodity_name, quantity + vol, limit, step, period, target) - \
+    vol * Lift(commodity_name, period, Market(commodity_name))
+
+
+GasStorage(Date('2011-6-1'), Date('2011-9-1'), 'GAS', 0, 0, 50000, TimeDelta('1m'), 'monthly')
+""",
     observation_date='2011-1-1',
     interest_rate=2.5,
     path_count=20000,
@@ -105,3 +147,116 @@ calc_and_plot(
 )
 
 ```
+
+Here's an evaluation of a power station.
+
+```python
+calc_and_plot(
+    title="Gas Storage",
+    source_code="""
+def PowerStation(start, end, gas, power, duration_off):
+    if (start < end):
+        Wait(start,
+            Choice(
+                ProfitFromRunning(gas, power, duration_off) + PowerStation(
+                    Tomorrow(start), end, gas, power, Running()
+                ),
+                PowerStation(
+                    Tomorrow(start), end, gas, power, Stopped(duration_off)
+                )
+            )
+        )
+    else:
+        return 0
+
+@inline
+def ProfitFromRunning(gas, power, duration_off):
+    if duration_off > 1:
+        return 0.75 * power - gas
+    elif duration_off == 1:
+        return 0.90 * power - gas
+    else:
+        return 1.00 * power - gas
+
+@inline
+def Running():
+    return 0
+
+@inline
+def Stopped(duration_off):
+    return duration_off + 1
+
+@inline
+def Tomorrow(today):
+    return today + TimeDelta('1d')
+
+PowerStation(Date('2012-01-01'), Date('2012-01-13'), Market('GAS'), Market('POWER'), Running())
+""",
+    observation_date='2011-1-1',
+    interest_rate=2.5,
+    path_count=20000,
+    perturbation_factor=0.01,
+    periodisation='monthly',
+    price_process={
+        'name': 'quantdsl.priceprocess.blackscholes.BlackScholesPriceProcess',
+        'market': ['GAS', 'POWER'],
+        'sigma': [0.5, 0.3],
+        'rho': [[1.0, 0.8], [0.8, 1.0]],
+        'curve': {
+            'GAS': [
+                ('2011-1-1', 13.5),
+                ('2011-2-1', 11.0),
+                ('2011-3-1', 10.0),
+                ('2011-4-1', 9.0),
+                ('2011-5-1', 7.5),
+                ('2011-6-1', 7.0),
+                ('2011-7-1', 6.5),
+                ('2011-8-1', 7.5),
+                ('2011-9-1', 8.5),
+                ('2011-10-1', 10.0),
+                ('2011-11-1', 11.5),
+                ('2011-12-1', 12.0),
+                ('2012-1-1', 13.5),
+                ('2012-2-1', 11.0),
+                ('2012-3-1', 10.0),
+                ('2012-4-1', 9.0),
+                ('2012-5-1', 7.5),
+                ('2012-6-1', 7.0),
+                ('2012-7-1', 6.5),
+                ('2012-8-1', 7.5),
+                ('2012-9-1', 8.5),
+                ('2012-10-1', 10.0),
+                ('2012-11-1', 11.5),
+                ('2012-12-1', 12.0)
+            ],
+            'POWER': [
+                ('2011-1-1', 13.5),
+                ('2011-2-1', 11.0),
+                ('2011-3-1', 10.0),
+                ('2011-4-1', 9.0),
+                ('2011-5-1', 7.5),
+                ('2011-6-1', 7.0),
+                ('2011-7-1', 6.5),
+                ('2011-8-1', 7.5),
+                ('2011-9-1', 8.5),
+                ('2011-10-1', 10.0),
+                ('2011-11-1', 11.5),
+                ('2011-12-1', 12.0),
+                ('2012-1-1', 13.5),
+                ('2012-2-1', 11.0),
+                ('2012-3-1', 10.0),
+                ('2012-4-1', 9.0),
+                ('2012-5-1', 7.5),
+                ('2012-6-1', 7.0),
+                ('2012-7-1', 6.5),
+                ('2012-8-1', 7.5),
+                ('2012-9-1', 8.5),
+                ('2012-10-1', 10.0),
+                ('2012-11-1', 11.5),
+                ('2012-12-1', 12.0)
+            ]
+        }
+    }
+)
+
+```

quantdsl/interfaces/calcandplot.py

@@ -61,7 +61,6 @@ def run(self, title, source_code, observation_date, interest_rate, path_count, p
                 if self.is_evaluation_ready.wait(timeout=2):
                     break
 
-            unsubscribe(self.is_evaluation_complete, self.on_evaluation_complete)
 
             fair_value_stderr, fair_value_mean, periods = self.read_results(app, evaluation, market_simulation,
                                                                             path_count)
@@ -75,6 +74,9 @@ def run(self, title, source_code, observation_date, interest_rate, path_count, p
             if not supress_plot and plt and len(periods) > 1:
                 self.plot_results(interest_rate, path_count, perturbation_factor, periods, title, periodisation)
 
+            unsubscribe(self.is_result_value_computed, self.count_result_values_computed)
+            unsubscribe(self.is_evaluation_complete, self.on_evaluation_complete)
+
     def calc_results(self, app, interest_rate, observation_date, path_count, perturbation_factor,
                      contract_specification, price_process_name, calibration_params):
         market_calibration = app.register_market_calibration(
@@ -97,7 +99,6 @@ def calc_results(self, app, interest_rate, observation_date, path_count, perturb
 
         subscribe(self.is_result_value_computed, self.count_result_values_computed)
         evaluation = app.evaluate(contract_specification, market_simulation)
-        unsubscribe(self.is_result_value_computed, self.count_result_values_computed)
         return evaluation, market_simulation
 
     def is_result_value_computed(self, event):