AqueousDiffusionReaction

{
  "reference": {
    "citation": "Seinfeld, J. H., & Pandis, S. N. (2006) Atmospheric Chemistry and Physics, 2nd ed., Wiley. Ch. 12, Eqs. 12.105 and 12.110. Migrated from Aerosol.jl:src/aqueous_transport.jl:14 (esm-x0e).",
    "notes": "source_ref: Aerosol.jl@c8c640bb636b22f2425460e9ed5cba50793bfee1:src/aqueous_transport.jl:14\n\nUpstream `AqueousDiffusionReaction` is an MTK `@component` with three `@parameters` (R_p, k_rxn, D_aq) and two `@variables` (q, Q), both algebraic (no `D()` time-derivative). The `eqs` array uses `~` (algebraic constraint) for both:\n  q ~ R_p * sqrt(k_rxn / D_aq)         [Eq. 12.105 — Hatta number]\n  Q ~ 3 * ((exp(q)+exp(-q))/(exp(q)-exp(-q))/q - 1/q^2)  [Eq. 12.110 — correction factor]\n\nPhysical interpretation: the Hatta number q = R_p√(k_rxn/D_aq) compares the reacto-diffusive length scale l_d = √(D_aq/k_rxn) to the droplet radius R_p. When q << 1 (slow reaction or small droplet), species diffuses freely throughout the interior before reacting and Q → 1 (no correction). When q >> 1 (fast reaction or large droplet), the reaction is confined near the droplet surface because the interior is depleted, and Q → 3/q → 0. The factor Q multiplies the bulk aqueous-phase reaction rate to give the effective rate; it appears as the `Q` parameter in the sibling `DropletMassBalance` component (src/aqueous_transport.jl:112).\n\nBoth variables become observed quantities after MTK `structural_simplify`: `q` is an explicit function of parameters only; `Q` is an explicit function of `q` (and hence also of parameters only). The ESM serializes both as `type: observed` following the established convention for purely-algebraic MTK components (cf. sibling `MassTransportLimitation.esm` and `AqueousEquilibrationTimescale.esm`). The `equations` array is empty.\n\nDependency ordering: `Q` references `q` in its expression. The ESM runtime topologically sorts observed variables, so `q` is evaluated first, then `Q`.\n\ncoth encoding: there is no native `coth` op in the ESS operator set, so coth q is expressed as the ratio of native hyperbolic ops cosh(q)/sinh(q). This equals the algebraic identity (exp(q)+exp(-q))/(exp(q)-exp(-q)) — the factor of ½ in each of cosh/sinh cancels in the ratio. The upstream Aerosol.jl code and the original ESM migration both wrote the exp-expansion directly; bead esm-ji2 substituted the native cosh/sinh ops (added to the ESS static operator set at origin/main ≥ 8b8d3056, bead ess-v9a.1). Besides removing the hand-rolled expansion, the native ops are better-conditioned at small q: the exp difference exp(q)-exp(-q) suffers catastrophic cancellation as q→0, whereas sinh(q) is accurate to machine precision (see the kinetically_controlled_slow_reaction test, where accuracy improves ~180×). For large q (fast reaction), the formula evaluates to 3*(1/q − 1/q²) ≈ 3/q → 0; cosh and sinh remain finite well past the atmospheric q range.\n\nTolerance rationale (model-level rel=1e-9): the model is a composition of sqrt and the native hyperbolic ops (cosh, sinh) over parameters, plus simple arithmetic. Float64 round-off is ~1e-15 relative; rel=1e-9 leaves six orders of headroom, matching sibling observed-only ESMs (mass_transport_limitation, aqueous_diffusion_timescale)."
  },
  "tolerance": {
    "rel": 1e-09
  },
  "variables": {
    "R_p": {
      "type": "parameter",
      "units": "m",
      "default": 1e-05,
      "description": "Droplet radius (Eq. 12.105 term R_p). The Hatta number scales as R_p, so larger droplets are more susceptible to intra-droplet diffusion limitation. Default 1e-5 m is a 10 μm cloud droplet."
    },
    "k_rxn": {
      "type": "parameter",
      "units": "1/s",
      "default": 1.0,
      "description": "First-order aqueous-phase reaction rate constant (Eq. 12.105 term k_rxn). Faster reactions increase q and reduce Q, driving the system toward surface-limited uptake."
    },
    "D_aq": {
      "type": "parameter",
      "units": "m^2/s",
      "default": 1e-09,
      "description": "Aqueous-phase molecular diffusivity of the dissolved species (Eq. 12.105 term D_aq). Higher diffusivity reduces q and increases Q toward 1. Default 1e-9 m²/s is a typical small-solute aqueous diffusivity."
    },
    "q": {
      "type": "unknown",
      "units": "1",
      "description": "Hatta number (dimensionless reaction-diffusion parameter, Eq. 12.105): q = R_p √(k_rxn/D_aq). The ratio of the droplet radius to the reacto-diffusive length scale l_d = √(D_aq/k_rxn). q < 1: kinetically controlled (well-mixed interior); q > 1: diffusion-limited (surface-confined reaction)."
    },
    "Q": {
      "type": "unknown",
      "units": "1",
      "description": "Aqueous diffusion correction factor (Eq. 12.110): Q = 3(q·coth q − 1)/q² = 3(cosh(q)/(sinh(q)·q) − 1/q²), with coth q evaluated via the native cosh/sinh ops. Ranges from 1 (q→0, kinetically controlled, no correction) to 0 (q→∞, strongly diffusion-limited). Applied as a multiplier on the aqueous-phase reaction rate in the sibling DropletMassBalance component."
    }
  },
  "equations": [
    {
      "lhs": "q",
      "rhs": {
        "op": "*",
        "args": [
          "R_p",
          {
            "op": "sqrt",
            "args": [
              {
                "op": "/",
                "args": [
                  "k_rxn",
                  "D_aq"
                ]
              }
            ]
          }
        ]
      }
    },
    {
      "lhs": "Q",
      "rhs": {
        "op": "*",
        "args": [
          3,
          {
            "op": "-",
            "args": [
              {
                "op": "/",
                "args": [
                  {
                    "op": "cosh",
                    "args": [
                      "q"
                    ]
                  },
                  {
                    "op": "*",
                    "args": [
                      {
                        "op": "sinh",
                        "args": [
                          "q"
                        ]
                      },
                      "q"
                    ]
                  }
                ]
              },
              {
                "op": "/",
                "args": [
                  1,
                  {
                    "op": "^",
                    "args": [
                      "q",
                      2
                    ]
                  }
                ]
              }
            ]
          }
        ]
      }
    }
  ],
  "tests": [
    {
      "id": "default_parameters",
      "description": "Default parameters (R_p=1e-5 m, k_rxn=1 s⁻¹, D_aq=1e-9 m²/s). Hatta number q = 1e-5 × √(1/1e-9) = 0.31623 (moderate diffusion-reaction coupling). Correction factor Q = 3(coth(0.31623)/0.31623 − 1/0.31623²) ≈ 0.9934 — a small (0.7%) reduction from unity, meaning diffusion barely limits the reaction at these default conditions.",
      "time_span": {
        "start": 0.0,
        "end": 1.0
      },
      "assertions": [
        {
          "variable": "q",
          "time": 0.0,
          "expected": 0.31622776601683794
        },
        {
          "variable": "Q",
          "time": 0.0,
          "expected": 0.993396196825245
        }
      ]
    },
    {
      "id": "diffusion_limited_fast_reaction",
      "description": "Fast reaction k_rxn=1e6 s⁻¹ (k_rxn increased 1e6-fold). Hatta number q = 1e-5 × √(1e6/1e-9) = 316.23 >> 1: the reacto-diffusive length l_d = √(D_aq/k_rxn) = √(1e-9/1e6) = 1e-7.5 m is ~0.03 μm, far smaller than the 10 μm droplet. Correction factor Q ≈ 3/q ≈ 0.00946 — the reaction is strongly surface-confined, reducing the effective rate by ~99%.",
      "time_span": {
        "start": 0.0,
        "end": 1.0
      },
      "parameter_overrides": {
        "k_rxn": 1000000.0
      },
      "assertions": [
        {
          "variable": "q",
          "time": 0.0,
          "expected": 316.22776601683796
        },
        {
          "variable": "Q",
          "time": 0.0,
          "expected": 0.009456832980505137
        }
      ]
    },
    {
      "id": "kinetically_controlled_slow_reaction",
      "description": "Slow reaction k_rxn=1e-4 s⁻¹ (k_rxn reduced 1e4-fold). Hatta number q = 0.003162 << 1: the reacto-diffusive length l_d >> R_p so the species fully equilibrates across the droplet interior before reacting. Q ≈ 1 − q²/15 ≈ 0.999999 — diffusion correction is negligible, well-mixed limit applies. ACCURACY NOTE (esm-ji2): this small-q point is the one place where swapping the hand-rolled (exp(q)±exp(−q)) coth standin for native cosh(q)/sinh(q) changes the value. The standin lost ~2.8e-9 relative precision to catastrophic cancellation in exp(q)−exp(−q), so the prior stored expected 0.9999993361125229 inherited that error. Native sinh/cosh are well-conditioned at small argument: the updated expected 0.9999993333185557 lies within 1.5e-11 of the analytic value 0.9999993333339683 (mpmath, 50 dps) — about 180× more accurate, and comfortably inside the model rel=1e-9 tolerance.",
      "time_span": {
        "start": 0.0,
        "end": 1.0
      },
      "parameter_overrides": {
        "k_rxn": 0.0001
      },
      "assertions": [
        {
          "variable": "q",
          "time": 0.0,
          "expected": 0.003162277660168379
        },
        {
          "variable": "Q",
          "time": 0.0,
          "expected": 0.9999993333185557
        }
      ]
    },
    {
      "id": "low_diffusivity_intermediate_hatta",
      "description": "Low aqueous-phase diffusivity D_aq=1e-10 m²/s (10× lower than default; representative of a larger or more complex solute molecule). Hatta number q = 1e-5 × √(1/1e-10) = 1.0 exactly (the droplet radius equals the reacto-diffusive length). Correction factor Q = 3(coth(1) − 1) = 3((e²+1)/(e²−1) − 1) ≈ 0.9391 — a notable 6% reduction.",
      "time_span": {
        "start": 0.0,
        "end": 1.0
      },
      "parameter_overrides": {
        "D_aq": 1e-10
      },
      "assertions": [
        {
          "variable": "q",
          "time": 0.0,
          "expected": 1.0
        },
        {
          "variable": "Q",
          "time": 0.0,
          "expected": 0.9391058564979942
        }
      ]
    }
  ],
  "analyses": [
    {
      "id": "correction_factor_vs_reaction_rate",
      "description": "Correction factor Q and Hatta number q as a function of aqueous-phase reaction rate constant k_rxn (log sweep 1e-4 to 1e6 s⁻¹) at default R_p=1e-5 m, D_aq=1e-9 m²/s. Q decreases from ~1 (slow reaction, kinetically controlled) toward 0 (fast reaction, diffusion-limited). The transition occurs near q=1, i.e. k_rxn ≈ D_aq/R_p² = 1e-9/(1e-5)² = 10 s⁻¹.",
      "time_span": {
        "start": 0.0,
        "end": 1.0
      },
      "parameter_sweep": {
        "type": "cartesian",
        "dimensions": [
          {
            "parameter": "k_rxn",
            "range": {
              "start": 0.0001,
              "stop": 1000000.0,
              "count": 40,
              "scale": "log"
            }
          }
        ]
      },
      "plots": [
        {
          "id": "Q_vs_k_rxn",
          "type": "line",
          "description": "Correction factor Q vs reaction rate k_rxn. Shows the transition from Q=1 (kinetically controlled) to Q→0 (diffusion-limited) near k_rxn ~ D_aq/R_p² = 10 s⁻¹.",
          "x": {
            "variable": "k_rxn",
            "label": "k_rxn (s⁻¹)"
          },
          "y": {
            "variable": "Q",
            "label": "Q (dimensionless)"
          }
        },
        {
          "id": "q_vs_k_rxn",
          "type": "line",
          "description": "Hatta number q vs reaction rate k_rxn (log-log slope 0.5, since q ∝ √k_rxn). The q=1 threshold marks the onset of significant diffusion limitation.",
          "x": {
            "variable": "k_rxn",
            "label": "k_rxn (s⁻¹)"
          },
          "y": {
            "variable": "q",
            "label": "q (dimensionless)"
          }
        }
      ]
    },
    {
      "id": "correction_factor_vs_droplet_radius",
      "description": "Correction factor Q as a function of droplet radius R_p (log sweep 0.1 μm to 100 μm) at default k_rxn=1 s⁻¹, D_aq=1e-9 m²/s. Larger droplets have larger q = R_p/l_d and lower Q — they are more diffusion-limited because the interior is harder to reach. Small droplets (R_p < l_d = √(D_aq/k_rxn) = 31.6 μm) are kinetically controlled (Q ≈ 1).",
      "time_span": {
        "start": 0.0,
        "end": 1.0
      },
      "parameter_sweep": {
        "type": "cartesian",
        "dimensions": [
          {
            "parameter": "R_p",
            "range": {
              "start": 1e-07,
              "stop": 0.0001,
              "count": 40,
              "scale": "log"
            }
          }
        ]
      },
      "plots": [
        {
          "id": "Q_vs_Rp",
          "type": "line",
          "description": "Correction factor Q vs droplet radius R_p. Q is near 1 for small droplets (R_p << reacto-diffusive length) and decreases for larger droplets.",
          "x": {
            "variable": "R_p",
            "label": "R_p (m)"
          },
          "y": {
            "variable": "Q",
            "label": "Q (dimensionless)"
          }
        }
      ]
    }
  ]
}