Your interactions of one 3-level categorical predictor with one 2-level categorical predictor (one combination of categories missing) and with a continuous predictor limits you to an intercept plus 9 coefficients. That has nothing to do with this being a mixed model. However you end up coding this, the results should all be identical in terms of outcome predictions at any set of predictor values when all of the interaction terms are taken into account. For example:

predict(mod,data.frame(fire_trt="re",fert_trt="fert",yr=3,subplot_id="41-1",plot="41"))
        1 
0.6235175 

predict(mod2,data.frame(fire_trt="re",fert_trt="fert",yr=3,subplot_id="41-1",plot="41"))
        1 
0.6235175

With respect to Option 3, creating a single 5-level categorical predictor, you say:

won't this change the lower order interactions and main effects...?

When there are interactions you should be wary of how you interpret regression coefficients for "main effects" or lower order interactions anyway. Those are the values only when all of the interacting predictors at higher levels are at 0 (continuous predictors) or at reference levels (categorical predictors). In your case, re-centering the yr values would change the intercept and the coefficients and some lower-level interactions (those that omit yr) of the other predictors. Try this to see:

mod3 <- lmer(biomass ~ fire_trt*fert_trt*I(yr-2.5) + (1|plot) (1|plot:subplot_id), data = my.df)

So don't be put off by that fact alone.

I might lean toward Option 3 but I don't see a harm in using mod2 if you wish. Be very careful in how you do your model comparisons, however. Your mod2 replaces the yr coefficient in mod1 with a numerically identical fire_trtcon:yr interaction coefficient. It might be safest to specify each of the "main effects" and interaction terms as defined in mod2 explicitly. That way, if you try to compare models, you can specify what models to compare rather than depend on how R happens to remove columns to fix rank deficiency.

Finally, if you want to compare mixed models having different fixed effects with likelihood-based metrics (likelihood-ratio tests, AIC), you should not be using the default REML option for fitting the model. See this page. Your random-effect estimates will change, but those don't seem to be of primary interest here.

Your interactions of one 3-level categorical predictor with one 2-level categorical predictor (one combination of categories missing) and with a continuous predictor limits you to an intercept plus 9 coefficients. That has nothing to do with this being a mixed model. However you end up coding this, the results should all be identical in terms of outcome predictions at any set of predictor values when all of the interaction terms are taken into account. For example:

predict(mod,data.frame(fire_trt="re",fert_trt="fert",yr=3,subplot_id="41-1",plot="41"))
        1 
0.6235175 

predict(mod2,data.frame(fire_trt="re",fert_trt="fert",yr=3,subplot_id="41-1",plot="41"))
        1 
0.6235175

With respect to Option 3, creating a single 5-level categorical predictor, you say:

won't this change the lower order interactions and main effects...?

When there are interactions you should be wary of how you interpret regression coefficients for "main effects" or lower order interactions anyway. Those are the values only when all of the interacting predictors at higher levels are at 0 (continuous predictors) or at reference levels (categorical predictors). In your case, re-centering the yr values would change the intercept and the coefficients and some lower-level interactions (those that omit yr) of the other predictors. Try this to see:

mod3 <- lmer(biomass ~ fire_trt*fert_trt*I(yr-2.5) + (1|plot) (1|plot:subplot_id), data = my.df)

So don't be put off by that fact alone.

I would lean toward Option 3, a 5-level categorical predictor, with post hoc tests of particular coefficient combinations of interest. Russ Lenth illustrates that approach with your data in another answer.

I don't see a harm in using mod2 if you wish. Be very careful in how you do your model comparisons, however. Your mod2 replaces the yr coefficient in mod1 with a numerically identical fire_trtcon:yr interaction coefficient. It might be safest to specify each of the "main effects" and interaction terms as defined in mod2 explicitly. That way, if you try to compare models, you can specify what models to compare rather than depend on how R happens to remove columns to fix rank deficiency.

Finally, if you want to compare mixed models having different fixed effects via likelihood-based metrics (likelihood-ratio tests on nested models, or the AIC as you indicate in the question), you should not be using the default REML option for fitting the models you compare. See this page. As Russ Lenth notes in a comment, once model comparisons are done it would be best to fit the final model with REML for further analysis.