Wombatlite april2026

[pearseb]

Hey @dougiesquire

These are small additions to WOMBAT-lite. Since I have been working on WOMBAT-mid a lot, some of the advances in WOMBAT-mid are transferable to WOMBAT-lite and probably should be the same. these include:

1. two ligand complexation with iron, not just one
2. dynamic prey switching by zooplankton
3. a control of colloidal shunt and two ligand complexation as logicals in the input.nml
4. a control of the dFe floor in the field_table

Sorry to add to your mountain of work but this should be quite simple and easy to merge with main.

[dougiesquire]

Thanks @pearseb. Some comments/questions below.

I was worried about the cost impact of doing 30 iterations for the two ligand bisection at every grid cell.

I ran 2 x 1 year runs of this configuration with do_two_ligands set to .true. and .false..

do_two_ligands = .false.

$ grep -r "payu_model_run_duration_seconds" payu_jobs/
payu_jobs/0/run/165895773.gadi-pbs.json:        "payu_model_run_duration_seconds": 1450.404177219607,
payu_jobs/1/run/165902070.gadi-pbs.json:        "payu_model_run_duration_seconds": 1430.7093261787668,

do_two_ligands = .true.

$ grep -r "payu_model_run_duration_seconds" payu_jobs/
payu_jobs/0/run/165895774.gadi-pbs.json:        "payu_model_run_duration_seconds": 1625.5848535923287,
payu_jobs/1/run/165903936.gadi-pbs.json:        "payu_model_run_duration_seconds": 1521.5471708243713,

So there is a pretty significant cost (~10%).

I think we could improve this using the Newton-Raphson method which I think suits the problem and will require much fewer iterations. I'll try to implement and test while you address my comments below.

[dougiesquire]

I rewrote and tested the two ligand bisection code to use Newton-Raphson. I tested across 10,000 sets of parameters across the ranges:

fe_sfe : 0.001  - 100  umol Fe/m3
K1     : 10  -  10000  m3/umol
K2     : K1 * 10^-1.5
L1     : 0.1   -  1.0  umol/m3
L2     : 0.5   -  3.0  umol/m3

where 

- K1/K2 = strong/weak ligand stability constants
- L1/L2 = strong/weak ligand concentrations

With 8 Newton-Raphson iterations, the solution for all sets of parameters produces a relative residual less than 1e-6, or an absolute residual that is less than that from 30 bisection iterations.

The timings (as above) for 8 Newton-Raphson iterations are:

do_two_ligands = .false.

$ grep -r "payu_model_run_duration_seconds" payu_jobs/
payu_jobs/0/run/165956251.gadi-pbs.json:        "payu_model_run_duration_seconds": 1425.3137839641422,

do_two_ligands = .true.

grep -r "payu_model_run_duration_seconds" payu_jobs/
payu_jobs/0/run/165956257.gadi-pbs.json:        "payu_model_run_duration_seconds": 1477.0212698340038,

So ~3% slower (as expected). As we discussed, let's make this change and default to do_two_ligands = .false..

I'll push these changes now.

[dougiesquire]

One small correction, and some formatting inconsistencies

[dougiesquire]

Thanks @pearseb!