CLM Canopy Air Space Project

At last week's meeting we decided that I repeat the tests that I conducted last summer with the "toy model." I saved last summer's work at /home/bluesky/slevis/sib3/oldset. I started the new tests at /home/bluesky/slevis/sib3. The README describes the tests. I have written up the complete derivation of the matrix coefficients at /home/slevis/papers/clm4/pcas1.doc.
Sam (slevis at 11 am on Jan 31, 2006)
– -
Today Gordon, Mariana, Keith, Forrest, and Sam went over eqs. 1 and 2 (sections 2.1 and 2.2 of pcas1.doc).
Questions that arose:
1) Solve for Dtheta and Dq (as sib does for the coupled case) or use sib's offline method (ie set Dtheta=0 and Dq=0) or drop the equations altogether?
2) Where does the atmosphere forcing appear in the matrix equations?
3) How should we initialize variables in arbitrary initialization?
Gordon, Keith, and Sam continued the conversation in the afternoon. Our answers:
1) Do not solve for Dtheta and Dq. Ask Reto Stockli at the LMWG meeting whether we need sib's offline method (that sets Dtheta=0 and Dq=0) or we can drop the first two rows and columns altogether.
2) The atm forcing appears in Hn (sens heat flux at timestep n).
3) We will discuss initialization at the next meeting. I think that we have what we need to initialize variables at timestep 1.
At the next meeting we will discuss section 2.3 of pcas1.doc. Next meeting on 2/23 probably at 10 am.
Sam (slevis 3 pm 2/13/06)
– -
On 2/23/06 we decided that we needed to know more clearly which variables were involved in the right hand side (RHS) of the matrix. So we decided to expand the heat flux terms on the RHS of all equations. Also we decided to remove from the matrix the equations that were solving for Dtheta_atm and Dq_atm.
I have completed both tasks twice: Once in the toy model (details in /home/bluesky/slevis/sib3/README) and once in the technical description (now at /home/slevis/papers/clm4/pcas2.doc). I have saved the document in pdf (see attachment at the bottom of this entry).
Assignment: At the next meeting (3/6 at 10 am) we will go over my changes. See sections 2.1 and 2.2.
Sam (slevis 5 pm 3/2/06)
– -
On 3/6/06 we decided on a few additional changes to the tech note:
1) Add a bar over theta_atm
2) Add a subscript j in fwet and fdry and
3) Write the terms C_Tv_j in an easier way to understand. In the end I simply removed the confusing brackets and added the words "for j=1,2,...,npft" next to the relevant equations.
We also decided to add a comment in the "toy model" that would remind us for the bare ground case that
1) (rb .ne. 0.) and
2) Now we have both rah1 and rah2 defined
I have made the changes listed above. This probably concludes the review of eqation 1, so the assignment for next Monday's meeting is to go over eqation 2 for discussion.
Sam (slevis 5 pm 3/7/06)
– -
On 3/13/06 we went over eq. 2. We discussed reminders for the bare ground case that should be entered in the 'toy model' and I have added those:

• (rb+rs .ne. 0)
• Both raw1 and raw2 must be defined
• I also changed the expression fdry/L to fdry/max(L,1.e-6).
  I added the following comments to the tech note:
  On 3/15/06 I met with Reto:
• He will send me the code that includes the water and energy limits. These limits are applied before solving the matrix.
• He offered to review the tech note after we finish reviewing it ourselves.
• He recommended leaving DTsoi in the tridiagonal matrix for computational efficiency.
  Assignment: Go over eq. 3 (section 2.3). Generally straightfwd, but I have left a couple of questions unanswered that one of you may help answer.
  Sam (slevis 11 am 3/16/06)
  • -
    On 3/20/06 we went over eq. 3 and eq. 4 and decided that I had chosen the wrong derivatives of longwave radiation. I have corrected these in the tech note and the toy model.
    Assignment: go over eq. 4 in more detail. Next meeting on April 3rd.
    We liked Reto's suggestions (see last week's entry).
    Sam (slevis 2 pm 3/22/06)
  • -
    On 4/17 we went over eq. 4 that solves for Tg. We decided that we would solve for Ts, qs, Tv, and Tg=Tsoisno(top) with the PCAS matrix and for remaining soil/snow temperatures with clm's existing tridiagonal matrix. Keith and I confirmed that clm reduces the heat capacity of the top layer by a fudge factor and that this factor is applied consistently through the code to conserve energy. We will keep this approach in solving for Tg in the PCAS matrix. Tg will now be an input rather than an output in clm's tridiagonal matrix.

For reference:

• In Sib3 Ian Baker solves for all the soil temperatures in the PCAS matrix. He considers Tg=Tsoisno(top), but he does not reduce the heat capacity of the top layer as we do in clm.
• In Sib2.5 Reto Stockli solves for soil temperatures separately the way we plan to, but he calculates separate Tsoisno(top) and Tg.
  I have created new versions of the tech note and the 'toy model' to reflect our decisions.
  Next Monday we will finalize these issues and in two weeks we will discuss what should come next.
  Sam (slevis 4 pm 4/17/06)
  • -
    Last Monday we clarified a few details about how:
    1) the pcas matrix would solve for Tg
    2) the tridiagonal matrix would solve for 'deep' soil temperatures
    3) the soil/snow water phase change would not affect the sensible and latent heat fluxes as it does in clm3
    Particularly for (1) and (2) we clarified that F(1+snl) at timestep n+1 is the lower boundary of the pcas matrix. The same term is then used as the upper boundary of the 9-row soil/snow temperature matrix.
    Gordon suggested that I include explicit definitions of the sensible, latent, and soil/snow heat flux terms in the tech note.
    I have made these changes. The document is still named pcas3.pdf. These changes did not require modification of the 'toy model.'
    Once we consider these issues resolved, we will discuss implementation of Reto's scheme of limiting the latent heat fluxes according to water availability.
    Sam (slevis 5 pm 4/28/06)
  • -
    On May 1st we decided to implement the water limits found in Reto's paper and in his code (sib2.5). I have done this in the "toy model." I have also written it all up in the tech note. At the next meeting we will discuss what I've done.
    Sam (slevis 11:30 am 5/12/06)
  • -
    I have had informal meetings with Gordon and Keith since May 1st. I have corrected the limits (discussed in the previous entry) by dividing them by dtime (the clm timestep) and by adding soil ice to the soil water. I have also modified the presentation of these limits in the tech note.
    Sam (slevis 3:15 pm 6/1/06)
  • -
    On June 5th I refamiliarized everyone with the latest developments (see previous entry). In limiting ground evap by the available soil moisture, we decided to normalize the pft-level values as in Oleson et al (2004), p. 76. The tech note is up-to-date. The 'toy model' is not because I received a few emails from Reto and P.L. Vidale suggesting a new way of applying the limits to evaporation. I have decided to hold off until we discuss this issue in Breckenridge.
    Sam (slevis 4:45 pm 6/13/06)
  • -
    Reto, Gordon, Keith, Forrest, Mariana, and I met in Breckenridge on June 20th and decided the following:
• Take out mods that pertain to limiting evapotranspiration according to Vidale & Stockli (2005).
• Instead, limit canopy evap at time step n+1 to half the max possible canopy evap at time step n. Do this by limiting the clm-calculated value of rb.
  The document pcas4.pdf includes these changes. Before I modify the toy model, I want you to look over my work. I also want to discuss the following questions, which I may do in an informal meeting with Gordon, but feel free to tell me your opinion:
• Should the limit to rb apply to the LHS of the equations only? I assume so, since the RHS represents what actually happened in time step n. If so, I will need to declare a new rb (say, rblim) that represents the limited version of rb.
• Should rblim apply to canopy evaporation only (Eq. 2c and 3c) or canopy sensible heat flux as well (Eq. 1c)?
  Sam (slevis 5 pm 6/26/06)
  • -
    Gordon and I answered the above questions in an informal meeting:
• Do not modify rb because it is also used in transpiration and sensible heat flux calculations
• Make a coefficient beta that equals the ratio of maximum possible canopy evap over actual canopy evap and ranges from 0 to 1
• Apply beta to canopy evaporation in the LHS of the matrix of equations
  I have updated the tech note and the toy model.
  Our next regular meeting will be on July 10th at 10 am.
  Reto has offered to look over the tech note by then, so we should know at that time whether we are ready to move on to the code. "Moving to the code" means going over the toy model first before implementing anything in clm.
  Sam (slevis 3:15 pm 6/28/06)
  • -

Gordon, Keith, Mariana, and I met.
We reviewed the changes mentioned in the last entry. We agreed that the scientific content of the tech note is sound.
We went over the main points in Reto's review of the tech note. We decided that I will reorganize the tech note according to his suggestions.
We will meet on Monday the 17th to discuss the next steps.
Sam (slevis 11:40 am 7/10/06)
– -
Today we assigned Forrest to look at the 'toy model' of the prog. canopy air space equations and compare to the tech note (latest posted version). Forrest will bring questions, concerns, feedback to the meeting next Monday.
Sam (slevis 11:40 am 7/17/06)
– -
After last Monday's meeting, I made the following changes to the toy model (slv16.F90):

• I changed Dz from max(4.,htop-hbot) to max(1.,htop) as Reto has done
• I removed Ca from the code and wrote it as forc_rho x cpair x Dz instead
• I added the terms (sabg-net_lw_grnd) to eq. 4, which I had omitted before
• I removed the 2 x dtime smoothing as Reto suggested
• I removed fracvegnosno from eqs. 3 and 4 where it was redundant. In CLM we will need to confirm that this works. Keith said that we should prob make sure that when L+S less than 0.05 we set L=0 and S=0 in the CLM implementation
• I corrected the usage of hvap and htvp in the toy model relative to lamda_vap and lamda_vapsub
• I made the tech note consistent with slv16.F90
• I reordered the tech note following Reto's advice: (1) Introduction (short) (2) CLM Heat Flux Equations, ie the clm3 diagnostic & iterative scheme (3) The Prognostic Equations (4) Numerical implementation, ie the time differencing scheme (5) CLM Implementation, ie multiple pfts, canopy evaporation limits, etc
  Please browse through pcas5.pdf to familiarize yourselves with the new chapters, but I hope you will also have the time to read it more closely. Forrest offered to compare the equations against the toy model a second time. Also Reto has offered to review the tech note a second time. I will wait for our group's comments first before I ask him to go over the document again.
  Based on our schedules, we will be meeting next on Monday September the 11th. Feel free to contact me in the mean time with comments and questions.
  Sam (slevis 1:30 pm 7/28/06)
  pcas5.pdf
  • -
    In September's meeting we asked Forrest to confirm again that the toy model and the tech note were consistent. In today's meeting Forrest said that he had compared the toy model and the tech note less thoroughly than before and did not find inconsistencies.
    In today's meeting we decided:
• Forrest would create a branch for the implementation of pcas in clm.
• Forrest would begin to introduce new parameters to the model based on what's labeled 'new' in the toy model. Such parameters should probably be defined in clm_varcon. At the meeting we said that some parameters could be defined locally, but now I'm not sure that's a good choice.
• The routines CanopyFluxes and BareGroundFluxes should not calculate temperatures, humidities, and fluxes anymore. However, they would continue to calculate other variables as they do now. From a quick search in BareGroundFluxes, I found no mention of glaciers and wetlands. This suggests that nothing special is done there for these two land units.
• The routine SoilTemperature should be modified to work with 9 instead of 10 soil layers.
  Testing will begin with timestep level and pft level comparisons of the old to the new in single timestep and one day runs.
  I will ask Reto to review the tech note. His name will appear in the list of co-authors.
  Next meeting: Thursday Oct 19th at 10 am.
  Sam (slevis 5:20 pm 10/12/06)
  • -
    Gordon, Forrest, Keith, and Reto met at NCAR on October 19 to discuss implementation of PCAS into CLM3. With regard to new variables, we decided:
• psb = 0.50 (BL mass depth [Pa]) would be added to clm_varcon
• leafh2o = 0.2 (leaf water [kg/m2]) would be added to clm_varcon
• Dz = max(1., htop) [m] would be a local variable
• Cc = (tla + tsai) x cpliq x leafh2o + cpliq x h2ocan [J/m2K] would be a local variable
• net_lw_veg would be a new global variable and output as FIRA_V on history files
• net_lw_gnd would be a new global variable and output as FIRA_G on history files
  We decided that biogeophys/CanopyAirSpaceMod.F90 would be an acceptable name for the new module, even though it includes calculations over bare ground, because the purpose of the module is clear.
  Gordon said he would go back and carefully read the draft tech note. Reto said he would provide comments on the tech note directly to Sam.
  Later the same day, Forrest and Mariana discussed where to branch off the CLM3 development line. There was no reason not to branch off the lastest tag, clm3_expa_74, and it was agreed the new branch name would be pcas_clm3_expa_74.
  On Tuesday, October 24, Forrest created the new branch and the initial branch tag in the NCAR subversion repository.
  Sam has set the next meeting for Monday, October 30 at 10:00 am MT.
  Forrest (forrest Tue Oct 24 23:01:49 EDT 2006)
  • -
    Forrest and I talked for 30 minutes this morning. Forrest asked for and I gave him further clarification of the method for calculating all that is necessary for the prog. canopy matrix of equations. Forrest brought up that restart files may require new information saved and the question was "what do we save?" I suggested that Forrest keep the existing equations in the code taking out the iterations. Then he should be able to trace back to the variables coming from the previous timestep that are used to calculate the Monin-Obukhov length and the resistance terms. It's possible that nothing new will be needed in restart files.
    Before ending the call, I asked him for a time-estimate for the completion of this project. Forrest mentioned that at this stage the coding felt more complex than he had anticipated. He said that once we passed the planning stage of how everything fits together in the CLM, he would begin to make faster progress, still limited by the 1/4-time allocation to this project and by other priorities (eg biogeochem) taking precedence. He and I decided that the June CCSM meeting was a reasonable time by which to have something concrete to show.
    No conference call on Dec 11th due to AGU. Next call on Dec 18th. We will continue in this format, where I call Forrest, and I gather the group if there is a larger issue to discuss.
    Sam (slevis 1.40 pm 12/4/06)