Skip to content

p4zlys

Overview

The p4zlys module implements the calcium carbonate (CaCO₃) dissolution and precipitation processes in the PISCES model. This module calculates the degree of calcite saturation in seawater and computes the resulting dissolution or precipitation rates and its loss to the sediment pool.

Module Name: p4zlys
Model: PISCES (NEMO 5.0)
Subroutines: p2z_lys, p4z_lys, p4z_lys_init

The four main steps modelled are:

  1. Calculation of chemical parameters : pH and H⁺ and CO₃²⁻ concentrations
  2. Degree of saturation (Ω): Calculation to determine whether the water is undersaturated or supersaturated with CaCO\(_3\)
  3. Dissolution rate: Two regimes depending on the degree of saturation, with different reaction orders
  4. Update of tracers: DIC, alkalinity and particulate CaCO₃

Physical and Chemical Background

Calcite Saturation State

The saturation state (Ω) of calcite in seawater is defined as: \(\Omega = {{\ce{[Ca^{2+}] \cdot [CO3^{2-}]}} \over Ksp}\)

Where:

  • [Ca²⁺]: Calcium ion concentration (mol/kg)
  • [CO₃²⁻]: Carbonate ion concentration (mol/kg)
  • Ksp: Stoichiometric solubility product of CaCO₃ (mol\(^2\)/kg\(^2\))

Saturation Conditions

  • Ω > 1: Supersaturation → Precipitation favored
  • Ω = 1: Saturation equilibrium
  • Ω < 1: Undersaturation → Dissolution occurs

Model Implementation

1. Chemical Concentration Calculations

Hydrogen Ion Concentration [H⁺]

The hydrogen ion concentration is computed iteratively using the solve_at_general subroutine, which solves the carbonate system equilibrium equations:

zhinit(ji,jj,jk) = hi(ji,jj,jk) / (zdens + rtrn)
CALL solve_at_general(zhinit, zhi, Kbb)
hi(ji,jj,jk) = zhi(ji,jj,jk) * zdens

Carbonate Ion Concentration [CO₃²⁻]

Calculated from dissolved inorganic carbon (DIC) and equilibrium constants:

zco3(ji,jj,jk) = tr(ji,jj,jk,jpdic,Kbb) * ak13(ji,jj,jk) * ak23(ji,jj,jk) / 
                 (zhi(ji,jj,jk)**2 + ak13(ji,jj,jk) * zhi(ji,jj,jk) + 
                  ak13(ji,jj,jk) * ak23(ji,jj,jk) + rtrn)

Where:

  • ak13: First dissociation constant (H₂CO₃ ⇌ H⁺ + HCO₃⁻)
  • ak23: Second dissociation constant (HCO₃⁻ ⇌ H⁺ + CO₃²⁻)

pH Calculation

pH = -LOG10(MAX(hi(ji,jj,jk), rtrn))

2. Saturation State Computation

Calcium Concentration

Salinity-dependent calcium concentration:

zcalcon = calcon * (salinprac(ji,jj,jk) / 35.0)

Where:

  • calcon = 1.03×10⁻² mol/kg: Mean calcium concentration in seawater

Saturation State (Ω)

zdens = rhop(ji,jj,jk) / 1000.0
zomegaca = (zcalcon * zco3(ji,jj,jk)) / (aksp(ji,jj,jk) * zdens + rtrn)

Excess Saturation

excess(ji,jj,jk) = 1.0 - zomegaca
zexcess0 = MAX(0.0, excess(ji,jj,jk))

3. Dissolution Rate Calculation

The dissolution rate depends on the saturation state, with two distinct kinetic regimes:

Regime A: Strong Undersaturation (Ω < 0.8)

High-order reaction kinetics:

zdispot = kdca * zexcess0**nca

Parameters:

  • kdca: Dissolution rate constant (month⁻¹)
  • nca: Reaction order (typically ~4.5)

Characteristics:

  • Fast dissolution rate
  • Dominant in deep, undersaturated waters

Regime B: Near-Saturation (Ω ≥ 0.8)

Reduced-order kinetics:

zkd = kdca * 0.2**(nca - 0.2)
zdispot = zkd * zexcess0**0.2

Characteristics:

  • Slower dissolution rate
  • Transition regime near saturation equilibrium
  • Prevents unrealistic dissolution at near-equilibrium conditions

4. Tracer Updates

The dissolution process affects three biogeochemical tracers:

Dissolved Inorganic Carbon (DIC)

Increases by the amount of CaCO₃ dissolved

Total Alkalinity (ALK)

Increases by twice the dissolution rate (2 equivalents per CaCO₃)

Reaction: 

CaCO₃(s) → Ca²⁺ + CO₃²⁻

Particulate CaCO₃ (CAL)

Decreases by the amount dissolved

Implementation

PISCES-Simple

p2z_lys: Flux-Based Approach

Key Features:

  • Computes CaCO₃ production not explictly resolved by PISCES-Simple from vertical flux
  • Integrates vertical transport via sinking velocity (wsbio4)
  • Calculates export fluxes at specific depths (e.g., 100m)
  • Tracks CaCO₃ sedimentation to seafloor

Vertical Transport: 

zwsbio = wsbio4(ji,jj,jk) / rday
zdepexp = exp(-ztra(ji,jj,jk) * e3t(ji,jj,jk,Kmm) / zwsbio)
zcaco3(ji,jj,jk) = prodcal(ji,jj,jk) * rfact2r / ztra(ji,jj,jk) * 
                   (1.0 - zdepexp) + zcaco3(ji,jj,jk-1) * zdepexp

Bottom Sink: 

sinkcalb(ji,jj) = wsbio4(ji,jj,mbkt(ji,jj)) * zcaco3(ji,jj,mbkt(ji,jj)) * 
                  rfact2 / rday

PISCES-Operational/Research

p4z_lys: Concentration-Based Approach

Key Features:

  • Direct dissolution based on local CaCO₃ concentration
  • Simpler implementation without explicit vertical transport
  • Local equilibrium approach

Dissolution: 

zdissol = zexcess * tr(ji,jj,jk,jpcal,Kbb) * rfact2 / rmtss

Namelist Parameters

Configuration via the nampiscal namelist block:

Parameters

Parameter Description Typical Value Units
kdca Dissolution rate constant 1.0-10.0 month⁻¹
nca Reaction order for dissolution 4.0-5.0 dimensionless

Internal Constants

Constant Value Description
calcon 1.03×10⁻² Mean [Ca²⁺] in seawater (mol/kg)
rmtss nyear_len × rday / 12 Seconds per month

Diagnostic Outputs

The module provides several diagnostic fields for analysis:

Available Diagnostics

Field Description Units
PH pH value -
CO3 Carbonate ion concentration mmol/m³
CO3sat Saturation carbonate concentration mmol/m³
DCAL Calcite dissolution rate mmol/m³/s
PCAL Calcite production rate mmol/m³/s
EPCAL100 Calcite export at 100m (p2z only) mmol/m²/s
EXPCAL Calcite export profile (p2z only) mmol/m³/s

Output Calculations

pH: 

PH = -1.0 * LOG10(MAX(hi(:,:,:), rtrn))

CO₃ concentration (mmol/m³): 

CO3 = zco3(:,:,:) * 1.0e+3

CO₃ saturation: 

CO3sat = aksp(:,:,:) * rhop(:,:,:) / 1000.0 / calcon * 
         (salinprac(:,:,:) / 35.0) * 1.0e+3

Scientific Background and References

Key Physical Processes

  1. Temperature dependence: Via equilibrium constants (ak13, ak23, aksp)
  2. Pressure dependence: Via solubility product (aksp)
  3. Salinity dependence: Via calcium concentration adjustment
  4. Density effects: Proper unit conversion using seawater density (rhop)

Model Assumptions

  • Instantaneous chemical equilibrium for carbonate system
  • Magnesium carbonate effects implicitly included
  • No kinetic barriers to dissolution at strong undersaturation
  • Smooth transition between kinetic regimes

Implementation Notes

Numerical Considerations

  1. Division by zero protection: Use of rtrn (small threshold value)
  2. Exponential stability: Careful handling of exp() calculations in p2z_lys
  3. Time stepping: Explicit forward Euler scheme
  4. Unit conversions: Proper handling of mol/kg vs mol/m³

Usage

Initialization

CALL p4z_lys_init()  ! Read namelist and initialize parameters

Time Integration

! For p4z version (concentration-based)
CALL p4z_lys(kt, knt, Kbb, Krhs)

! For p2z version (flux-based with vertical transport)
CALL p2z_lys(kt, knt, Kbb, Kmm, Krhs)

Typical Workflow

  1. Read model configuration and namelist parameters
  2. Initialize carbonate chemistry module
  3. At each time step:
  4. Compute chemical equilibrium
  5. Calculate saturation state
  6. Determine dissolution rate
  7. Update biogeochemical tracers
  8. Output diagnostics

Dependencies

Modules

  • oce_trc: Shared ocean-tracer variables
  • trc: Passive tracer common variables
  • sms_pisces: PISCES source-minus-sink variables
  • p4zche: Chemical model (carbonate equilibrium)
  • p4zsink: Particle sinking (for p2z_lys)
  • prtctl: Print control for debugging
  • iom: I/O manager for diagnostics

External Functions

  • solve_at_general: Iterative solver for carbonate system
  • timing_start/stop: Performance profiling
  • iom_use/put: Diagnostic output

Document Version: 1.0
Last Updated: October 2025
Source code: NEMO 5.0