Lagrangian particle histories provide a useful way to diagnose cloud-microphysical pathways in particle-based models, but their practical use in massively parallel super-droplet simulations depends on how particle identifiers, event records, output files, and post-processing workflows are implemented. This paper documents and evaluates a Lagrangian particle-tracking framework in SCALE-SDM 5.2.6-2.3.1. The framework includes forward tracking, backward tracking, and Two-Pass Hybrid Tracking (TPHT). Forward tracking saves the later states of super-droplets selected before or during a simulation. Backward tracking saves the rank and record location of each selected super-droplet's predecessor at the previous output time, so that post-processing can retrieve that record with an 𝒪(1) lookup instead of searching all super-droplets. Both methods use a common procedure to select all valid super-droplets, a random or height–radius-stratified sample, or an externally supplied target list. TPHT uses two integrations over the same physical time interval. The first records super-droplet identifiers when a chosen radius or coalescence condition is met. After duplicate identifiers are removed, the second starts from the same physical and random-number state and writes backward links only for this reduced set. Controlled computational tests are used to characterise runtime, peak random-access memory (RAM), storage (disk-output volume), the effect of output interval on file volume and mean write time, and sampling error and seed-to-seed variability. A three-dimensional DYCOMS-II RF02 stratocumulus case demonstrates TPHT. In this demonstration, the target criteria selected 19.3 % of the estimated valid super-droplet identities, corresponding to an estimated 5.2-fold reduction in retained backward output after temporary identifier files were removed. This factor is specific to the broad criteria used in this test; much larger reductions are possible for problems with sparser target populations. Agreement of the recorded coalescence events between the two integrations provided evidence of reproducibility under the tested configuration. The case study shows that TPHT can identify particles of interest and reconstruct their saved histories while reducing retained output. In future, the same approach could be applied to other history- or event-based questions, including entrainment and mixing, sedimentation, and ice-particle growth.
The selection of physical parameterization schemes is a major source of uncertainty in WRF simulations. A comprehensive evaluation of model performance requires the simultaneous consideration of multiple output variables, but simple multivariable averaging can dilute the influence of variables that are difficult to simulate or that respond strongly to alternative physical-process options. To address this issue, this study proposes a performance-based multivariable weighting framework for WRF parameterization-scheme evaluation. The framework combines two components: mean relative error, which represents variable-level simulation difficulty quantified by relative simulation error, and the coefficient of variation of absolute errors, which represents the scheme-dependent response of each variable within the evaluated physical-process group. The approach was applied and evaluated in assessing WRF parameterization schemes across two climatically distinct environments: the arid Northwest and the humid coastal Southeast of China. The results show that the method assigns greater weights to variables that are difficult to simulate, such as precipitation and 10 m wind, thereby helping identify physically interpretable and well-performing parameterization-scheme combinations under the prescribed evaluation framework. Compared with equal weighting and established statistical weighting methods, the proposed MRE × CV weighting helped identify parameterization-scheme combinations that achieved the lowest integrated cMIEI in both representative screening years and in the independent 2010 extreme-year validation. These results indicate that jointly considering variable-level simulation difficulty and scheme-dependent response improves the diagnostic relevance of multivariable WRF parameterization-scheme evaluation.
GPUs provide much higher memory bandwidth than CPUs. Data-parallel computations with regular memory access patterns exploit the GPU architecture efficiently. To solve the Stokes equations, we apply a preconditioned conjugate-gradient (PCG) method to the Schur complement system. The algorithm is dominated by sparse matrix-vector products, vector inner products, vector updates, and repeated velocity subproblem solves, all of which map onto GPU architectures. We implement the solver in GAUZZ v1.0.0 on GPUs by coupling CUDA-enabled PETSc with FEniCS. PETSc AIJCUSPARSE matrices with CUDA vectors keep operators and field vectors resident in VRAM, which limits CPU-GPU transfers. The P1 mass matrix, the gradient operator, and the divergence operator are reused in the pressure-correction L2-projection step. HYPRE BoomerAMG was adopted for preconditioning the velocity-only subproblem. Accuracy and performance are evaluated on a manufactured solution, the SolCx benchmark, the sticky-air benchmark, 2D Rayleigh–Taylor instabilities, 3D thermo-mechanical convection, and a 2D subduction model with nonlinear viscosity. We confirm that the GPU-based implementation reproduces CPU solutions. Performance is quantified by the execution time of the pressure-correction step and by throughput, defined as the number of Stokes degrees of freedom divided by the execution time of one Stokes solve. A single GPU delivers 6–11× the throughput of a 32-core CPU across the Stokes benchmarks. Two GPUs increase the throughput at sufficiently large problem sizes, reaching up to 2.1× the single-GPU throughput. The two-GPU configuration also extends the maximum problem size to 2.5×107 DOFs, beyond the 48 GB VRAM limit of a single GPU.
Spectropolarimetry is a powerful tool for characterising planetary atmospheres and surfaces. For the design and operation of spectro(polari)metric instrumentation, numerically simulated signals of the measured radiation are essential. We present MONKI (Monte Carlo KNMI), an efficient and accurate Monte Carlo radiative transfer code written in Fortran. MONKI computes total, linearly polarised, and circularly polarised radiances reflected by and transmitted through planetary atmospheres, fully accounting for the polarisation of light in all orders of scattering. Unlike classical Monte Carlo approaches that rely on post-sampling weight corrections and can therefore produce rare but disproportionately large photon weights, MONKI samples scattering directions from the full polarisation-dependent distribution, consistently in both forward and backward Monte Carlo modes. We provide the corresponding theoretical framework in the appendices. Because of this exact polarised sampling approach, the code converges robustly for optically thick and/or strongly polarising atmospheres. MONKI can handle horizontally homogeneous atmospheres as well as atmospheres with horizontal inhomogeneities, such as three-dimensional (3D) patchy clouds. We validate MONKI through comparisons with other radiative transfer codes. Finally, we present sample simulations of sunlight reflected by the Earth and Venus, and explain the total and polarised radiance features by analysing the altitudes at which photons are scattered. We conclude that MONKI is a versatile and accurate tool for simulations and detailed analyses of total and polarised light emerging from the atmospheres of the Earth, Venus, and, in principle, any other planet.
River plumes can strongly influence coastal dynamics but are poorly represented in global ocean models due to stringent resolution requirements. Following the idealized simulations presented in a companion paper (Hinson et al., 2026), this study presents submesoscale-permitting realistic simulations using unprecedented regional refinement with MPAS–O. We assess MPAS–O's ability to represent river plume dynamics by comparing with a validated ROMS configuration focused on the Mississippi–Atchafalaya (M–A) plume in the northern Gulf of Mexico (GoM). The variable-resolution mesh spans 1.4 km in the GoM to 100 km in the Pacific and Indian Oceans. We incorporate several improvements to the representation of river forcing in MPAS–O: the pseudo point source treatment of runoff, river temperature prescribed from air temperature, and passive tracers to track freshwater transport. MPAS–O generates vigorous submesoscales in the open GoM, as quantified by probability density functions of surface relative vorticity and divergence. While the model qualitatively captures the M–A plume's seasonal evolution, it does not realistically reproduce the summer submesoscale eddy field on the Texas–Louisiana shelf. A persistent brackish lens associated with reduced vertical mixing forms and is present throughout the M–A plume. We hypothesize that this arises from interacting parameterization and structural design choices within MPAS–O that affect plume stratification and instability growth under realistic forcing, rather than from deficiencies in numerics. Our results demonstrate MPAS–O's potential for simulating submesoscale processes in the open ocean while highlighting challenges in using regional refinement to model coastal dynamics.
