Skip to content

Record immersed-boundary forces and kinematics every time step - #1846

Open
sbryngelson wants to merge 2 commits into
masterfrom
feat/ib-force-history
Open

Record immersed-boundary forces and kinematics every time step#1846
sbryngelson wants to merge 2 commits into
masterfrom
feat/ib-force-history

Conversation

@sbryngelson

Copy link
Copy Markdown
Member

Why

ib_state_wrt records the force, torque and kinematic state of each immersed boundary only at snapshot intervals, so the force history is sampled at the field-output cadence — typically every few hundred steps. That is too coarse to compare a transient load against an experiment or a reference computation, or to drive a reduced-order model from a force signal. Every force comparison in the verification work behind this PR needed a per-step record.

What

One record per time step in D/ib<id>_forces.dat: time step, time, force, torque, velocity, angular velocity, angles, centroid. Gated by the existing ib_state_wrt, so no new switch to enable it.

Scaling

The obvious implementation — open, append, close, per body, per step — is a metadata operation per step on a parallel filesystem and does not scale. A particle bed of a thousand bodies over a long run would issue on the order of a hundred million of them. Instead records accumulate in a small fixed rank-local buffer (1024 rows) and are flushed in batches, with a final flush at shutdown. Each buffered row carries its own global body id, so a body changing owner mid-run needs no special handling: the previous owner's pending rows still reach the right file.

ib_force_stride (default 1) subsamples runs long enough for the record itself to become large. At the default the record is about 350 bytes per step per body — 25 MB for a 50k-step case, which is nothing next to the field output, but a thousand-body bed would want a stride.

Testing

Exercised on Frontier across 2D and 3D immersed-boundary cases, static and moving, on 4, 8 and 64 ranks. The per-step forces were checked against an independent control-volume momentum balance and, for a static cylinder at Re 40, against the published drag coefficient. ./mfc.sh precheck passes.

https://claude.ai/code/session_01HMJ7cycfo7kTFSFq5yhHLG

`ib_state_wrt` writes the force, torque and kinematic state of each immersed boundary only
at snapshot intervals, so the force history is sampled at the field-output cadence. That is
far too coarse to compare a transient load against an experiment or a reference
computation, or to drive a reduced-order model from a force signal: typical cases here
write a field every few hundred steps.

Write one record per time step to `D/ib<id>_forces.dat` (time step, time, force, torque,
velocity, angular velocity, angles, centroid), under the existing `ib_state_wrt` flag.

Records are buffered per rank and flushed in batches rather than opened per body per step:
opening a file per body per step is a metadata operation per step on a parallel filesystem
and does not scale, and a particle bed of a thousand bodies would issue on the order of a
hundred million of them over a long run. Each buffered row carries its own global body id,
so a body changing owner mid-run needs no special handling, and `ib_force_stride`
subsamples runs long enough for the record itself to become large.

Claude-Session: https://claude.ai/code/session_01HMJ7cycfo7kTFSFq5yhHLG
Copilot AI lite review requested due to automatic review settings September 11, 2026 14:42

Copilot AI left a comment

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Warning

Copilot couldn't run its full agentic review because it didn't start before the timeout. Make sure your repository has a runner available, or add a copilot-code-review.yml file specifying one with the runs-on attribute. See the docs for more details.

Pull request overview

Adds per-time-step immersed-boundary force/kinematics logging (buffered + flushed in batches) so IB force histories are recorded at high temporal resolution, with optional subsampling via ib_force_stride.

Changes:

  • Introduces ib_force_stride parameter (docs + toolchain registration + validation).
  • Writes buffered per-step IB force/torque/kinematics rows and flushes at shutdown.
  • Hooks per-step recording into the time-step loop when ib_state_wrt is enabled.

Reviewed changes

Copilot reviewed 8 out of 8 changed files in this pull request and generated 6 comments.

Show a summary per file
File Description
toolchain/mfc/params/descriptions.py Adds description string for ib_force_stride.
toolchain/mfc/params/definitions.py Registers ib_force_stride in parameter definitions.
toolchain/mfc/case_validator.py Validates ib_force_stride value in IBM checks.
src/simulation/m_time_steppers.fpp Calls per-step IB force logging routine.
src/simulation/m_start_up.fpp Flushes buffered IB force records at shutdown.
src/simulation/m_global_parameters.fpp Sets default ib_force_stride = 1.
src/simulation/m_data_output.fpp Implements buffered IB force record writing + flush to D/ib<id>_forces.dat.
docs/documentation/case.md Documents ib_force_stride and new per-step IB force record behavior.

💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.

Comment on lines +1160 to +1167
write (file_loc, '(A,I0,A)') '/D/ib', ib_id, '_forces.dat'
file_loc = trim(case_dir) // trim(file_loc)
inquire (file=trim(file_loc), exist=file_exist)
if (file_exist) then
open (newunit=file_unit, file=trim(file_loc), form='formatted', status='old', position='append')
else
open (newunit=file_unit, file=trim(file_loc), form='formatted', status='new')
write (file_unit, '(A)') '# t_step time Fx Fy Fz Tx Ty Tz vx vy vz wx wy wz ax ay az xc yc zc'
Comment on lines +46 to +48
!> Buffered immersed-boundary force records: (id, t_step, time, force, torque, vel, angular_vel, angles, centroid)
integer, parameter :: ib_force_buf_len = 1024
real(wp), dimension(21, ib_force_buf_len) :: ib_force_buf
Comment on lines +1127 to +1131
$:GPU_UPDATE(host='[patch_ib(1:num_ibs)]')

do i = 1, n_write
ib_idx = i
if (num_procs > 1) ib_idx = local_ib_patch_ids(i)
)
self.prohibit(not ib and num_ibs > 0, "num_ibs is set, but ib is not enabled")
self.prohibit(ib_state_wrt and not ib, "ib_state_wrt requires ib to be enabled")
ib_force_stride = self.get("ib_force_stride", 1) or 1
integer :: i, j, ib_id, file_unit
logical :: file_exist

do i = 1, ib_force_buf_n
write (file_unit, '(A)') '# t_step time Fx Fy Fz Tx Ty Tz vx vy vz wx wy wz ax ay az xc yc zc'
end if

do j = i, ib_force_buf_n ! all rows for this body, in time order
ES18.10E3 is exactly wide enough for a negative value with a three-digit exponent
(-1.2345678901E+003), so adjacent values were written with no space between them and the
file could not be read as whitespace-separated columns. Use an explicit 1X separator.

Claude-Session: https://claude.ai/code/session_01HMJ7cycfo7kTFSFq5yhHLG
@github-actions

Copy link
Copy Markdown

Claude Code Review

Head SHA: 0d79d0e

Files changed:

  • 8
  • docs/documentation/case.md
  • src/simulation/m_data_output.fpp
  • src/simulation/m_global_parameters.fpp
  • src/simulation/m_start_up.fpp
  • src/simulation/m_time_steppers.fpp
  • toolchain/mfc/case_validator.py
  • toolchain/mfc/params/definitions.py
  • toolchain/mfc/params/descriptions.py

Findings:

  • toolchain/mfc/case_validator.py (check_ibm): ib_force_stride = self.get("ib_force_stride", 1) or 1 silently coerces a user-supplied 0 into 1 before the very next line checks ib_force_stride < 1. Since 0 is falsy in Python, the or 1 fallback fires for the exact invalid value (0) this check exists to catch, so ib_force_stride = 0 passes validation instead of being prohibited. Negative values are still caught (they're truthy), so only 0 slips through silently. Use self.get("ib_force_stride", 1) without the or 1 fallback (or check for None explicitly) so 0 still hits the < 1 prohibit.
  • src/simulation/m_time_steppers.fpp: s_write_ib_force_files(t_step) is called unconditionally whenever ib_state_wrt is true (line ~480), but per the adjacent existing code (if (moving_immersed_boundary_flag) ... else if (ib_state_wrt) call s_compute_ib_forces(...), lines ~602-609) and the doc text this PR itself edits in case.md ("When no IBs are moving, it also triggers force and torque calculation..."), patch_ib(...)%force/%torque are only refreshed when IBs are not moving. When moving_immersed_boundary_flag is true, the new D/ib<id>_forces.dat writer will keep recording stale/unrefreshed force and torque values every step instead of the actual current-step forces, with no indication in the output that the data is invalid for that run configuration.

@github-actions

Copy link
Copy Markdown

Lines of Code

File Lines Diff
src/simulation/m_data_output.fpp 1422 +55
src/simulation/m_global_parameters.fpp 770 +1
src/simulation/m_start_up.fpp 1250 +1
src/simulation/m_time_steppers.fpp 838 +1
Directory Lines Diff
simulation 28095 +58
total 46320 +58

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Development

Successfully merging this pull request may close these issues.

2 participants