Poor performance of finite difference frozen soil temperature scheme

[tbohn]

Comparison with observations in Tibet indicates very large cold excursions during winter at depths between 1.5 and 3 m, when constant temperature specified at 7m depth.

Preliminary indications are that an incorrect term in the heat equation (in frozen_soil.c) is to blame; a different version of the code by Tara Troy (formerly of Eric Wood's lab) improves performance substantially. However, Tara's code only fixed the explicit scheme, with linear node spacing. This fix needs to be ported to the implicit scheme and to exponential node spacing before it can be officially incorporated into VIC for release. (and of course it should be tested at a sites outside of Tibet).

[tbohn]

Further testing by Lan Cuo on sites in Tibet indicate:

1. Tara's code reduces the discrepancies but not completely
2. Changing the soil moisture contents can fix the cold excursions during winter
3. A potential bug still remains: the middle depths of the soil column tend to become colder than the top and bottom boundaries (in annual average sense). This seems to happen for both implicit and explicit (with Tara Troy's fix) approaches.

Some possible culprits:

1. A while back, I implemented a "fall back" option for when the temperature iteration failed to converge, at either the surface or at one of the nodes. It sets the node temperature to the average of the neighboring nodes' temperatures (I don't remember, I may have also involved the previous step's temperature of the current node). A related fallback is when a node has a more extreme temperature than its neighbors (warmer or colder) and this difference grows over time rather than shrinking over time - this was to fix the "runaway cold nose" problem first encountered by Jenny Adam. Once again, in such cases we set the node temperature to the average of its neighbors. In either case, doing so doesn't conserve energy. So, this could conceivably lead to a long-term decline in temperature.
2. Energy could be more or less conserved (i.e. errors due to 1 are small), but heat could be going into phase changes, thus lowering the temperature. To check this, need to see if ice/water contents are changing over the long term. If this is happening, it would imply that longer spinup is needed.

So, in addition to porting Tara's code to the implicit scheme, this weird behavior needs to be investigated.

[tbohn]

I may have misunderstood what Lan was saying about temperature "drifting lower" - this might have been with depth instead of over time. So, this might not be a bug. I'm trying to get clarification.

[jhamman]

@tbohn , if possible, can you post a figure showing demonstrating the issue (or absence of an issue).

@bartnijssen and I will be interested to in making sure this is all working for the RASM project.

[tbohn]

Sure. Lan's simulations originally had a lot of discrepancies with observations. We've gotten a lot of the issues solved, through increasing the number of soil nodes, switching to Tara Troy's code, and tweaking soil moisture (Lan didn't give details on exactly what she changed about the soil parameters to change soil moisture contents). One example of where it is working well is: (Each line is a different depth in the soil column; unfortunately I don't know exactly which depths Lan is plotting, but the soil column goes to 7 m)

But apparently there are still some sites with problems - including a summertime low bias with depth:

At this site, the problems are not so bad, but there still appears to be a low bias at middle depths (approx 2-5 m depths)

We are now working on more soil moisture tweaks, since the soil moisture from the 3rd site doesn't match observations at middle depths (depths here are 0.4m, 1.4m, and 2.5m): It looks like having too much moisture in the middle layer might be reducing the seasonal temperature amplitude.

So, it may be that VIC is just fine, but to get the right answers you may need to tweak the soil parameters - hopefully it can be reconciled with streamflow calibrations.

[bartnijssen]

It's important that we distinguish between bugs (especially things that make the model crash) and imperfect representations. Seems to me that the current problem falls in the latter category.

On May 27, 2014, at 11:07 AM, Ted Bohn notifications@github.com wrote:

Sure. Lan's simulations originally had a lot of discrepancies with observations. We've gotten a lot of the issues solved, through increasing the number of soil nodes, switching to Tara Troy's code, and tweaking soil moisture (Lan didn't give details on exactly what she changed about the soil parameters to change soil moisture contents). One example of where it is working well is:

(Each line is a different depth in the soil column; unfortunately I don't know exactly which depths Lan is plotting, but the soil column goes to 7 m)

But apparently there are still some sites with problems - including a summertime low bias with depth:

At this site, the problems are not so bad, but there still appears to be a low bias at middle depths (approx 2-5 m depths)

We are now working on more soil moisture tweaks, since the soil moisture from the 3rd site doesn't match observations at middle depths (depths here are 0.4m, 1.4m, and 2.5m):

It looks like having too much moisture in the middle layer might be reducing the seasonal temperature amplitude.

So, it may be that VIC is just fine, but to get the right answers you may need to tweak the soil parameters - hopefully it can be reconciled with streamflow calibrations.

— Reply to this email directly or view it on GitHub.

[tbohn]

Yes. Certainly the port of Tara Troy's code to the implicit scheme seems to be a necessary change - consistent with calling this an enhancement. The other information is in flux - first thought the low bias at middle depths was a bug, then later it appears to be not a bug. I thought the information should go here... Should it go somewhere else?

[bartnijssen]

Yes - it should go here - I just want to make sure that we distinguish between these sources of error. A bug that crashes is critical and it will prevent us from using that code in a coupled model, an improvement is "merely" desirable ;)

On May 27, 2014, at 11:33 AM, Ted Bohn notifications@github.com wrote:

Yes. Certainly the port of Tara Troy's code to the implicit scheme seems to be a necessary change - consistent with calling this an enhancement. The other information is in flux - first thought the low bias at middle depths was a bug, then later it appears to be not a bug. I thought the information should go here... Should it go somewhere else?

— Reply to this email directly or view it on GitHub.

[tbohn]

So, after more testing, it appears that Lan's problems were solved via longer spin-up of at least 5 years (originally only a 1-year spinup was used) and switching from a 7m damping depth with constant temperature boundary condition to 25m with no-flux bottom boundary condition. I had thought that it was common knowledge that a spinup on the order of at least a decade was necessary for frozen soil / permafrost simulations.

Anyway, it looks like the only major to-do item here is to port the Tara Troy code into the implicit option. Perhaps it would also be worthwhile to point users to permafrost-modeling literature and recommend deep damping depths and no-flux boundary conditions.

[jhamman]

@tbhohn - thanks for the update and the good news. So if the major to-do is porting Tara's code, let's outline how/who/when to get it done.

[bartnijssen]

Also which version of the code to integrate it into ...

On May 28, 2014, at 1:34 PM, Joe Hamman notifications@github.com wrote:

@tbhohn - thanks for the update and the good news. So if the major to-do is porting Tara's code, let's outline how/who/when to get it done.

— Reply to this email directly or view it on GitHub.

[jhamman]

Since this is a science improvement and it probably will not make it into 4.2, I would think this is suited for 5.1.

Ted, is this something Lan Cuo is going to do? Maybe we could ask her to contribute.

[tbohn]

I don't know if we should count on her for it - sounds like she's working alone (she didn't have anyone to help her with the testing) and I don't get the sense that she's eager to do a lot of coding.

[jhamman]

An update via personal communication with Lan Cuo on 2014-12-17:

I have been using VIC 4.1.2. model to simulate the soil temperature at various layers down to 5 m below surface. Ted changed deep layer soil temperature simulation using Tara's codes to improve my simulation. Ted, please clarify if questions raised for this part. For some reason, VIC has problems in simulating soil temperature at both surface layer and deep layers... At my sites, soil temperature is basically underestimated at all sites. The first layer is ok during the cold season, but it tends to be underestimated during the warm and wet season at all sites. Deep layers are all underestimated regardless of seasonality.

[jhamman]

To add a bit more about Lan's simulations, here are the key soil temperature parameters from her global parameter file:

NLAYER              3
NODES               50
TIME_STEP           1
FULL_ENERGY         TRUE
FROZEN_SOIL         TRUE
QUICK_FLUX          FALSE
IMPLICIT            FALSE
QUICK_SOLVE         FALSE
NO_FLUX             TRUE
EXP_TRANS           FALSE
GRND_FLUX_TYPE      GF_410
TFALLBACK           TRUE