Issue 512 3D-1D coupling#563
Conversation
…larify time argument
…nTotalModels, improve comment
…el initialization
…set_flowrates() for DIR
…et_bc_neu_l (pressure)
…D), fixes svZeroD DIR
…e DOFs from linear solve For 1D Dirichlet coupling, `iBC_Dir` is cleared in read_files.cpp so that set_bc_cpl routes the BC correctly. However this caused fsi_ls_ini to skip registering the face with the FSILS linear solver as a Dirichlet (BC_TYPE_Dir) constraint. As a result, the preconditioner did not zero out those rows/columns, and the linear solver overwrote the velocity values applied by set_bc_dir with the NS solution, producing a lower-than-expected centerline velocity (~47 vs ~63.66). Fix: detect Coupled-DIR faces in fsi_ls_ini and register them as BC_TYPE_Dir so their DOFs are properly excluded from Ax=b.
- CouplingInterfaceParameters: add Coupling_ramp_steps (int) and Coupling_ramp_ref_pressure (double) XML parameters - CoupledBoundaryCondition: add oned_ramp_steps_ / oned_ramp_ref_pressure_ private members with set_oned_ramp() / get_oned_ramp_steps() / get_oned_ramp_ref_pressure() accessors; propagate through all copy/move constructors and assignment operators - read_files.cpp: read ramp params from CouplingInterfaceParameters and store them in lBc.coupled_bc via set_oned_ramp() - svOneD_subroutines.cpp: OneDModelState gains ramp_steps, ramp_ref_pressure, step_count fields; init_svOneD reads them from coupled_bc; calc_svOneD applies linear pressure ramp for DIR coupling; step_count is incremented on every committed (BCFlag=='L') step
Implements exponential smoothing of the pressure boundary condition sent to the 1D solver (DIR coupling only): P_sent = omega * P_target + (1 - omega) * P_prev_sent where P_target has already been through the ramp filter. - Parameters.h/cpp: add Coupling_dir_relax_factor parameter (default 1.0) - CoupledBoundaryCondition.h/cpp: add oned_relax_factor_ field, getter, setter, copy/move and distribute() broadcast - read_files.cpp: parse and apply the new parameter - svOneD_subroutines.cpp: add relax_factor / P_prev_sent_old/new to OneDModelState; apply under-relaxation after ramp; update history only on BCFlag=='L' committed steps; NEU coupling untouched
…at the same model
|
@taeoukkim Pleas request some reviews for this PR too |
There was a problem hiding this comment.
Pull request overview
Implements 3D–1D coupling via a dynamically loaded svOneDSolver interface and extends the existing 3D–0D coupling infrastructure to support mixed 0D/1D coupling per-face, including coexistence with RCR boundaries.
Changes:
- Added svOneDSolver dynamic-library interface and coupling subroutines, with MPI-parallel initialization/stepping per coupled face.
- Updated BC parsing/routing and BC application to support 1D coupling (Dirichlet + Neumann) and mixed 0D/1D configurations, including adjusted handling for DIR-coupled faces in the linear system/BC application.
- Extended CoupledBoundaryCondition to carry 1D metadata plus ramp/relax history; added an example coupling case under
tests/cases.
Reviewed changes
Copilot reviewed 18 out of 18 changed files in this pull request and generated 7 comments.
Show a summary per file
| File | Description |
|---|---|
| tests/cases/fluid/pipe_svOneD_2faces/solver.xml | Adds an example/test-case solver.xml for 3D–1D coupling on two faces. |
| Code/Source/solver/txt.cpp | Updates end-of-step coupling integration to support mixed svZeroD/svOneD and RCR coexistence. |
| Code/Source/solver/svZeroD_interface.cpp | Updates svZeroD coupled-face indexing (exclude 1D faces) and adds ramp/relax history handling + DIR-flow broadcast. |
| Code/Source/solver/svOneD_subroutines.h | Declares svOneD init/step entry points. |
| Code/Source/solver/svOneD_subroutines.cpp | Implements 3D–1D coupling orchestration, MPI ownership, stepping, and relaxation/ramping. |
| Code/Source/solver/svOneD_interface/OneDSolverInterface.h | Adds a dlopen/dlsym wrapper API for the svOneDSolver interface library. |
| Code/Source/solver/svOneD_interface/OneDSolverInterface.cpp | Implements dynamic symbol loading and wrapper calls into the 1D shared library. |
| Code/Source/solver/set_bc.cpp | Updates coupled BC derivative evaluation, mixed solver calling, RCR integration indexing, and DIR-coupled handling in Dir/Neu application. |
| Code/Source/solver/read_files.cpp | Extends XML parsing/routing for svOneD per-face input files and enables mixed 0D/1D coupling configuration rules. |
| Code/Source/solver/Parameters.h | Adds svOneDSolver interface parameter list + coupling-interface parameters for per-face 1D input and ramp/relax knobs. |
| Code/Source/solver/Parameters.cpp | Implements XML parameter parsing for svOneD interface + coupling-interface additions. |
| Code/Source/solver/distribute.cpp | Broadcasts svOneD interface metadata and ensures cplBC state is distributed correctly with RCR coexistence. |
| Code/Source/solver/CoupledBoundaryCondition.h | Extends CoupledBoundaryCondition with svOneD input file + ramp/relax state and adds new cap-surface exception types. |
| Code/Source/solver/CoupledBoundaryCondition.cpp | Implements distribution of new fields and adds additional cap-surface robustness checks + DIR flow broadcast helper. |
| Code/Source/solver/ComMod.h | Adds svOneD interface type and new flags/fields in coupling-related structs. |
| Code/Source/solver/ComMod.cpp | Implements svOneDSolverInterfaceType::set_data. |
| Code/Source/solver/CMakeLists.txt | Adds new svOneD sources to the solver build. |
| Code/Source/solver/baf_ini.cpp | Initializes svOneD coupling during solver init and tags RCR faces for integration. |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
|
|
||
| // Coupled BC class | ||
| CoupledBoundaryCondition coupled_bc; | ||
|
|
||
| // svOneD: per-face 1D solver input file path (set when Time_dependence=Coupled | ||
| // and svOneDSolver_interface is active). | ||
| std::string oned_input_file; | ||
| }; |
There was a problem hiding this comment.
I kept these members for now because they are used as intermediate storage
There was a problem hiding this comment.
I'm not sure I understand this: the AI comment says that these members are written to but never read from. Do you mean to say that they are indeed read from, and that AI missed it? Or that you plan on reading from them in further developments of this work?
Codecov Report❌ Patch coverage is Additional details and impacted files@@ Coverage Diff @@
## main #563 +/- ##
==========================================
- Coverage 69.44% 68.66% -0.79%
==========================================
Files 181 184 +3
Lines 34072 34630 +558
Branches 5930 6018 +88
==========================================
+ Hits 23662 23778 +116
- Misses 10273 10715 +442
Partials 137 137 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
aabrown100-git
requested review from
aabrown100-git,
javijv4,
ktbolt and
michelebucelli
ktbolt
left a comment
ktbolt
left a comment
There was a problem hiding this comment.
@taeoukkim A couple of general comments
- C++ does not have subroutines so rename
svOneD_subroutines*files tosvOneD_interface* - I see a lot of
sv1Dprefixes in the code; change these tosvOneD
michelebucelli
left a comment
michelebucelli
left a comment
There was a problem hiding this comment.
Thank you @taeoukkim, I have left some comments. I do apologize for the slight pedantry 😅 , I am happy to further discuss any of the comments I left, even over a call if necessary.
A couple of general comments:
-
as I pointed out in a few of my suggestions, I think that it would be useful to document a bit better how this works. This can be done by adding some details in the Doxygen comments and/or expanding those comments themselves. In some cases it might be enough to just point to the paper where the precise model/method used is implemented, as long as there is no ambiguity about it.
-
I think it might be useful to make an effort to further encapsulate the 1D-solver-related stuff into a coherent whole, and try to condense them into as few files as possible (and/or a few files placed in a dedicated folder). For example, some of the functions reside in the
svOneDnamespace, but many others do not, and I see no reason why they shouldn't.
| /// @param system_size Output: total number of DOFs (nodes * 2: flow + area). | ||
| /// @param coupling_type "NEU" or "DIR" coupling direction. | ||
| void initialize(const std::string& input_file, int& problem_id, | ||
| int& system_size, const std::string& coupling_type); |
There was a problem hiding this comment.
I would recommend against using std::string for representing an object that has only a finite set of values. A dedicated enum class would be more appropriate (e.g. for type safety).
There was a problem hiding this comment.
I agree that an enum class would provide better type safety here. I kept the string representation in this PR to minimize changes and remain consistent with the existing interface.
There was a problem hiding this comment.
@taeoukkim I don't think we need to maintain consistency between the 0D and 1D interface code; better to improve the 1D interface code for readability and robustness.
There was a problem hiding this comment.
The svOneDSolver interface currently expects the coupling type as a string ("NEU" or "DIR"), so I kept std::string in the wrapper to match that interface.
There was a problem hiding this comment.
I see your point, but I would argue that the use of std::string to select this is a slightly bad design on the svOneDSolverSide, and while it might be worth improving there, I think it's a good idea to hide it behind a more robust interface on the svmp side in the meantime.
An implementation might look like this:
// ...somewhere in the header file...
/// @brief Type of coupling with svOneDSolver.
enum class svOneDCoupling {
/// Documentation on what this type of coupling means.
Dirichlet,
/// Documentation on what this type of coupling means.
Neumann
}
// ...in the cpp file...
void initialize(const std::string& input_file, int& problem_id,
int& system_size, const svOneDCoupling coupling_type) {
std::string coupling_type_str;
if (coupling_type == svOneDCouplingType::Dirichlet) {
coupling_type_str = "DIR";
} else if (coupling_type == svOneDCouplingType::Neumann) {
coupling_type_str = "NEU";
} else {
svmp::raise<svmp::FE::NotImplementedException>(
SVMP_HERE,
"Unsupported type of coupling to svOneDSolver");
}
// ...do whatever comes next...
}This code is way more type safe and robust to change than the version using std::string, which makes the overall implementation easier to maintain and build upon.
| const std::string& coupling_type, double* params, | ||
| double* solution, double& cpl_value, char last_flag, |
There was a problem hiding this comment.
Also similar to before, I suggest that both coupling_type and last_flag be turned into enum classes. The same might be done for error_code as well.
There was a problem hiding this comment.
I agree that an enum class would provide better type safety here. I kept the string representation in this PR to minimize changes and remain consistent with the existing interface
There was a problem hiding this comment.
Similar comment as above.
| // A Coupled face belongs to svZeroD only when it is not a 1D face. | ||
| // The is_sv1d_face() helper encapsulates the routing invariant: | ||
| // 1D face → oned_input_file_ non-empty, block_name_ empty | ||
| // 0D face → block_name_ non-empty, oned_input_file_ empty | ||
| if (utils::btest(bc.bType, consts::iBC_Coupled) && | ||
| !bc.coupled_bc.is_sv1d_face()) { |
There was a problem hiding this comment.
Perhaps due to my unfamiliarity with this part of the code, I'm not sure I understand this comment: is there a situation where the user might intentionally specify the same face to be coupled to both the 0D and 1D solvers?
There was a problem hiding this comment.
Correct, the intended use case is not to couple the same face to both the 0D and 1D solvers. A single face should be routed to either the 0D interface or the 1D interface.
What this logic allows is a mixed model where different faces of the same 3D model can be coupled to different reduced-order models. For example, one outlet face can be coupled to svZeroD, while another outlet face can be coupled to svOneD. The is_svOneD_face() check prevents the 1D-coupled faces from also being handled by the svZeroD interface.
There was a problem hiding this comment.
I see, thanks for explaining.
So if I understand correctly the user now could assign (most likely by mistake) the same face to both the 0D and the 1D interface, and when this line is reached the 1D takes precedence over the 0D, and the 0D coupling that the user declared is silently disregarded. Is this accurate?
If so, I do not think that this is a desirable behavior. I think what should happen is that if the user does accidentally assign the same face to 0D and 1D coupling, an error should be raised, with a clear message (and ideally immediately after parsing the XML file).
If instead there already is something that prevents the user from assigning two coupling BCs to the same face, the check !bc.coupled_bc.is_sv1d_face() is redundant and can probably be removed.
Is there something else that I am missing?
javijv4
left a comment
javijv4
left a comment
There was a problem hiding this comment.
Hey @taeoukkim, thanks for working on this! I left some comments.
One important thing, the test case for the 1D implementation does not have all the necessary files to run. Can you update that?
Also, can you show tests for the Dirichlet implementation? I'm not sure whether they should be CI tests, but I think it would be helpful to show how to use the new functionality in the respective GitHub issues.
|
The svZeroD coupling code is not a good object-oriented implementation; had hoped the svOneD code could be implemented with a better design, perhaps implement an external solver abstract class of some sort, and then refactor the svZeroD code. @taeoukkim has used the current implementation as a template which was the most straight-forward thing to do. But now to implement a more object-oriented code we would need to refactor two codes; this should be a separate PR. |
ktbolt
left a comment
ktbolt
left a comment
There was a problem hiding this comment.
@taeoukkim Thank you for squashing your commits; much easier to review !
| /// @param system_size Output: total number of DOFs (nodes * 2: flow + area). | ||
| /// @param coupling_type "NEU" or "DIR" coupling direction. | ||
| void initialize(const std::string& input_file, int& problem_id, | ||
| int& system_size, const std::string& coupling_type); |
There was a problem hiding this comment.
@taeoukkim I don't think we need to maintain consistency between the 0D and 1D interface code; better to improve the 1D interface code for readability and robustness.
There was a problem hiding this comment.
Thanks for the many improvements @taeoukkim! I have done another pass (and also reopened some threads from previous pass)
| /// (pressure for NEU, flow for DIR). | ||
| /// @param last_flag 'L' for the final (committed) iteration, 'D' for | ||
| /// derivative / predictor steps. | ||
| /// @param error_code Output: non-zero on failure. |
There was a problem hiding this comment.
Se my other comment for my thoughts on using enum class. Having said that, even if you want to keep the int error code here (actually, especially if you want to keep it), it would be a good idea to document the meaning of the possible error codes.
|
|
||
| // Coupled BC class | ||
| CoupledBoundaryCondition coupled_bc; | ||
|
|
||
| // svOneD: per-face 1D solver input file path (set when Time_dependence=Coupled | ||
| // and svOneDSolver_interface is active). | ||
| std::string oned_input_file; | ||
| }; |
There was a problem hiding this comment.
I'm not sure I understand this: the AI comment says that these members are written to but never read from. Do you mean to say that they are indeed read from, and that AI missed it? Or that you plan on reading from them in further developments of this work?
| } | ||
| } | ||
|
|
||
| void OneDSolverInterface::load_library(const std::string& interface_lib) |
There was a problem hiding this comment.
Yes, I agree that svOneDSolver should remain optional (same as svZeroDSolver). CMake configuration would provide a way to make this optional dependency more systematic (and possibly static instead of dynamic, but that's a bigger issue).
I do agree that it's maybe a bigger refactoring than this PR is concerned with, but I do think it's something worth considering.
| throw svmp::CoreException( | ||
| std::string("[OneDSolverInterface] Could not load shared library '") + | ||
| interface_lib + "': " + dlerror(), | ||
| svmp::StatusCode::DependencyError, | ||
| __FILE__, | ||
| __LINE__, | ||
| __func__); |
There was a problem hiding this comment.
The exception throwing syntax that @zasexton set up uses the svmp::raise function. With it, this would look like
svmp::raise<svmp::DependencyException>(
SVMP_HERE,
"[OneDSolverInterface] Could not load shared library" + interface_lib + ": " + dlerror());which is a bit more compact, readable, and better adheres to the exception conventions established in #537.
I suggest using this syntax here and elsewhere.
| // For DIR Coupled BCs (svZeroD or svOneD), the downstream solver always returns a | ||
| // volumetric flow rate Q [m³/s], not a velocity [m/s]. Without bType_flx, bc_ini | ||
| // sets gx(a) = 1, and set_bc_dir_l applies velocity = Q * nV which has the wrong | ||
| // dimensions. With bType_flx, bc_ini normalises gx(a) = 1/area, so the applied | ||
| // velocity = (Q/area) * nV [m/s] is physically correct. Enforce this automatically | ||
| // regardless of the user's <impose_flux> setting. |
There was a problem hiding this comment.
Similar to before, I am a bit wary of overriding the user's settings without telling them. There might be a situation where the user is convinced they're doing something, while in fact they're doing something else without ever knowing it (and maybe end up publishing something incorrect, or unintended, by accident).
I think it is better to throw an error message if the user chooses a configuration that is incorrect, with the error message telling them what to do.
In this case, this would translate to checking whether bType_flx has been set, and throwing an exception if it hasn't saying that it needs to be set (and why).
Current situation
resolves #512
Release Notes
New features:
(1) 3D–1D coupling for both Dirichlet and Neumann boundary conditions
(2) RCR boundary condition support with 0D/1D coupling interfaces
(3) Dirichlet-type 3D–0D coupling implemented
(4) Enabled simultaneous use of 1D and 0D coupling in a single model
Documentation
I will also update documentation.
Testing
Successfully compiled on macOS.
Tested with 3D–1D coupling cases.
Tested with mixed 3D–1D and 3D–0D coupling configurations.
Code of Conduct & Contributing Guidelines