Investigate PRICE_COMMODITY spikes with non-uniform duration_period

[Wegatriespython]

What happened?

For the baseline scenario under model SSP_SSP2_v6.1, price commodity for electricity shows a spike at 2060.

This has important downstream implications as any technology portfolio with electricity inputs in some small neighbourhood will see a basis flip to the cheaper option under the price spike barring other constraints. This happens for instance in the water model, where groundwater extraction is massively reduced in favour of a pre-dominantly surface water strategy which is problematic.

Now from a purely economic perspective such a price shock should be arbitraged away by the central planner under perfect foresight mode. Since the model doesn't do this, the only implication is that its cheaper for the model to have the price shock than to arbitrage it away. This however would invite the external validity Lucas critique as real world agents would be able to exploit an arbitrage the model cannot, breaking the policy invariance of the model.

The fact that this happens at the shift in the model timestep from 5-10 year periods is also suspicious.

• This means that either there is some misspecification in how the model is discounting Present Value at the period length change at 2060
• Or there are some discrete constraints creating a lumpy choice that forces this behaviour.

Code Sample

What did you expect to happen?

A strict 0 arbitrage price trajectory.

Versions

ixmp:              3.11.2.dev44+g1730afd24
              git: 1730afd (HEAD -> main, origin/main, origin/HEAD) Merge pull request #587 from iiasa/fix/gamsmodel-init
message_ix:        3.11.2.dev67+g00bec749
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
message_ix_models: (not installed)
message_data:      (not installed)

click:             8.2.1
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
dask:              2025.5.1
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
genno:             1.28.2
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
graphviz:          0.21
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
ixmp4:             0.10.0
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
jpype:             1.5.2
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
     Java VM path: /usr/lib/jvm/java-21-openjdk-amd64/lib/server/libjvm.so
openpyxl:          3.1.5
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
pandas:            2.2.3
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
pint:              0.24.4
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
xarray:            2025.4.0
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
yaml:              6.0.2
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock

iam_units:         2023.9.12
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
jupyter:           5.8.1
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
matplotlib:        3.10.3
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
plotnine:          0.15.0
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock
pyam:              3.0.0
              git: 00bec749 (HEAD -> main, origin/main, origin/HEAD) Merge pull request #451 from iiasa/material_stock

GAMS:              46.5.0
       system dir: /home/raghunathan/opt/gams/gams46.5_linux_x64_64_sfx

python:            3.12.3 (main, Jun 18 2025, 17:59:45) [GCC 13.3.0]
python-bits:       64
OS:                Linux
OS-release:        6.6.87.1-microsoft-standard-WSL2
machine:           x86_64
processor:         x86_64
byteorder:         little
LC_ALL:            None
LANG:              C.UTF-8
LOCALE:            ('C', 'UTF-8')

Additional Context

No response

[khaeru]

Thanks for filing this. I've removed the 'timeslice' label, as that pertains to "sub-annual time slices", and this appears to occur in a scenario that does not use them.

In Slack, someone suggested that the culprit may be df_period, which is used in the calculation of PRICE_COMMODITY:

message_ix/message_ix/model/MESSAGE/model_solve.gms

Lines 53 to 55 in cede75c

	PRICE_COMMODITY.l(node,commodity,level,year,time) = (
	- COMMODITY_BALANCE_AUX.M(node,commodity,level,year,time) / df_period(year)
	);

and computed here:

message_ix/message_ix/model/includes/period_parameter_assignment.gms

Lines 90 to 97 in 3742426

	df_period(year_all) =
	df_year(year_all) * (
	* multiply the per-year discount factor by the geometric series of over the duration of the period
	( ( POWER( 1 + interestrate(year_all) , duration_period(year_all) ) - 1 )
	/ interestrate(year_all) )$( interestrate(year_all) )
	* if interest rate = 0, multiply by the number of years in that period
	+ ( duration_period(year_all) )$( interestrate(year_all) eq 0 ) )
	;

@Wegatriespython I believe you shared specific df_period values in Slack. Can you please also include them in the issue description, and ideally as text rather than a screenshot image?

[Wegatriespython]

The value calculation for df_period is correct as per the formula. Initially I had thought the value for 2060 was suspicious as it was lower than the one for 2070 but that is a consequence of 2070 having a longer period. At 5% the length of the period outweighs the later time.

Which means that df_period does not have an implementation issue. The problem is that then the change in duration period affects the actual economic logic of the model.

For the central planner, then the decision to invest in capacity in 2060 is weighted by the fact that waiting one more period gives 10 years of returns vs 5 years of return for the investment in 2060. Hence the incentive to allow for a price shock in 2060 because the overall ROI is higher by waiting till 2070 for that capital.

This can be tested by setting up a model, one with uniform periods, and the the other with different period lengths. Since the specific price shock comes from the full model a more direct test of change in investment schedules (investment deferral) can proxy the dominant mechanism of issue here.

However since the 5->10year periods are there as a computational trick to improve performance, this behaviour is undesirable(for all the reasons mentioned in the first post on thread). Worth thinking of different approaches to address this.

year_all	df_year	df_period	ratio
1950	1.00000e+00	5.52563e+00	5.52563
1955	7.83526e-01	4.32948e+00	5.52563
1960	6.13913e-01	3.39226e+00	5.52563
1965	4.81017e-01	2.65792e+00	5.52563
1970	3.76889e-01	2.08255e+00	5.52563
1975	2.95303e-01	1.63173e+00	5.52563
1980	2.31377e-01	1.27851e+00	5.52563
1985	1.81292e-01	1.00174e+00	5.52563
1990	1.42046e-01	7.84892e-01	5.52563
1995	1.11297e-01	6.14983e-01	5.52563
2000	8.72037e-02	4.81856e-01	5.52563
2005	6.83264e-02	3.77547e-01	5.52563
2010	5.35355e-02	2.95818e-01	5.52563
2015	4.19465e-02	2.31781e-01	5.52563
2020	3.28662e-02	1.81606e-01	5.52563
2025	2.57515e-02	1.42293e-01	5.52563
2030	2.01772e-02	1.11491e-01	5.52563
2035	1.58092e-02	8.73557e-02	5.52563
2040	1.23869e-02	6.84455e-02	5.52563
2045	9.70547e-03	5.36289e-02	5.52563
2050	7.60449e-03	4.20196e-02	5.52563
2055	5.95832e-03	3.29235e-02	5.52563
2060	4.66850e-03	2.57964e-02	5.52563
2070	2.86605e-03	3.60449e-02	12.5751
2080	1.75951e-03	2.21309e-02	12.5751
2090	1.08018e-03	1.35865e-02	12.5751
2100	6.63144e-04	8.34090e-03	12.5751
2110	4.07114e-04	5.12059e-03	12.5751

[khaeru]

OK, thanks for the update.

Which means that df_period does not have an implementation issue. The problem is that then the change in duration period affects the actual economic logic of the model.

IMO that sensitivity means there is an implementation issue. The model solution should be stable regardless of how time is discretized. In order for this to be true, sometimes it is necessary to use a more 'complicated' expression. Often a 'simple' expression of intended logic will give correct values in some cases, but not all—for instance, at a change in duration_period.

A recent example is #924; see https://gist.github.com/khaeru/59e362941f3eb84afcf08d66d3f36c72. Cell 8 in the notebook shows an expression which is more 'complicated' than what was previously in the GAMS code, but which actually gives correct results for changing period durations.

I think if we went through a similar (careful, first-principles) derivation of df_period and PRICE_COMMODITY, we would be able to reach a symbolic expression and then confirm whether the current code aligns.

[Wegatriespython]

OK, thanks for the update.

Which means that df_period does not have an implementation issue. The problem is that then the change in duration period affects the actual economic logic of the model.

IMO that sensitivity means there is an implementation issue. The model solution should be stable regardless of how time is discretized. In order for this to be true, sometimes it is necessary to use a more 'complicated' expression. Often a 'simple' expression of intended logic will give correct values in some cases, but not all—for instance, at a change in duration_period.

A recent example is #924; see https://gist.github.com/khaeru/59e362941f3eb84afcf08d66d3f36c72. Cell 8 in the notebook shows an expression which is more 'complicated' than what was previously in the GAMS code, but which actually gives correct results for changing period durations.

I think if we went through a similar (careful, first-principles) derivation of df_period and PRICE_COMMODITY, we would be able to reach a symbolic expression and then confirm whether the current code aligns.

Wouldn't that imply we can discretize at no information loss ?

Agreed though with the goal of making the model robust, my intuition is that more that just the duration_period will need to be touched for this.

[khaeru]

Wouldn't that imply we can discretize at no information loss ?

To be more precise with the example in the notebook I linked, the code could not (and I guess does not) claim to model any arbitrary "continuous function $n(t)$ which is additions to capacity per unit time". The derived expression there only works in the case where $n(t)$ is piecewise exponential in form (drawn in red on all the plots). We can imagine other functions that are continuous (splines, sinusoids) for which the formulation would not capture the actual ratio of $N_b$ to $N_a$.

In that case (growth_activity_up), our documentation should remind the user that the expression provides a consistent discretization of an upper bound at a certain (piecewise) rate of exponential growth. And the examples in the notebook show that the new formulation gives the same results for that functional form and a few alternate discretizations of time.

So here, by similar reasoning (and to be clear I don't want to do the whole "(careful, first-principles) derivation" in the comments, just sketching), we could imagine a bunch of cases:

1. 10 periods of uniform duration_period = 1 year. (Without using the "subannual time slice" feature, this is the finest resolution that the MESSAGE formulation can handle.)
2. 1 period of duration_period = 10 years. (This is the coarsest resolution covering the same model horizon.)
3. 2 periods of uniform duration_period = 5 years.
4. 5 periods of uniform duration_period = 2 years.
5. 2 periods of non-uniform duration_period = (2, 8) years.
6. Alternate versions of (5).
7. 3 periods of non-uniform duration_period = (1, 8, 1) years.
8. etc.

We could then:

• Express what we expect to be invariant across theses cases—values of CAP, ACT, etc.
• Observe and understand how the actual marginals of COMMODITY_BALANCE_AUX relate to those.
• Think about the notion of "commodity price" (in general) and decide how we believe those values should look: in general and in each case.

This should tell us how to derive a consistent expression (including any intermediate quantities) for PRICE_COMMODITY.

[khaeru]

Also, mulling over this I realized why the term "price shock" wasn't sitting right. This is a term used in some kinds of economic modeling, e.g. CGE modeling. Crucially, in such models, prices are endogenous: forcing or constraining a price to be higher or lower than a reference value is a "shock" in the sense that it forces the model to find a different equilibrium.

In MESSAGE, PRICE_COMMODITY is a reporting variable; it does not enter into the objective function and is computed only after the model is solved (per the first code snippet linked in my earlier #973 (comment) —observe the solve statement several lines earlier). So (apparent) jumps and drops in prices are not the only or chief influence on the model solution. If there is some constraint other than COMMODITY_BALANCE with a larger marginal, the model will still choose to consume more at given indices (n,c,l,y,h) than at other indices (e.g. y) where PRICE_COMMODITY(n,c,l,y,h) is lower.

For the central planner, then the decision to invest in capacity in 2060 is weighted by the fact that waiting one more period gives 10 years of returns vs 5 years of return for the investment in 2060.

As far as I understand, remaining_capacity as derived from technical_lifetime is supposed to prevent this, by ensuring that the number of years of returns a given vintage gets is not affected by the duration_period of the period in which it's constructed. For instance, something constructed with technical_lifetime = 10 years:

• as of the start of period '2070' = 2061-01-01 to 2070-12-31 —should only be active within that one period.
• as of the start of period '2060' = 2056-01-01 to 2060-12-31 —should last until 2065-12-31, i.e. halfway through the period '2070', thus should have remaining_capacity(n,t,2065,2070) = 0.5.

Again, this is what I understand to be the intended representation in the GAMS code. If it has some implementation error, the solver could choose to exploit it.

The example cases in my previous comment would be the best way to investigate this.

[Wegatriespython]

https://github.com/Wegatriespython/message_ix/blob/price_spike/%23973/message_ix/tests/test_nonuniform_periods.py

Got a simple reproduction of the issue going.

Structure:

Single region, single commodity (electricity), single technology (gas power plant)

• Growing demand: 100 * (1.05^t_elapsed) GWa for t_elapsed = year_i - start_year.
• Gas plant: 800 USD/GWa investment cost, 50 USD/GWa variable cost, 30-year lifetime
• Capacity factor: 1.0
• Historical year = 2015

In general price in excess capacity regimes (CAP_NEW = 0): Price = Variable Cost.
In capacity investment regimes: Price = Variable cost + Periodized Capital Cost.
For the special case TL =1, the periodized capital cost is the full investment cost.
Therefore the price spikes are coming from the periodization of capital costs over different period lengths.

Observe the starting price of the all 10 year periods at interest rate 5 &10% (Cases 2 and 3).
Now compare them to the 2040 prices in the mixed periods (Case 4 & 5) at the respective rates. They are identical. This makes sense as the first period in all the 10 year cases is exactly after a historical period 5 years in length.

Given the model treats investment cost as a lump sum amount, the periodization of capital costs is grid dependent on the discretization schema used.

To make this time invariant, we could consider using a capital recovery factor based approach. This can give us a grid invariant annualized capital cost.

---

Difference in Period Structure and Interest Rates

1. All 5-year periods (5%)

• Years: [2020, 2025, 2030, 2035, 2040, 2045, 2050, 2055, 2060, 2065, 2070]
• Interest Rate: 0.05

Year	ACT	CAP_NEW	CAP	Demand	Price	Dual_COMM_BAL_AUX
2020	100.000000	20.000000	100.000000	100.000000	95.142130	-525.720324
2025	127.628156	5.525631	127.628156	127.628156	95.093740	-411.706130
2030	162.889463	7.052261	162.889463	162.889463	95.016513	-322.320552
2035	207.892818	9.000671	207.892818	207.892818	94.893931	-252.220771
2040	265.329771	11.487391	265.329771	265.329771	94.700677	-197.219113
2045	338.635494	14.661145	338.635494	338.635494	95.538169	-155.892900
2050	432.194238	38.711749	432.194238	432.194238	95.407693	-121.979352
2055	551.601537	29.407091	551.601537	551.601537	95.198558	-95.364514
2060	703.998871	37.531728	703.998871	703.998871	94.864787	-74.458618
2065	898.500779	47.901053	898.500779	898.500779	94.334992	-58.014461
2070	1146.739979	61.135230	1146.739979	1146.739979	93.499760	-45.053387

2. All 10-year periods (5%)

• Years: [2020, 2030, 2040, 2050, 2060, 2070]
• Interest Rate: 0.05

Year	ACT	CAP_NEW	CAP	Demand	Price	Dual_COMM_BAL_AUX
2020	100.000000	20.000000	100.000000	100.000000	108.885432	-1369.549264
2030	162.889463	6.288946	162.889463	162.889463	90.182678	-696.366736
2040	265.329771	10.244031	265.329771	265.329771	89.959040	-426.448619
2050	432.194238	21.686447	432.194238	432.194238	90.558142	-263.545991
2060	703.998871	38.469410	703.998871	703.998871	90.175646	-161.110995
2070	1146.739979	54.518142	1146.739979	1146.739979	89.209096	-97.848025

3. All 10-year periods (10%)

• Years: [2020, 2030, 2040, 2050, 2060, 2070]
• Interest Rate: 0.1

Year	ACT	CAP_NEW	CAP	Demand	Price
2020	100.000000	20.000000	100.000000	100.000000	133.403381
2030	162.889463	6.288946	162.889463	162.889463	102.143920
2040	265.329771	10.244031	265.329771	265.329771	102.044335
2050	432.194238	21.686447	432.194238	432.194238	102.412176
2060	703.998871	38.469410	703.998871	703.998871	102.128049
2070	1146.739979	54.518142	1146.739979	1146.739979	100.390349

4. Mixed 5 and 10-year periods (5%)

• Years: [2020, 2025, 2030, 2035, 2040, 2050, 2060, 2070]
• Interest Rate: 0.05

Year	ACT	CAP_NEW	CAP	Demand	Price	Dual_COMM_BAL_AUX
2020	100.000000	20.000000	100.000000	100.000000	95.316399	-526.683273
2025	127.628156	5.525631	127.628156	127.628156	92.925410	-402.318394
2030	162.889463	7.052261	162.889463	162.889463	95.215672	-322.996151
2035	207.892818	9.000671	207.892818	207.892818	92.725978	-246.458519
2040	265.329771	11.487391	265.329771	265.329771	108.530174	-226.019764
2050	432.194238	21.686447	432.194238	432.194238	90.558142	-263.545991
2060	703.998871	36.706344	703.998871	703.998871	90.175646	-161.110995
2070	1146.739979	53.409359	1146.739979	1146.739979	89.209096	-97.848025

5. Mixed 5 and 10-year periods (10%)

• Years: [2020, 2025, 2030, 2035, 2040, 2050, 2060, 2070]
• Interest Rate: 0.1

Year	ACT	CAP_NEW	CAP	Demand	Price	Dual_COMM_BAL_AUX
2020	100.000000	20.000000	100.000000	100.000000	114.573138	-699.480467
2025	127.628156	5.525631	127.628156	127.628156	113.073745	-428.638456
2030	162.889463	7.052261	162.889463	162.889463	114.551885	-269.629974
2035	207.892818	9.000671	207.892818	207.892818	112.993035	-165.140722
2040	265.329771	11.487391	265.329771	265.329771	133.252144	-120.924217
2050	432.194238	21.686447	432.194238	432.194238	102.412176	-93.538238
2060	703.998871	36.706344	703.998871	703.998871	102.128049	-35.962988
2070	1146.739979	53.409359	1146.739979	1146.739979	100.390349	-13.629372

6. Mixed 5 and 10-year periods (0%)

• Years: [2020, 2025, 2030, 2035, 2040, 2050, 2060, 2070]
• Interest Rate: 0.0

Year	ACT	CAP_NEW	CAP	Demand	Price	Dual_COMM_BAL_AUX
2020	100.000000	20.000000	100.000000	100.000000	76.666667	-383.333333
2025	127.628156	5.525631	127.628156	127.628156	76.666667	-383.333333
2030	162.889463	7.052261	162.889463	162.889463	76.666667	-383.333333
2035	207.892818	9.000671	207.892818	207.892818	76.666667	-383.333333
2040	265.329771	11.487391	265.329771	265.329771	76.666667	-383.333333
2050	432.194238	21.686447	432.194238	432.194238	76.666667	-766.666667
2060	703.998871	36.706344	703.998871	703.998871	76.666667	-766.666667
2070	1146.739979	53.409359	1146.739979	1146.739979	76.666667	-766.666667

---

Large Initial CAP to Neutralize CAP_NEW at Kink Year

This scenario assumes 5000 historical new capacity with a technical lifetime of 30. New capacity (CAP_NEW) is only required post-kink year.

7. Mixed 5 and 10-year periods (5%)

• Years: [2020, 2025, 2030, 2035, 2040, 2050, 2060, 2070]
• Interest Rate: 0.05

Year	ACT	CAP_NEW	CAP	Demand	Price	Dual_COMM_BAL_AUX
2020	100.000000	0.000000	5000.000000	100.000000	50.000000	-276.281563
2025	127.628156	0.000000	5000.000000	127.628156	50.000000	-216.473834
2030	162.889463	0.000000	5000.000000	162.889463	50.000000	-169.612913
2035	207.892818	0.000000	5000.000000	207.892818	50.000000	-132.896155
2040	265.329771	0.000000	5000.000000	265.329771	50.000000	-104.127615
2050	432.194238	43.219424	432.194238	432.194238	90.558142	-263.545991
2060	703.998871	27.180463	703.998871	703.998871	90.175646	-161.110995
2070	1146.739979	44.274111	1146.739979	1146.739979	89.209096	-97.848025

8. All 5-year periods (5%) - Counterfactual

• Years: [2020, 2025, 2030, 2035, 2040, 2045, 2050, 2055, 2060, 2065, 2070]
• Interest Rate: 0.05

Year	ACT	CAP_NEW	CAP	Demand	Price	Dual_COMM_BAL_AUX
2020	100.000000	0.000000	5000.000000	100.000000	50.000000	-276.281563
2025	127.628156	0.000000	5000.000000	127.628156	50.000000	-216.473834
2030	162.889463	0.000000	5000.000000	162.889463	50.000000	-169.612913
2035	207.892818	0.000000	5000.000000	207.892818	50.000000	-132.896155
2040	265.329771	0.000000	5000.000000	265.329771	50.000000	-104.127615
2045	338.635494	67.727099	338.635494	338.635494	95.538169	-155.892900
2050	432.194238	18.711749	432.194238	432.194238	95.407693	-121.979352
2055	551.601537	23.881460	551.601537	551.601537	95.198558	-95.364514
2060	703.998871	30.479467	703.998871	703.998871	94.864787	-74.458618
2065	898.500779	38.900382	898.500779	898.500779	94.334992	-58.014461
2070	1146.739979	49.647840	1146.739979	1146.739979	93.499760	-45.053387

---

Technical Lifetime = 1

This scenario sets the technical lifetime to 1, forcing per-period CAP_NEW to equal DEMAND.

9. Mixed 5 and 10-year periods (5%)

• Years: [2020, 2025, 2030, 2035, 2040, 2050, 2060, 2070]
• Interest Rate: 0.05

Year	ACT	CAP_NEW	CAP	Demand	Price	Dual_COMM_BAL_AUX
2020	100.000000	100.000000	100.000000	100.000000	850.000000	-4696.786563
2025	127.628156	127.628156	127.628156	127.628156	850.000000	-3680.055170
2030	162.889463	162.889463	162.889463	162.889463	850.000000	-2883.419520
2035	207.892818	207.892818	207.892818	207.892818	850.000000	-2259.234643
2040	265.329771	265.329771	265.329771	265.329771	850.000000	-1770.169459
2050	432.194238	432.194238	432.194238	432.194238	850.000000	-2473.704582
2060	703.998871	703.998871	703.998871	703.998871	850.000000	-1518.640028
2070	1146.739979	1146.739979	1146.739979	1146.739979	850.000000	-932.313241

10. All 5-year periods (5%)

---

• Years: [2020, 2025, 2030, 2035, 2040, 2045, 2050, 2055, 2060, 2065, 2070]
• Interest Rate: 0.05

Year	ACT	CAP_NEW	CAP	Demand	Price	Dual_COMM_BAL_AUX
2020	100.000000	100.000000	100.000000	100.000000	850.000000	-4696.786563
2025	127.628156	127.628156	127.628156	127.628156	850.000000	-3680.055170
2030	162.889463	162.889463	162.889463	162.889463	850.000000	-2883.419520
2035	207.892818	207.892818	207.892818	207.892818	850.000000	-2259.234643
2040	265.329771	265.329771	265.329771	265.329771	850.000000	-1770.169459
2045	338.635494	338.635494	338.635494	338.635494	850.000000	-1386.974090
2050	432.194238	432.194238	432.194238	432.194238	850.000000	-1086.730492
2055	551.601537	551.601537	551.601537	551.601537	850.000000	-851.481776
2060	703.998871	703.998871	703.998871	703.998871	850.000000	-667.158252
2065	898.500779	898.500779	898.500779	898.500779	850.000000	-522.735948
2070	1146.739979	1146.739979	1146.739979	1146.739979	850.000000	-409.577293

[Wegatriespython]

Using the CRF formula :
for r >0
$$\text{CRF}(r, n) = \frac{r (1+r)^n}{(1+r)^n - 1}$$

r = 5%, n = 30,

CRF ~= 0.0654

For Inv cost = 800, annualized cost = 800 * 0.0654 = 52.3

For the same scenario set inv cost = eps (so that membership in inv techs is not lost) and set fixed cost = 52.3

Years: [2020, 2025, 2030, 2035, 2040, 2050, 2060, 2070], Interest Rate: 0.05

Year	ACT	CAP_NEW	CAP	Demand	Price	Dual_COMM_BAL_AUX
2020	100.000000	20.000000	100.000000	100.000000	102.041148	-563.841757
2025	127.628156	5.525631	127.628156	127.628156	102.041148	-441.784771
2030	162.889463	7.052261	162.889463	162.889463	102.041148	-346.149927
2035	207.892818	9.000671	207.892818	207.892818	102.041148	-271.217526
2040	265.329771	46.013693	265.329771	265.329771	102.041148	-483.723554
2050	432.194238	39.693293	432.194238	432.194238	102.041148	-296.964301
2060	703.998871	70.399887	703.998871	703.998871	102.041148	-182.310320
2070	1146.739979	44.274111	1146.739979	1146.739979	102.041148	-111.922722

Notice that the model is now effectively retiring previous capacity every period, and choosing to start fresh. Exploiting an early retirement arbitrage.

Deactivating early retirement we get :

Years: [2020, 2025, 2030, 2035, 2040, 2050, 2060, 2070], Interest Rate: 0.05

Year	ACT	CAP_NEW	CAP	Demand	Price	Dual_COMM_BAL_AUX
2020	100.000000	20.000000	100.000000	100.000000	102.041148	-563.841757
2025	127.628156	5.525631	127.628156	127.628156	102.041148	-441.784770
2030	162.889463	7.052261	162.889463	162.889463	102.041148	-346.149927
2035	207.892818	9.000671	207.892818	207.892818	102.041148	-271.217526
2040	265.329771	11.487391	265.329771	265.329771	102.041148	-483.723554
2050	432.194238	21.686447	432.194238	432.194238	102.041148	-296.964301
2060	703.998871	36.706344	703.998871	703.998871	102.041148	-182.310320
2070	1146.739979	53.409359	1146.739979	1146.739979	102.041148	-111.922722

Which is virtually the same new capacity investment schedule as case 4.

So effectively the issue is resolved by removing the lump sum INV cost from model objective, replacing it with the annuity A = CRF(r,n) * INV, and adding an early retirement penalty = PV sum of the remaining annual payments per CRF.

Pre-liminary results -> the price spike at 2060 vanishes.

[khaeru]

This is a solid, systematic investigation of the issue and a potential solution, and I think worth presenting and discussing during the MESSAGE weekly meeting some week.

So effectively the issue is resolved by removing the lump sum INV cost from model objective, replacing it with the annuity A = CRF(r,n) * INV, and adding an early retirement penalty = PV sum of the remaining annual payments per CRF.

Just to confirm I've understood the diagnosis, here's what I hope is the same thing stated in different terms:

• The objective function value is computed as a sum over the year dimension of COST_NODAL:

  message_ix/message_ix/model/MESSAGE/model_core.gms

  Lines 340 to 341 in 33edac8

  	OBJECTIVE..
  	OBJ =E= SUM( (node,year), df_period(year) * COST_NODAL(node,year) ) ;

• COST_NODAL is defined here:

  message_ix/message_ix/model/MESSAGE/model_core.gms

  Lines 389 to 470 in 33edac8

  	COST_ACCOUNTING_NODAL(node, year)..
  	COST_NODAL(node, year) =E=
  	* resource extraction costs
  	SUM((commodity,grade)$( map_resource(node,commodity,grade,year) ),
  	resource_cost(node,commodity,grade,year) * EXT(node,commodity,grade,year) )
  	* technology capacity investment, maintainance, operational cost
  	+ SUM((tec)$( map_tec(node,tec,year) ),
  	( inv_cost(node,tec,year) * construction_time_factor(node,tec,year)
  	* end_of_horizon_factor(node,tec,year) * CAP_NEW(node,tec,year)
  	+ SUM(vintage$( map_tec_lifetime(node,tec,vintage,year) ),
  	fix_cost(node,tec,vintage,year) * CAP(node,tec,vintage,year) ) )$( inv_tec(tec) )
  	+ SUM((vintage,mode,time)$( map_tec_lifetime(node,tec,vintage,year) AND map_tec_act(node,tec,year,mode,time) ),
  	var_cost(node,tec,vintage,year,mode,time) * ACT(node,tec,vintage,year,mode,time) )
  	)
  	* additional cost terms (penalty) for relaxation of 'soft' dynamic new capacity constraints
  	+ SUM((inv_tec)$( map_tec(node,inv_tec,year) ),
  	SUM((mode,time)$map_tec_act(node,inv_tec,year,mode,time),
  	( ( abs_cost_new_capacity_soft_up(node,inv_tec,year)
  	+ level_cost_new_capacity_soft_up(node,inv_tec,year) * inv_cost(node,inv_tec,year)
  	) * CAP_NEW_UP(node,inv_tec,year) )$( soft_new_capacity_up(node,inv_tec,year) )
  	+ ( ( abs_cost_new_capacity_soft_lo(node,inv_tec,year)
  	+ level_cost_new_capacity_soft_lo(node,inv_tec,year) * inv_cost(node,inv_tec,year)
  	) * CAP_NEW_LO(node,inv_tec,year) )$( soft_new_capacity_lo(node,inv_tec,year) )
  	)
  	)
  	* additional cost terms (penalty) for relaxation of 'soft' dynamic activity constraints
  	+ SUM((tec)$( map_tec(node,tec,year) ),
  	SUM(time$( map_tec_time(node,tec,year,time) ),
  	( ( abs_cost_activity_soft_up(node,tec,year,time)
  	+ level_cost_activity_soft_up(node,tec,year,time) * levelized_cost(node,tec,year,time)
  	) * ACT_UP(node,tec,year,time) )$( soft_activity_up(node,tec,year,time) )
  	+ ( ( abs_cost_activity_soft_lo(node,tec,year,time)
  	+ level_cost_activity_soft_lo(node,tec,year,time) * levelized_cost(node,tec,year,time)
  	) * ACT_LO(node,tec,year,time) )$( soft_activity_lo(node,tec,year,time) )
  	)
  	)
  	* emission taxes (by parent node, type of technology, type of year and type of emission)
  	+ SUM((type_emission,emission,type_tec,type_year)$( emission_scaling(type_emission,emission)
  	AND cat_year(type_year,year) ),
  	emission_scaling(type_emission,emission)
  	* tax_emission(node,type_emission,type_tec,type_year)
  	* EMISS(node,emission,type_tec,year) )
  	* cost terms from land-use model emulator (only includes valid node-land_scenario-year combinations)
  	+ SUM(land_scenario$( land_cost(node,land_scenario,year) ),
  	land_cost(node,land_scenario,year) * LAND(node,land_scenario,year) )
  	* cost terms associated with linear relations
  	+ SUM(relation$( relation_cost(relation,node,year) ),
  	relation_cost(relation,node,year) * REL(relation,node,year) )
  	* implementation of slack variables for constraints to aid in debugging
  	+ SUM((commodity,level,time)$( map_commodity(node,commodity,level,year,time) ), ( 0
  	%SLACK_COMMODITY_EQUIVALENCE% + SLACK_COMMODITY_EQUIVALENCE_UP(node,commodity,level,year,time)
  	%SLACK_COMMODITY_EQUIVALENCE% + SLACK_COMMODITY_EQUIVALENCE_LO(node,commodity,level,year,time)
  	) * 1e6 )
  	+ SUM((tec)$( map_tec(node,tec,year) ), ( 0
  	%SLACK_CAP_NEW_BOUND_UP% + 10 * SLACK_CAP_NEW_BOUND_UP(node,tec,year)
  	%SLACK_CAP_NEW_BOUND_LO% + 10 * SLACK_CAP_NEW_BOUND_LO(node,tec,year)
  	%SLACK_CAP_NEW_DYNAMIC_UP% + 10 * SLACK_CAP_NEW_DYNAMIC_UP(node,tec,year)
  	%SLACK_CAP_NEW_DYNAMIC_LO% + 10 * SLACK_CAP_NEW_DYNAMIC_LO(node,tec,year)
  	%SLACK_CAP_TOTAL_BOUND_UP% + 10 * SLACK_CAP_TOTAL_BOUND_UP(node,tec,year)
  	%SLACK_CAP_TOTAL_BOUND_LO% + 10 * SLACK_CAP_TOTAL_BOUND_LO(node,tec,year)
  	) * ABS( 1000 + inv_cost(node,tec,year) ) )
  	+ SUM((tec,time)$( map_tec_time(node,tec,year,time) ), ( 0
  	%SLACK_ACT_BOUND_UP% + 10 * SUM(mode$( map_tec_act(node,tec,year,mode,time) ), SLACK_ACT_BOUND_UP(node,tec,year,mode,time) )
  	%SLACK_ACT_BOUND_LO% + 10 * SUM(mode$( map_tec_act(node,tec,year,mode,time) ), SLACK_ACT_BOUND_LO(node,tec,year,mode,time) )
  	%SLACK_ACT_DYNAMIC_UP% + 10 * SLACK_ACT_DYNAMIC_UP(node,tec,year,time)
  	%SLACK_ACT_DYNAMIC_LO% + 10 * SLACK_ACT_DYNAMIC_LO(node,tec,year,time)
  	) * ( 1e8
  	+ ABS( SUM(mode$map_tec_act(node,tec,year,mode,time), var_cost(node,tec,year,year,mode,time) ) )
  	+ fix_cost(node,tec,year,year) ) )
  	+ SUM(land_scenario, 0
  	%SLACK_LAND_SCEN_UP% + 1e6 * SLACK_LAND_SCEN_UP(node,land_scenario,year)
  	%SLACK_LAND_SCEN_LO% + 1e6 * SLACK_LAND_SCEN_LO(node,land_scenario,year)
  	)
  	+ SUM(land_type, 0
  	%SLACK_LAND_TYPE_UP% + 1e6 * SLACK_LAND_TYPE_UP(node,year,land_type)
  	%SLACK_LAND_TYPE_LO% + 1e6 * SLACK_LAND_TYPE_LO(node,year,land_type)
  	)
  	+ SUM((relation), 0
  	%SLACK_RELATION_BOUND_UP% + 1e6 * SLACK_RELATION_BOUND_UP(relation,node,year)$( is_relation_upper(relation,node,year) )
  	%SLACK_RELATION_BOUND_LO% + 1e6 * SLACK_RELATION_BOUND_LO(relation,node,year)$( is_relation_lower(relation,node,year) )
  	)
  	;

• By "the [current] lump sum inv_cost [in the] model objective", you refer to the fact that inv_cost here is multiplied by CAP_NEW in the same/single period in which a particular technology vintage is constructed (the "vintage period", index year).
• If meeting demand in later periods requires increasing a single value CAP_NEW(node, tec, year), this affects OBJ via an increase in COST_NODAL(node, year) for that vintage period. Conversely, there is no connection at all between CAP_NEW(…, year)/inv_cost(…, year) and future periods' COST_NODAL(…, year+1) etc.
• The marginal of COMMODITY_BALANCE_AUX.M(node,commodity,level,year,time) (used to compute PRICE_COMMODITY, same dimensionality) reflects the potential reduction in OBJ via that COST_NODAL(…, year), only for that vintage period (year) and not for future periods (year+1 et seq) in which output from that same technology vintage will also be used.
• Thus the prices in those future periods are artificially low, and in the vintage period artificially high.

[Wegatriespython]

I'm not sure, so the picture as I see it is at the optimal solution the zero profit condition holds and the model balances all its costs across the horizon(allowing our use of crf). This implicitly requires a capital recovery schedule. However since cap new is a lump sum via the objective and the objective is on an uneven grid, the model's capital recovery schedule inherits the approximation error of the grid.

So for instance if we set the resolution at 1 year and same everything else we get an annual price series closer to the calculated CRF value of 102.

Now I can't for sure say why there is a that small difference between the model's endogenous calculation and 102. But the idea is we are giving the model a nicer pre computed method for investment cost annualization than its actual grid. This is what removes the artefact.

Years: [2020, 2021, 2022, 2023, 2024, 2025, 2026, 2027, 2028, 2029, 2030, 2031, 2032, 2033, 2034, 2035, 2036, 2037, 2038, 2039, 2040, 2041, 2042, 2043, 2044, 2045, 2046, 2047, 2048, 2049, 2050, 2051, 2052, 2053, 2054, 2055, 2056, 2057, 2058, 2059, 2060, 2061, 2062, 2063, 2064, 2065, 2066, 2067, 2068, 2069, 2070], Interest Rate: 0.05

Year	ACT	CAP_NEW	CAP	Demand	Price	Dual_COMM_BAL_AUX
2020	100.000000	100.000000	100.000000	100.000000	99.632543	-99.632543
2021	105.000000	5.000000	105.000000	105.000000	99.621985	-94.878081
2022	110.250000	5.250000	110.250000	110.250000	99.610387	-90.349558
2023	115.762500	5.512500	115.762500	115.762500	99.597649	-86.036194
2024	121.550625	5.788125	121.550625	121.550625	99.583662	-81.927725
2025	127.628156	6.077531	127.628156	127.628156	99.568305	-78.014372
..
2067	990.597109	58.085664	990.597109	990.597109	98.061152	-9.899196
2068	1040.126965	60.989947	1040.126965	1040.126965	97.842668	-9.406800
2069	1092.133313	64.039444	1092.133313	1092.133313	97.603731	-8.936980
2070	1146.739979	67.241417	1146.739979	1146.739979	97.342524	-8.488631