Legacy Fortran Solver Revived for Laminar Flow Benchmarks

A recovered 1995 Fortran finite-volume solver for laminar contraction flow has been recompiled, documented, and exercised as a reusable research-software object. The study focuses on a two-dimensional incompressible flow problem in a channel with a symmetric step contraction, a geometry that generates sharp gradients and recirculation zones useful for testing numerical methods. By reviving this specific code, the authors aim to demonstrate how older scientific software can be inspected, validated, and repurposed for modern reproducibility workflows.

The contraction flow problem is deceptively simple in its geometry but complex in its numerical behavior. A channel with a symmetric step contraction forces the computed flow to adjust through acceleration, pressure redistribution, and the possible formation of separated regions near the step. This setup tests whether a code can accurately represent the qualitative response of a viscous flow to a localized geometric singularity. The step corners create localized regions where gradients can become large, making the benchmark more revealing than a straight-channel test because it stresses the interaction among geometry, pressure correction, convection, and boundary conditions.

The original dissertation implemented a projection time-stepping method for the two-dimensional Stokes and Navier–Stokes equations. It compared three classical convection discretizations: central differencing, upwinding, and QUICK. The current reproduction revisits this framework to answer whether the recovered code can produce finite, nearly divergence-free fields on the original-style grid. The authors also investigate whether the three convective schemes reproduce the qualitative hierarchy reported in the 1995 dissertation and where modern reproducibility checks reveal behavior that might be missed from console output alone.

The software package is written in Fortran and uses a finite-volume spatial discretization with a staggered arrangement, where velocity components are stored on cell faces and pressure is stored at cell centers. The pressure is handled by a projection step involving a dense LU factorization through LAPACK routines. The convective term in the Navier–Stokes runs is approximated using one of the three schemes: central differencing, upwind differencing, or QUICK. The code is compact enough to be inspected directly, consisting of a global variable module, input parsing, array allocation, and assembly of the pressure Poisson matrix. This design favors conceptual directness over the sparse iterative solvers expected in modern large-scale CFD packages.

The reproduction campaign used a consistent 120 by 30 pressure grid and the same contraction ratio for all runs. The authors recorded the maximum velocity-change residual printed by the code and computed the maximum absolute discrete divergence from the output file. In valid reproduced runs, the maximum absolute divergence is typically of order 10 to the power of negative 12, indicating that the pressure-correction step performs its intended role even when the velocity residual has not reached the nominal tolerance. This divergence diagnostic serves as an independent check on the incompressibility constraint, which is central to the interpretation of the results.

At a Reynolds number of 20, all three convection schemes produced finite, nearly divergence-free fields. Their streamfunctions were close, with QUICK being especially close to central differencing, showing a mean absolute streamfunction difference of approximately 7.9 times 10 to the power of negative 5. The corresponding upwind-central difference was approximately 5.1 times 10 to the power of negative 4. This result is consistent with the expected behavior that QUICK is less diffusive than first-order upwinding. The residual histories for all three schemes decayed to the same practical residual band, although none reached the nominal tolerance of 10 to the power of negative 8.

At a Reynolds number of 500, QUICK remained stable and continued to reduce the residual through 10,000 time steps. The streamfunction changed very little between the 7,000-step and 10,000-step outputs, with a mean absolute change of approximately 1.4 times 10 to the power of negative 6. Although the printed residual did not reach the nominal tolerance, the field itself was effectively steady for qualitative plotting and comparison. This demonstrates the distinction between residual convergence and field stationarity, a nuance important for interpreting legacy code outputs.

At a Reynolds number of 2000, the schemes separated more clearly. Upwinding converged robustly to a residual plateau near 2 times 10 to the power of negative 8. QUICK remained finite and the projection step continued to produce small divergence, but the residual did not approach a steady-state tolerance in the tested runs. Central differencing failed due to overflow and NaNs under the tested settings. These observations reproduce the qualitative hierarchy reported in the original dissertation: upwinding is robust, QUICK is less diffusive but more demanding, and central differencing becomes fragile under higher Reynolds numbers.

One numerical-safety issue was identified during reproduction: if a run overflows and produces NaNs, the convergence function may report false convergence because comparisons with NaNs do not update the residual maximum. Such cases are identified as failures, even if the final console line says “converged.” Valid fields are required to contain no NaNs and to have a small divergence diagnostic. This finding highlights why reproducibility packages should test for finite-valued fields rather than relying only on a final convergence message, a lesson applicable to many legacy scientific codes.

The software is released under the MIT License and includes the preserved numerical source, an immutable recovered 1995 copy, archived example configurations, validation and post-processing scripts, automated tests, and documentation. The executable workflow is built with make and run with ./run.sh from the repository root. The input file step.dat controls all principal run parameters, including the method, Reynolds number, time step, convergence tolerance, and grid size. This simple text interface makes the code useful as a reproducibility object, allowing individual experiments to be reconstructed from a small input file and a console log.

The recovered solver is not intended to compete with modern production CFD packages but serves as a small, inspectable, and transparent benchmark. It is useful in teaching settings where students need to connect finite-volume discretization, pressure projection, convergence histories, and flow visualization without first learning a large software framework. The package can also serve as a lightweight reference case for researchers building automated workflows around legacy scientific codes, where build commands, input files, run logs, and derived figures must be archived together.

The results support the original numerical narrative with a more explicit reproducibility framing. The solver can be compiled and run with modern Fortran tooling, and its valid outputs are strongly divergence-free. Differences among convection schemes appear mainly in their stability, convergence rate, and numerical diffusion. The exercise shows why resurrecting legacy scientific codes requires more than successful compilation, necessitating explicit tests for finite residuals, termination on NaNs or infinities, and machine-readable run metadata. The contribution is a documented reconstruction of a classical finite-volume contraction-flow computation, highlighting the broader computational-science lesson that fast modern workflows still require conservative baseline tests and reconstructable run evidence.

Subject of Research: Step-contraction finite-volume solver: A reproducible fortran solver for laminar contraction flow

Article Title: Step-contraction finite-volume solver: A reproducible fortran solver for laminar contraction flow

Article References: Chaisi, M. (2026). Step-contraction finite-volume solver: A reproducible fortran solver for laminar contraction flow. SoftwareX, 36, Article 103199. https://doi.org/10.1016/j.softx.2026.103199

Image Credits: AI Generated

DOI: 10.1016/j.softx.2026.103199

Keywords: Step-contraction, finite-volume, solver, reproducible, fortran, laminar, contraction, flow, scientific research