-
Notifications
You must be signed in to change notification settings - Fork 15
Add RCR-Windkessel model #935
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
base: main
Are you sure you want to change the base?
Conversation
8096a2d to
bfbf3b6
Compare
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Pull Request Overview
This PR introduces the RCR-Windkessel lumped parameter model for simulating vascular systems at open boundaries. The model captures the relationship between pressure and flow in pulsatile systems using three components: characteristic resistance (R1), peripheral resistance (R2), and compliance (C).
Key changes:
- Added
RCRWindkesselModelstruct and associated calculation methods - Integrated pressure model support into the
OpenBoundarySystemandBoundaryZoneinfrastructure - Added comprehensive validation tests using ferret vascular system data
Reviewed Changes
Copilot reviewed 13 out of 13 changed files in this pull request and generated 4 comments.
Show a summary per file
| File | Description |
|---|---|
src/schemes/boundary/open_boundary/pressure_model.jl |
Core implementation of the RCR-Windkessel model including pressure and flow rate calculations |
src/schemes/boundary/open_boundary/system.jl |
Extended OpenBoundarySystem to support pressure models with new fields and initialization logic |
src/schemes/boundary/open_boundary/boundary_zones.jl |
Added pressure_model field to BoundaryZone with validation logic |
src/schemes/boundary/open_boundary/dynamical_pressure.jl |
Modified pressure reference to support pressure model computations |
src/schemes/boundary/open_boundary/mirroring.jl |
Updated boundary quantity handling to work with pressure models |
src/callbacks/update.jl |
Modified update callback to pass timestep to open boundary updates |
test/schemes/boundary/open_boundary/pressure_model.jl |
Comprehensive tests including pulsatile flow validation against ODE solution |
docs/src/systems/boundary.md |
Added documentation section for pressure models |
docs/src/refs.bib |
Added references for Gasser (2021) and Westerhof (2008) |
| Other files | Minor updates for file inclusion and exports |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
|
/run-gpu-tests |
| end | ||
| end | ||
|
|
||
| function extract_pressure_models(boundary_zones) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I don't understand why this has to be so complicated.
It looks like this model can only be used if reference_pressure = nothing is passed to the BoundaryZone, so it is essentially just a complex function for reference_pressure. It is also only used where a reference_pressure function would be used (imposed_pressure) plus it is updated. These are the only two occurrences.
So why do you need this extra struct and extra field in the boundary zone? Why can't you just make this a functor? Basically, just add a field pressure to the struct and make it a functor as
function (model::RCRWindkesselModel)(x, t)
return model.pressure[]
endAnd then the only thing you have to change in the open boundary system and boundary zones files is that you add an update. This update function updates reference_pressure for each boundary zone. If it's ::RCRWindkesselModel, do the update, otherwise do nothing.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'm not sure I understood correctly. What happens if some BoundaryZones use a pressure_model while others have a prescribed pressure (as a function or scalar)?
Is this then type stable?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Also, IIRC, storing RefValues in a GPU array is not allowed.
| previous_flow_rate = flow_rate[] | ||
| flow_rate[] = current_flow_rate | ||
|
|
||
| # Calculate new pressure according to eq. 22 in Zhang et al. (2025) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'm confused. What do Zhang et al. have to do with this? You listed different references.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
When I implemented this, I based it on this literature because it was also mentioned in the SPH context. It's a simple ODE. I don't really care which literature I cite here...
The other references are purely in a biomechanical context where Zhang is in a SPH context.
| # dp[1] = -p[1] / (R_2 * C) + R_1 * dq_dt + q_func(t) * (R_2 + R_1) / (R_2 * C) | ||
| # end | ||
| # | ||
| # The solution is obtained using simple Euler integration over the time interval `tspan`. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Why Euler? You could just use the default, which is an adaptive method that doesn't even need a time step.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
and then? It doesn't makes a difference here, since I only needed this for the pressure reference values
Lumped parameter model to simulate vascular systems