Change renormalization of landunit areas in mksurfdata_map to be consistent with raw datasets

[olyson]

Description of changes

This PR is based on code branched off of ctsm5.1.dev062.
Per Issue #1555 the following problems are corrected:

1. Urban no longer adjusts bare ground or PCT_NAT_PFT in general. Urban is now treated the same as other special landunits.
2. The areas of all landunits (including crop) other than natural veg are now treated symmetrically; in places where we assure sums don't exceed 100% and in places where we adjust all other landunits proportionally. This is accomplished by replacing the suma rescaling code that is done just before calling normalizencheck_landuse. This rewritten section of code is also inserted in the time loop for transient landcover, also just before calling normalizencheck_landuse.
3. For determining landunit areas, PCT_NATVEG is now set to 100 minus the sum of the other landunits. This replaces the current rescaling code in normalizencheck_landuse.

Specific notes

There is an existing subroutine called adjust_total_veg_area that adjusts the total vegetated area (natural veg & crop) to a new specified total, which is no longer appropriate since crops are treated symmetrically as with other landunits. Instead of modifying this subroutine, I created a new one, adjust_total_natveg_area that adjusts only the natural vegetated area to a new specified total. I haven’t removed the old subroutine for now even though it is no longer used.
There are two errors checks I’ve left in (marked with “!KO”) that may no longer be useful or necessary or correct. It may be that they need to be changed with respect to the crop landunit, but I’m not quite sure what the intended purpose of these are. I’ve left these in unmodified for now since they are simply error checks and don’t actually change any of the landcover variables. And they don’t seem to be triggered in my testing.

Contributors other than yourself, if any: Bill Sacks

CTSM Issues Fixed (include github issue #):
Resolves #1555

Are answers expected to change (and if so in what way)? Yes, surface datasets will be different

Any User Interface Changes (namelist or namelist defaults changes)? No

Testing performed, if any:

For initial testing, I’ve created a 0.9x1.25 surface dataset and landuse file from ctsm5.1.dev062 and then the same set of files from the modified mksurfdata_map. The new surface dataset has 7 variables that are different: PCT_NATVEG, PCT_CROP, PCT_NAT_PFT, PCT_CFT, PCT_LAKE, PCT_GLACIER, PCT_URBAN. The MAXimum (increase) and MINimum (decrease) of these differences are:

(0) MAX pct_natveg: 2.842170943040401e-14
(0) MIN pct_natveg: -10.67630993642207
(0) MAX pct_crop: 10.67630993642206
(0) MIN pct_crop: -1.4210854715202e-14
(0) MAX pct_lake: 1.55750967678614e-11
(0) MIN pct_lake: -0.2633625082721949
(0) MAX pct_glacier: 4.760636329592671e-13
(0) MIN pct_glacier: -7.105427357601002e-15
(0) MAX pct_urban: 6.350475700855895e-13
(0) MIN pct_urban: -0.03106996610432056
(0) MAX pct_wetland: 0
(0) MIN pct_wetland: 0
(0) MAX pct_nat_pft: 90.21452585120576
(0) MIN pct_nat_pft: -52.42969544487663
(0) MAX pct_cft: 2.81071874515417e-09
(0) MIN pct_cft: -9.416396551387152e-10
The largest differences are in PCT_NAT_PFT. The largest increase is associated with bare soil because urban is no longer replacing bare soil preferentially. Glacier has roundoff differences. Lake and Urban show a small decrease, presumably because Crop is now included in the initial normalization. The largest increase is crop is equivalent to the largest decrease in natural vegetation. I might not expect PCT_CFT to change but the largest changes are on the order of e-09 (roundoff?).
The new landuse timeseries has 6 variables that are different: PCT_CROP_MAX, PCT_NAT_PFT_MAX, PCT_CFT_MAX, PCT_CROP, PCT_NAT_PFT, PCT_CFT. The MAXimum (increase) and MINimum (decrease) of these differences are:
(0) MAX pct_crop: 16.4328269096183
(0) MIN pct_crop: -5.684341886080801e-14
(0) MAX pct_crop_max: 16.4328269096183
(0) MIN pct_crop_max: -2.842170943040401e-14
(0) MAX pct_nat_pft: 98.50930869418386
(0) MIN pct_nat_pft: -77.092800323726
(0) MAX pct_nat_pft_max: 90.6516908188529
(0) MIN pct_nat_pft_max: -77.092800323726
(0) MAX pct_cft: 0.2290990088236286
(0) MIN pct_cft: -0.2965845809809196
(0) MAX pct_cft_max: 0.172655463450238
(0) MIN pct_cft_max: -0.2965845427981293

I ran a short historical simulation with the new files successfully.
Are there any mksurdata unit tests or other tests that should be run to check functionality?

[billsacks]

As @ekluzek pointed out at today's ctsm software meeting, we should wait to merge this until we're ready to make new surface datasets, so that mksurfdata_map stays more in sync with our out-of-the-box surface datasets. I'll still go ahead and review it.

[billsacks]

Thanks a lot for working on this, and for documenting all of the testing you have done! I'm sorry it has taken me a couple of weeks to finish reviewing it.

I have some comments below. They are mostly pretty minor, but there is one more substantial one relating to the possible need for merging this with #1579 and #1073, because I'm now seeing that it's hard to bring this one in without some of the changes from #1579 (specifically related to the storing of the "_orig" variables).

[billsacks]

• I think I understand why you added this here (because now you potentially need to redo the renormalization each year, since crop areas may have changed). However, I think that doing it right requires the addition of "_orig" variables that are used to reset non-changing landunit areas each year, as is done now by @keerzhang1 in Enable the mksurfdata_map to generate landuse_timeseries for dynamic urban and lake #1579 . I also like the solution in Enable the mksurfdata_map to generate landuse_timeseries for dynamic urban and lake #1579 of moving this code into "normalizencheck_landuse" to avoid duplication.

For this reason and maybe others, I'm starting to think that it's not going to work well to do my original suggestion of keeping this PR separate from #1579 , and instead we should do a single set of changes that combines the changes from this PR with #1579 as well as those in #1073, since there is so much overlap between these three sets of changes. @olyson and @keerzhang1 what do you think about that? Alternatively, if we want to keep the behavior change in this PR separate from the feature additions of #1579 and #1073, I think the best way to do it would be to (1) combine #1579 with #1073, then (2) apply the changes in this PR in a new PR in their own (but probably new) PR on top of (1). It comes down to how much we care about keeping feature additions separate from behavior/answer changes in mksurfdata_map; @ekluzek do you have any feeling about how important it is to keep these separate?

[olyson]

After reviewing the changes made by @keerzhang1 last week, I was beginning to come to the same conclusion. My preferred approach would be to combine these three PRs into one PR (#1579 ), rather than in two steps. But will leave it up to you and Erik to decide. As you've noted, there is also a fourth PR forthcoming from issue #1572 , code that reads in HASURBAN and initializes landunits in memory (which will now need to be changed to PCT_URBAN_MAX I guess). I think that would come in as part of the feature additions if we separated out feature from behavior.

[keerzhang1]

I also like the idea to combine these three PRs into one PR (#1579) by the two steps proposed by Bill. I think perhaps I can combine #1073 with #1579 first. Then after we finalize #1564, we can apply the changes in #1564 to #1579 (or a new PR) or keep them separate.

[olyson]

Thanks @keerzhang1 for offering to combine #1073 with #1579 , that sounds good to me. Let us know if you have any questions. I'll work on responding to the review of #1564

[billsacks]

Thanks a lot @keerzhang1 ! I'm thinking it will be easiest if you use #1579 as a starting point, since that's already in good shape, and then bring in the relevant changes from #1073 into that. Once we look that over and are happy with that merge, then @olyson can bring his changes from this PR (#1564 ) on top of that.

A few notes when integrating the changes from #1073

• There's no need to bring in the namelist changes from that PR
• For the text file that connects the perl script to the mksurfdata_map fortran code, you can use the approach you implemented for urban: using the same text file as for transient PFTs but adding an extra line for each year.
• You don't need to update the perl code to add these lake files for now, unless that will help you with your testing.

Please let us know if you have questions or run into any problems.

[billsacks]

Closing this, since it has been superseded by #1587