MFlowCode/MFC

Euler–Lagrange two-way coupling gives rank-count-dependent results — investigated: round-off seed amplified by chaos, not a decomposition bug Summary Moving Lagrangian bubbles with two-way coupling (lag_params%vel_model = 2, lag_params%solver_approach = 2) produce visibly different results when the same case is run with different MPI rank counts. A systematic bisection campaign (attached tarball) traced this to a one-ulp (2.2e-16) summation-order difference in the two-way smearing deposition wherever a bubble's kernel straddles a subdomain boundary — the halo-accumulate reduce in s_mpi_reduce_beta_variables_buffers sums the same contributions in a decomposition-dependent order. The strongly chaotic coupled-bubble dynamics then amplify that seed exponentially (measured e-folding time ≈ 0.38 nondim time units) to O(1) trajectory differences within a few hundred steps. Verdict: there is no decomposition bug in the moving-bubble / two-way-coupling path. Every mechanism tested in is

Summary s_compute_viscous_stress_tensor (added in #1614 for higher-order central-difference IB drag) performs an out-of-bounds heap read whenever an immersed-boundary body's finite-difference stencil is centered on a cell adjacent to a domain boundary. The bug has been latent in master since #1614 (2026-06-24); it is silent in optimized builds and only manifests as a hard crash under bounds-checking (./mfc.sh ... --debug/reldebug) builds. cc @danieljvickers (author of #1614 and #1618) Root cause In src/simulation/m_viscous.fpp, the velocity-gradient stencil indexes the FD-coefficient arrays by cell: velocity_gradient_tensor(l, 2) = velocity_gradient_tensor(l, 2) & + fd_coeff_y(r, j)*q_prim_vf(...)%sf(i, j + r, k) But fd_coeff_x/y/z are allocated over the interior only — e.g. fd_coeff_y(-fd_number:fd_number, 0:n) (in m_derived_variables.fpp). The IB drag path (s_compute_ib_forces in m_ibm.fpp) evaluates this stencil at body ghost-point cells, which for a body sitting against a bou

Summary Add a helper (toolchain utility and/or smarter case-file parameterization) to make hyperbolic-tangent grid-stretching parameters easier to choose. Today users have to hand-tune the raw stretch parameters (a_x, x_a, x_b, loops_x, etc.) by trial and error, which is unintuitive and error-prone. Motivation This came out of the discussion on #1625, which fixed MPI-rank-dependent divergent grid stretching. Now that stretching is rank-invariant, it's a good time to make setup friendlier. As @danieljvickers suggested: I think it may also be neat to add a grid-stretching helper to the toolchain or a smarter parameter grid stretching in the case files to make getting your parameters more obvious. Something to help people with parameter selection based on their desired domain size, uniform region, uniform spacing, etc. Proposed idea Let users specify what they actually care about, and derive the raw stretch parameters from it: desired total domain size extent of the uniform (unstretche

In the recent addition of higher-order finite differencing for immersed boundaries from #1614, the fd_order case file parameter has begun to be overloaded in terms of its use. 4 different case file settings all use the same finite difference order. This means that we are somewhat less flexible and the term appears more ambiguous to users. It would be helpful to consider a refactor of the finite difference scheme for paramters like probes, hypo-elasticity, immersed boundaries, MHD, post-processing, etc. In particular, the finite difference order used in all of these terms should be disambiguated, The goal seems to be separate the use cases of this in pre_process, simulation, post_process to be distinct

The IB patches have just been refactored to significantly reduce the number of lines of code required in the module. This was done by making use of shared subroutines and functions that are common to all IB patch set up, with mild if/then trees for various patches. There is arguably more overlapping/redundant code in the ICPP patches than there was in the IB patch code, and similar improvements can be made. The product of this should be unifying as many subroutines as possible to use common initial set up and geometry checking. The result should be a significant reduction in the number of lines of code in the module.