Skip to content

Implements partitioned Fluid-Structure Interaction (FSI) coupling in svMultiPhysics#554

Draft
Eleven7825 wants to merge 115 commits into
SimVascular:mainfrom
Eleven7825:feature/issue-431-partitioned-fsi
Draft

Implements partitioned Fluid-Structure Interaction (FSI) coupling in svMultiPhysics#554
Eleven7825 wants to merge 115 commits into
SimVascular:mainfrom
Eleven7825:feature/issue-431-partitioned-fsi

Conversation

@Eleven7825
Copy link
Copy Markdown
Collaborator

This PR implements partitioned Fluid-Structure Interaction (FSI) coupling in svMultiPhysics, building on the work in mrp089:feature/issue-431-partitioned-fsi. It introduces a PartitionedFSI class that orchestrates three independently solved sub-simulations (fluid, solid, mesh) with a Dirichlet-Neumann coupling loop and Aitken relaxation.

This is the first step toward the broader partitioned FSI + Growth and Remodeling framework described in issue #528, which proposes InterfaceData, an abstract CouplingAlgorithm strategy, and a PartitionedSimulation driver capable of supporting multi-rate time integration and IQN/JFNK acceleration. The architecture introduced here — three independent Simulation objects with explicit interface transfer — is designed to extend cleanly toward those goals.

What's included

  • PartitionedFSI class: drives the coupling loop between independently initialized fluid, solid, and mesh sub-simulations
  • Dirichlet-Neumann coupling: solid displacement → mesh solve → fluid ALE solve → fluid traction → solid solve
  • Interface relaxation: constant-factor and Aitken Δ² acceleration (Küttler & Wall 2008)
  • MPI multi-core support: all three sub-simulations are fully distributed; interface arrays are globally replicated via MPI_Allgatherv so convergence checks and Aitken updates are consistent across all ranks
  • initialize_partitioned_fsi in Simulation: broadcasts config from rank 0 to all ranks so slave ranks enter the coupling constructor (fixes deadlock on >1 MPI process)
  • Test case: tests/cases/fsi/pipe_3d_partitioned/ — 3D pipe flow with a compliant wall, ramp and steady-state XMLs
  • Coupling convergence log written to coupling.dat each time step

Eleven7825 and others added 30 commits October 13, 2025 01:57
- Add CMake-generated files (CMakeCache.txt, CMakeFiles/, Makefile, etc.)
- Add build directories (*-build/, *-prefix/)
- Add compiled outputs (Code/bin/, Code/lib/)
- Add generated headers
- Add vim swap files (*.swp, *.swo, *~)

Related to SimVascular#442
- Create integrator.h and integrator.cpp with Integrator class
- Encapsulates solution variables (Ag, Yg, Dg)
- Handles Newton iteration loop
- Manages linear system assembly and solve
- Applies boundary conditions
- Update CMakeLists.txt to include new source files

This is the first step toward implementing partitioned FSI as described in SimVascular#442
- Replace Newton iteration loop with Integrator::step() call
- Remove local Ag, Yg, Dg arrays (now managed by Integrator)
- Simplify iterate_solution() function
- Add integrator.h include

Addresses SimVascular#442: Encapsulate Newton iteration in main.cpp
The iEqOld variable is needed for output::output_result() calls
after the Newton iteration completes. This variable tracks which
equation was solved in the last iteration.
- Rename integrator.h to Integrator.h
- Rename integrator.cpp to Integrator.cpp
- Update all #include statements in main.cpp and Integrator.cpp
- Update CMakeLists.txt reference

Addresses PR SimVascular#450 review feedback on file naming convention.
- Use JavaDoc-style /** */ comments throughout Integrator.h
- Add @brief, @param, and @return tags following svZeroDSolver conventions
- Document all public methods and private members
- Improve descriptions for solution variables (Ag, Yg, Dg) and helper methods

Addresses PR SimVascular#450 review feedback on Doxygen documentation format.
- Add istr_ as private member variable in Integrator class
- Update step() method to set istr_ instead of using local variable
- Replace all istr references with istr_ throughout the code
- Simplify debug write statements in apply_boundary_conditions()

This change enables individual helper functions to access istr_ for
their own debug outputs, preparing for the next step of moving write
statements into individual functions.

Addresses PR SimVascular#450 review feedback on variable management.
- Move all .write() debug statements from step() to their respective helper functions
- initiator_step() now handles its own debug output (Ag, Yg, Dg, Yn)
- allocate_linear_system() now handles Val output
- set_body_forces() now handles Val output
- assemble_equations() now handles R and Val output
- apply_boundary_conditions() now handles Val, R, Yg, Dg output
- solve_linear_system() now handles Val and R output
- corrector_and_check_convergence() now handles Yn output

This makes step() cleaner and more focused on orchestration, while
each helper function manages its own debug logging. All functions
now use the istr_ class member variable for consistent debug naming.

Addresses PR SimVascular#450 review feedback on code organization.
…oop'

- Rename inner_count_ to newton_count_ throughout Integrator class
- Update comments to use 'Newton iteration' instead of 'inner loop'
- Update debug messages to show 'Newton Iteration' instead of 'Inner Loop'

This makes the terminology more precise and consistent with standard
computational mechanics nomenclature, clarifying that we're referring
to the Newton iteration within each time step (as opposed to the outer
time loop).

Addresses PR SimVascular#450 review feedback on terminology consistency.
- Move debug print statements from step() into individual helper functions
- Each function now handles its own debug output:
  * initiator_step()
  * allocate_linear_system()
  * set_body_forces()
  * assemble_equations()
  * apply_boundary_conditions()
  * update_residual_arrays()
  * solve_linear_system()
  * corrector_and_check_convergence()
- Makes step() method cleaner and more focused on orchestration
- Each helper function is now self-contained with its own debug logging

This addresses PR SimVascular#450 review feedback on code organization.
Moved all pic namespace functions (picp, pici, picc, pic_eth) into the
Integrator class as member functions:
- picp -> predictor() (public method, called from main.cpp)
- pici -> pici() (private method)
- picc -> picc() (private method)
- pic_eth -> pic_eth() (private method)

Updated all call sites to use the new member functions.
Removed pic.h and pic.cpp files and updated CMakeLists.txt.

Related to GitHub issue SimVascular#459.
Move time integration variables An (acceleration), Dn (displacement), and
Yn (velocity) from the global ComMod structure into the Integrator class.
These variables are now passed as function parameters throughout the
codebase instead of accessing them through ComMod.

Changes include:
- Add An_, Dn_, Yn_ members to Integrator with public getters
- Initialize from Ao, Do, Yo in Integrator::initialize_arrays()
- Update 19 files to pass these variables as Array<double>& parameters
- Remove An, Dn, Yn declarations from ComMod.h
This commit completes the extraction of time-level variables from the global
ComMod structure to the Integrator class. Ao, Do, Yo (old acceleration,
displacement, velocity at time level n) join An, Dn, Yn as Integrator members.

Key changes:
- ComMod.h: Removed Ao, Do, Yo member variables
- Integrator.h/cpp: Added Ao_, Do_, Yo_ members with move constructor and getters
- Simulation.h/cpp: Added initialize_integrator() to transfer ownership
- initialize.cpp: Made Ao, Do, Yo local variables, moved to Integrator
- Updated 25+ files with function signatures to pass Ao/Do/Yo parameters
- Fixed default parameter placement (declarations only, not definitions)

All time integration state is now encapsulated in the Integrator class,
improving modularity and preparing for future multi-physics enhancements.
…nd integration functions

Following the extraction of Ao/Do/Yo from ComMod to Integrator class, this commit
completes the parameter threading work by updating all functions that perform face
integration or boundary condition operations to accept and pass Do/Dn parameters.

Changes include:

Boundary Condition Functions:
- bc_ini(): Add Do parameter, update 3 integ() calls for parabolic profiles and flux normalization
- set_bc_cmm(), set_bc_cmm_l(): Add Do parameter for CMM boundary conditions
- set_bc_neu(), set_bc_neu_l(): Add Do parameter for Neumann boundary conditions
- set_bc_dir_w(), set_bc_dir_wl(): Add Do parameter for Dirichlet wall conditions
- set_bc_trac_l(): Add Do parameter for traction boundary conditions
- set_bc_rbnl(): Add Do parameter for Robin boundary conditions
- set_bc_cpl(): Update to pass Do to bc_ini() for coupled BC area recalculation

FSI and Mesh Functions:
- fsi_ls_upd(): Add both Dn and Do parameters for level set updates
- fsi_ls_ini(): Add Do parameter for FSI initialization
- b_assem_neu_bc(), b_neu_folw_p(): Add Do parameter for Neumann BC assembly

RIS (Reduced Immersed Surface) Functions:
- ris_meanq(): Add Do parameter, update 2 integ() calls for pressure/flow computation
- ris0d_status(): Add Do parameter, update 3 integ() calls for 0D coupling
- ris_resbc(), setbc_ris(), ris0d_bc(): Add Do parameter for RIS boundary conditions

CMM (Coupled Momentum Method):
- cmm_b(): Add Do parameter, update nn::gnnb() call with full parameters

Integration Functions:
- All face-based integ() calls updated to pass Do parameter via pointer
- Updated calls to use proper signatures with pFlag, cfg, Dn, and Do parameters
- Enhanced error reporting in nn::gnnb() to include face name, mesh name, and element

Main Integration Loop:
- Integrator::iterate(): Updated all set_bc calls to pass Do_ parameter
- main.cpp: Updated ris_meanq() and ris0d_status() calls to pass Do
This commit fixes systematic issues with displacement parameter passing throughout the solver:

1. Corrected parameter confusion in set_bc.cpp where An/Ao were incorrectly used instead of Dn/Do
2. Replaced nullptr with proper Dn/Do parameters in 14 integ() calls across ris.cpp, txt.cpp, baf_ini.cpp, and set_bc.cpp
3. Fixed critical bug in txt.cpp:520 where missing pFlag parameter caused &Do to be misinterpreted as boolean
4. Enhanced error messages in all_fun.cpp to identify exact failure locations

These fixes resolve 12 test failures in moving mesh simulations (genBC and RIS tests).
…and `solution.current.A`. However, getter method for each variables (`get_Ao` for example) still exists for backward compatibility purpose.
mrp089 and others added 25 commits April 4, 2026 23:13
…ection

- Add Coupling_method XML parameter: "constant", "aitken" (default), "iqn-ils"
- Refactor relax_interface into separate methods per algorithm
- Implement IQN-ILS (Degroote et al. 2009) with QR decomposition,
  NaN safeguards, and Aitken fallback on instability
- Add NaN detection in Integrator::step_equation Newton loop to abort
  immediately on divergence instead of running max iterations
- Add NaN/divergence check after interface relaxation

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
- Use correct update formula: x_new = x_tilde + W*c (not (W-V)*c)
- Persist V/W columns across time steps (trimmed to 20 max)
- QR filtering with eps=0.1 threshold (modified Gram-Schmidt)
- Aitken warm-up for first time step and early iterations
- Track x_tilde and residual history per time step for building
  difference vectors
- IQN-ILS converges for TS1 but still needs debugging for TS2+

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
- Add XML parameters: IQN_ILS_q (max columns, default 10),
  IQN_ILS_eps (QR filter tolerance, default 1e-2),
  IQN_ILS_warmup (Aitken warm-up iterations, default 5)
- Use V_cols_.size() < warmup for Aitken/IQN-ILS switching
- Persist V/W columns across time steps, clear per-TS history only
- Fix NaN detection in step_equation: only check iNorm, skip first
  Newton iteration where norms are uninitialized
- IQN-ILS oscillates for this strongly-coupled problem; Aitken
  remains default

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
….dat

Swap output destinations so the tabular coupling.dat format (with |disp|)
goes to stdout for easy logging, while the compact solver-style format
goes to histor.dat.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
…e steps

The mesh equation solves for total displacement from the original reference,
but the old code used += to add it to x_ref (which already carried accumulated
displacements from prior steps). After N time steps the fluid mesh was
deformed ~N times too much.

Fix: apply only the incremental displacement (Dn - Do) relative to x_ref,
giving fluid_x = x_ref + (Dn - Do) = x_original + Dn regardless of the
time step number.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Only coupling iteration data is printed to screen. Sub-field Newton
iteration output (NS/ST/MS lines) still goes to each sub-simulation's
history file.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Use the same numbering scheme as NS/ST/MS output lines, with CP as the
equation name and s suffix indicating saved time steps.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
- Use (W - V) * c instead of W * c for the IQN-ILS update, matching
  the secant condition J^{-1} * ΔR = ΔX where ΔX = ΔX̃ - ΔR.
- Compute vel_prev_ from relaxed disp_prev_ via Newmark relationship
  instead of independent relaxation. Applies to all coupling methods.
- Remove spurious omega_ diagnostic from relax_iqn_ils.
- Use constant relaxation (not Aitken) during IQN-ILS warmup.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
- Filter criterion: |R[i,i]| < eps * ||R||_F (L2 norm of R matrix),
  matching svFSGe QRfiltering. More aggressive than the previous
  |R[i,i]| < eps * max|R[j,j]| which caused batch removal.
- Modified back-substitution: set c[i]=0 when |R[i,i]| is small
  instead of dividing by near-zero (svFSGe solve_upper_triangular_mod).
- Newest columns prepended so QR prioritizes recent information.
- Remove debug output.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
- Fix out-of-bounds Dg access in extract_fluid_traction: use com_mod.x
  directly (already deformed in partitioned FSI) instead of x + Dg
  which accesses non-existent mesh DOFs in the fluid sub-sim.
- Remove unused vel_current parameter from all relaxation methods
  (velocity is now computed from relaxed displacement via Newmark).
- Remove unused variables in apply_velocity_on_fluid (dt, gam, v_old,
  a_old, Yo, Ao).
- Remove unused cm/cm_mod in build_node_maps.
- Remove dead stub methods (create_sub_simulations, generate_sub_xml,
  init_sub_simulation).
- Remove legacy use_aitken parameter (replaced by coupling_method).
- Fix al array oversizing in construct_fluid (only yl needs ALE rows).
- Add comments for x_ref semantics and extract_fluid_traction coordinate
  handling.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
- Add newmark::state_from_displacement/velocity to Integrator.h and use
  them in fsi_coupling and PartitionedFSI (eliminates 4 duplicated
  Newmark computations).
- Move extract_fluid_traction to post::compute_face_traction (traction
  computation belongs with post-processing, alongside bpost).
- Move enforce_dirichlet_on_face and enforce_dirichlet_dofs_on_face to
  set_bc namespace (BC enforcement in the assembled system belongs with
  set_bc, alongside set_bc_undef_neu_l).
- Remove extract_solid_velocity (velocity is now always computed from
  displacement via Newmark).
- Clean up fsi_coupling includes (no longer needs fluid.h, fs.h, nn.h,
  all_fun.h).

fsi_coupling now only contains the thin BC-application functions
(apply_velocity_on_fluid, apply_traction_on_solid,
apply_displacement_on_mesh, extract_solid_displacement).

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
IQN-ILS was experimental and not working correctly. Remove the
implementation, parameters, and all references. Only constant and
Aitken relaxation remain as coupling methods.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Extract from the ~280-line step() function:
- compute_interface_velocity(): Newmark velocity from disp_prev_
- solve_fluid(): fluid solve with interface velocity and ALE
- solve_solid(): traction extraction and solid solve
- solve_mesh(): mesh solve and fluid mesh deformation

step() is now a readable coupling loop (~80 lines) that calls
these helpers in sequence.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
PartitionedFSI now sets ale_mesh_velocity via Integrator::set_ale_mesh_velocity()
instead of writing directly to ComMod. The Integrator copies it to ComMod
before assembly and clears it after, so the field in ComMod is only used
as a transient transport mechanism for the assembly functions.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Restore content lost from Simulation.h, Simulation.cpp, and main.cpp due
to a regex bug that swallowed code across conflict boundaries. Update
compute_face_traction to remove explicit Yg/Dg parameters, matching
main's convention where these are retrieved from SolutionStates::intermediate.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
- Broadcast partitioned FSI config from rank 0 to all ranks in
  initialize_partitioned_fsi so slave ranks enter the constructor
- Gather interface face coordinates globally before nearest-neighbor
  node mapping so nodes on remote ranks are found correctly
- Replicate all interface arrays (disp_prev_, vel_prev_, node maps)
  on every rank via MPI_Allgatherv to ensure consistent convergence
  checks and Aitken relaxation across ranks
- Fix Use_Aitken XML element to Coupling_method in solver.xml and
  solver_ramp.xml
The function was only visible to svmultiphysics (main.cpp is excluded
from the run_all_unit_tests build). Moving the definition to
Simulation.cpp and the declaration to Simulation.h makes it available
to both binaries without introducing a new translation unit.
The static local copy conflicted with the external definition now in
Simulation.cpp. The function is declared in Simulation.h so the call
site in setup_simulation resolves correctly without the local copy.
cmType methods (id, np, com) are not marked const, so passing const cmType&
to the static helpers caused a compile error under clang's stricter
const-correctness enforcement. Drop const from the cmType parameters in
compute_face_global_info, gather_face_data, and gather_global_map.
@codecov
Copy link
Copy Markdown

codecov Bot commented May 22, 2026

Codecov Report

❌ Patch coverage is 57.84314% with 602 lines in your changes missing coverage. Please review.
✅ Project coverage is 68.36%. Comparing base (ed4c544) to head (318aed0).

Files with missing lines Patch % Lines
Code/Source/solver/PartitionedFSI.cpp 0.00% 461 Missing ⚠️
Code/Source/solver/Simulation.cpp 31.14% 42 Missing ⚠️
Code/Source/solver/fsi_coupling.cpp 60.37% 21 Missing ⚠️
Code/Source/solver/set_bc.cpp 55.55% 20 Missing ⚠️
...nitTests/integrator_tests/test_partitioned_fsi.cpp 96.52% 11 Missing ⚠️
Code/Source/solver/Parameters.cpp 61.53% 10 Missing ⚠️
.../unitTests/integrator_tests/test_step_equation.cpp 94.00% 9 Missing ⚠️
Code/Source/solver/post.cpp 92.30% 7 Missing ⚠️
Code/Source/solver/Integrator.cpp 92.30% 4 Missing ⚠️
Code/Source/solver/Integrator.h 33.33% 4 Missing ⚠️
... and 6 more
Additional details and impacted files
@@            Coverage Diff             @@
##             main     #554      +/-   ##
==========================================
- Coverage   68.80%   68.36%   -0.45%     
==========================================
  Files         181      186       +5     
  Lines       34157    35559    +1402     
  Branches     5903     6170     +267     
==========================================
+ Hits        23503    24310     +807     
- Misses      10517    11112     +595     
  Partials      137      137              

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:
  • ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

@Eleven7825 Eleven7825 requested a review from mrp089 May 22, 2026 23:43
shiyi added 3 commits May 23, 2026 00:43
Users now write a single solver XML with role="partitioned_fluid/solid/mesh"
on each Add_equation. PartitionedFSI automatically extracts each tagged
equation and its domain-matched meshes into a minimal sub-simulation XML
(written to a temp file) and feeds it to init_sub_sim.

- Parameters.h/cpp: add role attribute to EquationParameters; remove
  Fluid_xml/Solid_xml/Mesh_xml from PartitionedCouplingParameters
- PartitionedFSIConfig: remove fluid_xml/solid_xml/mesh_xml fields
- PartitionedFSI::build_sub_xml: new static helper that parses the main
  XML with tinyxml2, clones the target equation + matching meshes into a
  sub-document, and writes a temp file (cleaned up in destructor)
- Simulation::initialize_partitioned_fsi: active check now looks for any
  equation with a partitioned role; remove sub-XML path broadcasts
- solver.xml/solver_ramp.xml: add role attributes, add Domain block to
  mesh equation, remove sub-XML path elements from Partitioned_coupling
- Override Name_prefix_of_saved_VTK_files per sub-sim in build_sub_xml
  so fluid/solid/mesh each write distinct result_fluid_*, result_solid_*,
  result_mesh_* files instead of overwriting a shared result_* prefix
- Add minimal <Partitioned_coupling> block in mesh sub-XML to set mvMsh flag
- Extend conftest.py run_by_name/run_with_reference with name_result
  parameter to support non-default VTK output filenames
- Add test_pipe_3d_partitioned_fluid and test_pipe_3d_partitioned_solid
  to test_fsi.py comparing against committed reference VTUs
- Relax solver.xml coupling settings (tol 1e-4, relax 0.5) for convergence
- Commit reference result_fluid_010.vtu and result_solid_010.vtu
The standalone mesh sub-sim includes a <Partitioned_coupling> stub so
that read_files accepts the mesh equation (which requires mvMsh=true).
However baf_ini assumes full FSI DOF layout when mvMsh=true and accesses
Do(i+nsd+1,Ac) — row 4 — which is out-of-bounds when tDof=3 for the
mesh-only sub-sim.  This crash was caught by -DENABLE_ARRAY_INDEX_CHECKING
in CI, causing all partitioned FSI tests to fail with no output.

Reset mvMsh=false after read_files in init_sub_sim when the sub-sim
contains only the mesh equation.  Validation (requiring mvMsh for the mesh
equation) already passed at that point, so resetting is safe.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants