Ideas for optimisation in water use model

[itowers1]

Just noting some ideas to optimise speed of water use model:

• [x] in ff16w strategy, optimum stem water potential/ci is identified, and the profit is recalculated at this value. Instead, could profit be stored iteratively as optimum position is calculated, and then recalled in FF16w strategy
  • Update: Solved by storing previous value
• [x] recalling optimum stem water potential/ci between nodes/timesteps breaks due to required transpiration stream exceeding what is theoretically possible using Newton's method of differentiation, but is likely associated with the fact that the last node of a given time step has very different environmental conditions to the first node of the next time step once the forest has matured a little bit. Thus, storing the recalled value between nodes, rather than across calculations, could potentially solve this issue and unlock Newtons method which is proposed to be faster than GSS.
  • Update: Solved by switching to GSS which is more stable.
• [x] calculation of certain values at a high level without having to input them as arguments. For example, the ratio of maximum profits to maximum costs (lambda in sperry's method) is constant for a given node in a given timestep, meaning that it need only be calculated once for a given call to the leaf model. However, it seems unnecessary to have this as another argument akin to k_l_max in calc_profit, thus I am proposing some other way of having this as a pre-computed variable in ff16w strategy which can be passed to leaf model.
  • Update: No longer relevant

[aornugent]

Awesome - looks like you're really getting to the finer points. We discussed some potential solutions today and I think there's a way to achieve all three and significantly clean up the code.

Here's a sketch of some of the core ideas.

The environment and vars objects are passed to the leaf as pointers. set_physiology updates the leaf state, pre-calculating useful values as you suggest in point 3. When we optimise psi_stem the output is saved in the vars object for access in the next timestep as you suggest in point 2. @dfalster described this approach as 'handing over your checkbook' - what this means is that the side-effects of leaf.optimise_psi_stem is not known from FF16w_Strategy - we just trust that it does what we want. optimise_psi_stem then sets profit in the leaf rather than returning the value as you suggest in point 1.

@dfalster @devmitch - any major drawbacks for optimise_psi_stem having the side-effect of setting vars? I like the clean code but don't have a strong sense of modern C++ style.

FF16w_Strategy.cpp

FF16w_Strategy::net_mass_production_dt(...) {
   ...

   // update leaf to reflect env. conditions
   leaf.set_physiology(environment, height);

   // track stem hydrology as individual auxiliary state
   leaf.optimise_psi_stem(vars);

   double assimilation = leaf.profit * area_leaf_; 
   // or
   // double assimilation = leaf.assimilation_per_area * area_leaf_;

   ...
}

FF16w_Strategy.h

std::vector<std::string> aux_names() {
  std::vector<std::string> ret(
      {"competition_effect", "net_mass_production_dt",  "opt_ci"});
  // add the associated computation to compute_rates and compute there
  if (collect_all_auxiliary) {
    ret.push_back("area_sapwood");
  }
  return ret;
}

leaf_model.h

class Leaf {
   ...
   void set_physiology(const FF16_Environment &environment, double height);
   void optimise_psi_stem(Internals &vars);
}

leaf_model.cpp

// set up leaf internal state based on environmental conditions
void Leaf::set_physiology(const FF16_Environment &environment, double height) {
   psi_soil = environment.get_psi_soil() / 1000000;
   PPFD = environment.PPFD;
   huber_value = ...;
   k_l_max = K_s * huber_value / (height * eta_c);
   lambda = ...;
}

void Leaf::optimise_psi_stem(Internals &vars) {
   // start search at prev. node optima
   double psi_stem_initial = vars.aux(aux_index.at("opt_ci"));

   // calculate profit at initial psi_stem
   y_0 = calc_profit(PPFD, psi_soil, psi_stem_initial, k_l_max)

   // optimise
   ...

   // set new node optima for next timestep
   vars.set_aux(aux_index.at("opt_ci"), psi_stem_next)  // might have variable names confused.

   // set profit at optimal psi_stem
   leaf.profit = y_0;

  // return nothing because void 
}

[itowers1]

Looks good, thanks for writing that up Andrew. I'll address these now and close when ready.

[a-orn]

I suppose it wouldn't be much extra code to avoid passing vars and setting opt_ci as a side-effect. This has the benefit of being able to test from R without having to create an Individual or Internals object.

FF16w_Strategy.cpp

FF16w_Strategy::net_mass_production_dt(...) {
   // update leaf to reflect env. conditions
   leaf.set_physiology(environment, height);

   // optimise for profit
   leaf.optimise_psi_stem(vars.aux(aux_index.at("opt_ci")));

   // track stem hydrology as individual auxiliary state
   vars.set_aux(aux_index.at("opt_ci"), leaf.opt_ci)

   double assimilation = leaf.profit * area_leaf_; 
   // or
   // double assimilation = leaf.assimilation_per_area * area_leaf_;
}

leaf_model.cpp

void Leaf::optimise_psi_stem(double psi_stem_initial) {
   // optimise for profit starting at psi_stem_initial aka. opt_ci at prev. timestep
   y_0 = calc_profit(PPFD, psi_soil, psi_stem_initial, k_l_max)

   // set new node optima for next timestep
   opt_ci = psi_stem_optimised; // might have variable names confused.

   // set profit at optimal psi_stem
   profit = y_0;

  // return nothing because void 
}

[dfalster]

Thanks for the suggested code Andrew. Yes, I agree with not passing in the vars variable, as this keeps clearer boundary between the different parts. Isaac and I caught up about that and I explained a bit more about vars.

[dfalster]

Closed as these have been addressed