Differentiating the operating point

One surface, two consumers

The package differentiates the operating point one way, and two kinds of caller read it.

gradient::at is a thin reader of the surface: it seats the parameters, solves, makes one pass through the three calls below, and reads the requested lanes. Its status is the kind, coarsened:

ProfitMax with no objective (no light, or no stream to normalise over) is no-gradient whatever its kind. prescribed, refused and no-gradient are the statuses of the path that takes psi: refused means the collar fell outside the feasible interval, and no-gradient that the leaf potential it implies lies in the stem curve’s tail, where every row would be rounding. The rest of this vignette is about the surface itself.

Terms

• row: one output’s derivatives against a list of inputs.
• supplied row: a row computed here and handed to the caller as data, with nothing recorded behind it. The caller’s tape gets one statement per output, however long the leaf’s arithmetic is.
• kind: Leaf::OperatingPointKind, which condition placed the operating point. There are eight.
• draw: SupplyDraw<S>, the soil’s water draw at one collar potential: the total uptake E_up and its slope in the collar, and each layer’s uptake. It records the collar it was taken at.
• the three scalars: the leaf’s derivative surface is a set of member templates on a scalar S. At double they are an ordinary evaluation; at odelia::ode::tangent_scalar<double> they carry a directional derivative; at odelia::ode::adjoint_scalar<double> they record onto the caller’s reverse tape. One source, not three copies.
• the pack: leaf_pars<S>, the differentiable parameters as one array of n_leaf_pars indexed by par_*. All but the last are what a calibration fits (gradient::n_theta); the last, par_PPFD, is the light the leaf stands in, seated per observation.

The three calls

After the ordinary solve at double, a consumer at scalar S makes three calls. K is the cost curve the solve used.

l.find_root_collar_psi();                                   // the solve, at double

phylloptim::leaf_pars<S> pars = /* the seated traits, at S, with directions set as needed */;
phylloptim::SupplyInputs<S> supply{psi_soil, r_R_H_min, r_R_V_sum,
                                   pars[phylloptim::par_root_P50],
                                   pars[phylloptim::par_root_c]};

auto draw = l.supply_draw_at<S>(S(l.opt_root_psi_), supply);
phylloptim::Leaf::HeldCoords<S> held;
S    psi  = l.operating_collar<K, S>(draw, pars, curvature, &held);
auto out  = l.billed_outputs_at<K, S>(psi, draw, pars, &held);
// out.profit, out.carbon, out.uptake[j] and out.residual_flux carry their
// derivatives in pars and supply; physical_outputs_at returns the same and
// out.A, out.E, out.psi_leaf, out.ci and out.gc

curvature is l.profit_curvature_at_collar<K>() at an interior point and is not read at any other kind. profit_row_at_current_point in tests/cpp/test_leaf.cpp is this sequence in full, compiled.

profit is the objective that placed the point. carbon is what a plant grows on: A less the part of the curve’s cost that is carbon spent, which each curve’s Curve<K> states (decision D6). Only TF24’s hydraulic cost is carbon, so carbon is profit on TF24 and A on every curve whose cost is a price or a risk term. physical_outputs_at adds the point’s own outputs at the live collar, for a model fitted to what a leaf is observed to do; it costs more, because on a curve whose price reads only the flux (TF24, TF24_floor, CF77) billed_outputs_at bills carbon from profit and the price at the soil’s flux, and never forms the collar’s channel at an interior point.

held is optional. At an interior point, operating_collar solves the leaf potential and the intercellular CO2 at the collar on its way to the marginal, and held carries them to billed_outputs_at, which otherwise solves them again. Passing it moves no number: the forward tangent is bit-identical and the reverse rows agree to rounding. held lives on the caller’s tape, like the draw, so it must not outlive the recording it was made on. At double and at every other kind operating_collar leaves it unplaced and billed_outputs_at solves its own.

The calls assume, and mostly do not check:

• the draw is taken at the solved collar (operating_collar refuses another);
• held, if passed, comes from operating_collar at this point (billed_outputs_at refuses one placed at another collar or kind);
• nothing moves the leaf between the solve and billed_outputs_at;
• pars and supply hold, in value, the leaf’s seated traits (passive_pars()) and its seated soil state and root network;
• supply has the shape Leaf::held_supply() gives on the leaf’s supply path: each layer’s potential and resistances on the multi-layer path, or the one soil potential and the series resistance (r_R_V_sum[0], the pack’s par_resistance) on the single-potential path;
• on ProfitMax, whose normalisers draw the supply again at its reference collar, the draw keeps its own copy of supply for that; no other curve’s does;
• with the thermal-damage layer on, its inputs are the pack’s slots par_A to par_c_repair_flux: the capacities on the pack apply its shifted optimum, its functional fraction and its maintenance costs, as the solve’s members do.

With the energy balance on, the surface follows it on every curve: leaf temperature moves with the stem’s flux and with the light that heats the leaf, and the capacities, the deficit, the cuticle and the first-order condition’s temperature term carry its rows. ProfitMax’s |A|max, which can then peak inside the stream rather than at psi_crit, takes its rows at that peak.

The solve path, and the eight kinds

The solve runs at double and is never differentiated. It brackets the collar between the wet bound, where uptake from the soil equals the cuticle’s residual flux (zero at g_min = 0), and the dry bound, where the stem reaches its critical potential or the root reaches its own; then it root-finds the first-order condition dprofit/dcollar = 0 inside the bracket, or reports that a bound binds.

  prepare_collar_solve()               kind := unsolved
        |
        +-- no interval to optimise over --------> hydraulic-shutdown
        +-- interval closed or inverted ---------> hydraulic-shutdown
        v
  root-find on dprofit/dcollar = 0      plain double, off any tape
        |
        +--> interior              profit stationary in the collar
        +--> boundary-soil         the wet bound binds
        +--> boundary-dry          the stem's dry end binds
        +--> solver-refused        the stationary point is a minimum
        +--> non-finite-gradient   no usable gradient at an endpoint
        v
  supply_draw_at -> operating_collar -> billed_outputs_at

What operating_collar does at each kind:

A dry pin is decided at the end itself. The solve probes the marginal a step inside the dry end and, where it rises there, solves the end to collar precision and reads the marginal there: falling, the maximum is interior however close to the end it sits, and its curvature is then taken from a difference in the logarithm of its distance to the end, where the marginal is close to linear; rising, the leaf is pinned. The two answers are one continuous function of the drivers, which a consumer differencing a whole run depends on: a pin decided a step inside the end sat past its maximum, and a run crossing that decision jumped there.

A narrow interval is polished, and a closed one is shut. When the feasible interval is narrower than bracket_tol_at(bound_b), both of its ends are solved to collar precision before anything is decided on it; one no wider than collar_root_tol is hydraulic-shutdown, as is a soil from which the stem can carry nothing (past the end of its curve). Every open leaf therefore takes its stem from continuity. The collar route loses precision as the root resistance goes to zero (#176).

The leaf route is the same surface with a different decision variable. There the stem’s upstream end is held at the soil, so the leaf potential is what the solve places, and decision_psi() returns it; supply_draw_at carries the soil potential as an input of its own and walks no supply. The solve records the same kinds, read off which candidate the one solver returned: an optimum at the wet end, where the stem carries the residual flux from the soil, is boundary-soil; one at the stem’s dry end is boundary-dry; a root of dprofit/dpsi_leaf = 0 is interior; the one feasible potential, and ProfitMax with no objective, are hydraulic-shutdown; an argmax resolved only to a width is solver-refused. operating_collar places each by its condition in the leaf potential, and the coordinates follow the leaf potential one for one. On the stem curve’s far tail the leaf route’s profit is flat to ~1e-17 per MPa, so its argmax is undetermined there; operating_collar refuses such a point as flatter than interior_curvature_floor.

A shut leaf is the ordinary evaluator at its wet bound (#182): E is the cuticle’s residual flux, the conductance to CO2 is the cuticle’s rho g_min, and ci solves the same supply-demand balance as an open leaf, above ca where the leaf is a net source. Nothing is written by hand, so the rows differentiate the value the solve reports and the shut leaf is continuous with the open one next to it. At g_min = 0 there is no conductance and the shut leaf takes the g_min -> 0 limit: the root of A(ci) = 0, or A as ci -> infinity where there is none. Light and heat shut no leaf: a leaf in the dark or too hot to fix carbon is an ordinary optimum, and the kind shade-death is gone. The shut stem sits at the collar rather than being held at a critical potential, so it keeps following a drying soil past psi_crit with no kink in the stem or root traits.

Why the implicit function theorem

Putting a trait on the tape and running the root-find live gives a derivative back through every iteration: the sensitivity of where that particular sequence of iterations stopped. It depends on the starting guess, the tolerance and the iteration count, converges to the right answer only as the tolerance tightens, and its error belongs to the solver’s settings rather than to the model.

The implicit function theorem gives the model’s derivative directly. With the point p* defined by a residual F(p*, theta) = 0,

    dp*/dtheta  =  -(dF/dp)^-1  dF/dtheta        at p*

Both partials are of the residual, which is a closed form, so neither needs anything the solver did. The kind decides which residual: stationarity of profit at an interior point, zero flux at the wet bound, stem continuity at the stem’s critical potential. One statement covers all of them. For profit at the optimum,

    dPi*/du  =  dPi/du |_p  -  sum_j mu_j dc_j/du |_p

At an interior optimum every multiplier mu_j is zero and this is the envelope theorem: differentiate the objective holding the optimum fixed. At a bound the multiplier term is exactly the correction the envelope theorem leaves out, and billed_outputs_at supplies it by letting profit see the collar move through the bound’s own implicit function. Only the objective is stationary in the collar, so only profit is taken with the collar held: carbon and the physical outputs carry the collar’s movement at every kind.

Why the kind is recorded, not inferred

dprofit_at_collar_psi returns a literal 0.0 on its shutdown and reversed-gradient exits. A test of the obvious form, marginal profit near zero means an interior optimum, reads a shut leaf as stationary, and a curvature formed by differencing the same function agrees with it, so the two checks confirm each other and are both wrong. The kind is therefore written by the branch the solve exited and never inferred afterwards. prepare_collar_solve resets it to unsolved at the start of every solve, so a branch that forgets reports unsolved rather than the previous solve’s kind.

A consumer branches on operating_point_kind() and never on a returned value. A consumer that records a solve and replays it later puts the point back with replay_operating_point(collar, kind), which hands over exactly the rows the solve did: a pinned solve places its collar at a polished bound and corrects its rows to that bound’s root, and the replay re-takes both.

Inside the surface

The draw is one walk. The soil layers are in parallel, coupled only through the collar potential, so the draw’s slope in the collar is a sum of per-layer terms and falls out of the same loop that accumulates the draw. The collar is stored passive in the draw, because the draw is the supply at a point; where the point moves is the supplied row’s business, and taking the collar live would record the search that placed it.

Second derivatives without a nested tape. kernel_slope_at returns a kernel’s slope as its value and the kernel’s mixed second partials as its rows. It evaluates the kernel once at a tangent-over-tangent scalar, both levels at double, so nothing reaches a tape; supply_derivatives then puts the result on the consumer’s tape as one statement. An adjoint nested inside a tangent would cost many times the kernel’s own statements, and odelia refuses that nesting at compile time.

Value from the table, rows from the curve. Both vulnerability curves read their value from a pre-integrated table, because the solve ran on the table, and their trait rows from the closed form, because a table carries none. closed_form_rows.hpp states which query slope each uses. The two differ by about 1e-4 relative at the default 400 knots, which bounds how far a row can sit from a finite difference of the table.

Extending it

for_each_active is a silent-zero hazard. odelia walks an aggregate’s active members through its for_each_active, and does nothing for an aggregate that has none. PhotoCapacity, SupplyDraw, HeldCoords and Normalisers declare it. A member added to either without being visited contributes no rows: its columns read exactly zero, every number stays finite, and nothing is raised. Visit every scalar member.

Another curve is its Curve<K>’s scalar forms (the SurfaceForms concept: its factor, its cost slope, the carbon part of its cost and, where they apply, its upstream slope and its price), beside the double forms the solve calls and held to them bit for bit by a test. Nothing in the surface is written per curve; a curve without its forms is refused at compile time. A curve that seeds normalisers before its search, as ProfitMax does, also needs them on the pack (normalisers_at).

A new fitted parameter goes before par_PPFD, because gradient::n_theta is par_PPFD; a slot appended last falls outside every caller’s theta.

What is tested

• test_supplied_rows: each layer’s supplied collar slope against a central difference of the uptake.
• test_transpose: a forward tangent against the reverse sweep through the same rows, on every curve at every kind it reaches, for every output of both calls. It cannot see a row that is wrong both ways.
• test_leaf: each curve’s scalar forms against its double forms, bit for bit; finite rows at every answered kind on a grid; a replay hands over its solve’s rows bit for bit; TF24_floor at a zero price is TF24 on the surface, bit for bit, in every value and row; every output of both calls equal to the solve’s own numbers at double; profit continuous, and every open leaf transpiring what its roots supply, across hydraulic shutdown and the light-compensation point; the collar, profit and output rows of every curve with the energy balance on, against a difference of the solve, at every kind each reaches; a second-difference sweep over soil potential, light, deficit and air temperature, with the energy balance off and on, one and three layers and the thermal-damage layer, finding no slope jump in profit, A, E, psi_leaf, lambda, or (with the energy balance off) two derivative rows at g_min > 0; the water budget E = E_s + r and E_up = E; the residual conductance’s rows against a difference of the solve; the collar, profit, carbon, A, E, psi_leaf, ci, gc, uptake and residual-flux rows of every curve at every kind it reaches against a difference of the solve, and a shut leaf whose stem is seated at its collar billing its root traits nothing; the root curve’s rows far along it against the closed form.

Not yet tested: boundary-dry driven through operating_collar on a tape beyond test_transpose’s one fixture.

Deliberately not built

A fully analytic block, supplying every row with nothing recorded anywhere, is possible in principle. It needs eight closed forms that do not exist: four mixed second partials of the assimilation kernel and four of the cost. The supply side stays recorded either way, so the statements walked per leaf would fall by about half, not to the handful the kernels alone cost, and the transpose test that referees the current rows could not referee those, because a wrong closed form makes both routes wrong in the same way.