svZeroDSolver
Loading...
Searching...
No Matches
ClosedLoopHeartPulmonary.h
Go to the documentation of this file.
1// Copyright (c) Stanford University, The Regents of the University of
2// California, and others.
3//
4// All Rights Reserved.
5//
6// See Copyright-SimVascular.txt for additional details.
7//
8// Permission is hereby granted, free of charge, to any person obtaining
9// a copy of this software and associated documentation files (the
10// "Software"), to deal in the Software without restriction, including
11// without limitation the rights to use, copy, modify, merge, publish,
12// distribute, sublicense, and/or sell copies of the Software, and to
13// permit persons to whom the Software is furnished to do so, subject
14// to the following conditions:
15//
16// The above copyright notice and this permission notice shall be included
17// in all copies or substantial portions of the Software.
18//
19// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS
20// IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED
21// TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
22// PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
23// OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
24// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
25// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
26// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
27// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
28// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
29// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
30/**
31 * @file ClosedLoopHeartPulmonary.h
32 * @brief model::ClosedLoopHeartPulmonary source file
33 */
34#ifndef SVZERODSOLVER_MODEL_CLOSEDLOOPHEARTPULMONARY_HPP_
35#define SVZERODSOLVER_MODEL_CLOSEDLOOPHEARTPULMONARY_HPP_
36
37#include "Block.h"
38#include "SparseSystem.h"
39
40/**
41 * @brief Define math constants (for M_PI)
42 *
43 */
44#define _USE_MATH_DEFINES
45
46#include <cmath>
47
48/**
49 * @brief Heart and pulmonary circulation model
50 *
51 * Models the mechanics of the 4 heart chambers and pulmonary circulation
52 *
53 * References: \cite sankaran2012patient and \cite menon2023predictors
54 *
55 * TODO: Equations and circuit diagram
56 *
57 * ### Parameters
58 *
59 * Parameter sequence for constructing this block
60 *
61 * * `0` Atrial systole time fraction
62 * * `1` Time for P-wave
63 * * `2` Scaling for right ventricle elastance
64 * * `3` Scaling for left ventricle elastance
65 * * `4` Scaling for intramyocardial pressure (left coronaries)
66 * * `5` Scaling for intramyocardial pressure (right coronaries)
67 * * `6` Right atrium inductance
68 * * `7` Right atrium outflow resistance
69 * * `8` Right ventricle inductance
70 * * `9` Right ventricle outflow resistance
71 * * `10` Left atrium inductance
72 * * `11` Left atrium outflow resistance
73 * * `12` Left ventricle inductance
74 * * `13` Left ventricle outflow resistance
75 * * `14` Right ventricle unstressed volume
76 * * `15` Left ventricle unstressed volume
77 * * `16` Pulmonary resistance
78 * * `17` Pulmonary capacitance
79 * * `18` Aortic capacitance
80 * * `19` Right atrium pressure scaling
81 * * `20` Right atrium volume scaling
82 * * `21` Left atrium pressure scaling
83 * * `22` Left atrium volume scaling
84 * * `23` Right atrium elastance
85 * * `24` Left atrium elastance
86 * * `25` Right atrium resting volume
87 * * `26` Left atrium resting volume
88 *
89 * ### Internal variables
90 *
91 * Names of internal variables in this block's output:
92 *
93 * * `V_RA`: Right atrium volume
94 * * `Q_RA`: Right atrium outflow
95 * * `P_RV`: Right ventricle pressure
96 * * `V_RV`: Right ventricle volume
97 * * `Q_RV`: Right ventricle outflow
98 * * `P_pul`: Pulmonary pressure
99 * * `P_LA`: Left atrium pressure
100 * * `V_LA`: Left atrium volume
101 * * `Q_LA`: Left atrium outflow
102 * * `P_LV`: Left ventricle pressure
103 * * `V_LV`: Left ventricle volume
104 * * `Q_LV`: Left ventricle outflow
105 *
106 */
108 public:
109 /**
110 * @brief Construct a new ClosedLoopHeartPulmonary object
111 *
112 * @param id Global ID of the block
113 * @param model The model to which the block belongs
114 */
116 : Block(id, model, BlockType::closed_loop_heart_pulmonary,
117 BlockClass::closed_loop,
118 {{"Tsa", InputParameter()}, {"tpwave", InputParameter()},
119 {"Erv_s", InputParameter()}, {"Elv_s", InputParameter()},
120 {"iml", InputParameter()}, {"imr", InputParameter()},
121 {"Lra_v", InputParameter()}, {"Rra_v", InputParameter()},
122 {"Lrv_a", InputParameter()}, {"Rrv_a", InputParameter()},
123 {"Lla_v", InputParameter()}, {"Rla_v", InputParameter()},
124 {"Llv_a", InputParameter()}, {"Rlv_ao", InputParameter()},
125 {"Vrv_u", InputParameter()}, {"Vlv_u", InputParameter()},
126 {"Rpd", InputParameter()}, {"Cp", InputParameter()},
127 {"Cpa", InputParameter()}, {"Kxp_ra", InputParameter()},
128 {"Kxv_ra", InputParameter()}, {"Kxp_la", InputParameter()},
129 {"Kxv_la", InputParameter()}, {"Emax_ra", InputParameter()},
130 {"Emax_la", InputParameter()}, {"Vaso_ra", InputParameter()},
131 {"Vaso_la", InputParameter()}}) {}
132
133 /**
134 * @brief Local IDs of the parameters
135 *
136 */
137 enum ParamId {
138 TSA = 0, ///< Fractions of cardiac cycle (not sure)
139 TPWAVE = 1, ///< Fraction of cardiac cycle (P-wave)
140 ERV_S = 2, ///< Scaling for right ventricle elastance
141 ELV_S = 3, ///< Scaling for left ventricle elastance
142 IML = 4, ///< Scaling for intramyocardial pressure (left coronaries)
143 IMR = 5, ///< Scaling for intramyocardial pressure (right coronaries)
144 LRA_V = 6, ///< Right atrium inductance
145 RRA_V = 7, ///< Right atrium outflow resistance
146 LRV_A = 8, ///< Right ventricle inductance
147 RRV_A = 9, ///< Right ventricle outflow resistance
148 LLA_V = 10, ///< Left atrium inductance
149 RLA_V = 11, ///< Left atrium outflow resistance
150 LLV_A = 12, ///< Left ventricle inductance
151 RLV_AO = 13, ///< Left ventricle outflow resistance
152 VRV_U = 14, ///< Right ventricle unstressed volume
153 VLV_U = 15, ///< Left ventricle unstressed volume
154 RPD = 16, ///< Pulmonary resistance
155 CP = 17, ///< Pulmonary capacitance
156 CPA = 18, ///< Aortic capacitance
157 KXP_RA = 19, ///< Right atrium pressure-volume relationship (?)
158 KXV_RA = 20, ///< Right atrium pressure-volume relationship (?)
159 KXP_LA = 21, ///< Left atrium pressure-volume relationship (?)
160 KXV_LA = 22, ///< Left atrium pressure-volume relationship (?)
161 EMAX_RA = 23, ///< Right atrium elastance (?)
162 EMAX_LA = 24, ///< Left atrium elastance (?)
163 VASO_RA = 25, ///< Right atrium rest volume (?)
164 VASO_LA = 26, ///< Left atrium rest volume (?)
165 };
166
167 /**
168 * @brief Set up the degrees of freedom (DOF) of the block
169 *
170 * Set \ref global_var_ids and \ref global_eqn_ids of the element based on the
171 * number of equations and the number of internal variables of the
172 * element.
173 *
174 * @param dofhandler Degree-of-freedom handler to register variables and
175 * equations at
176 */
177 void setup_dofs(DOFHandler &dofhandler);
178
179 /**
180 * @brief Update the constant contributions of the element in a sparse
181 system
182 *
183 * @param system System to update contributions at
184 * @param parameters Parameters of the model
185 */
186 void update_constant(SparseSystem &system, std::vector<double> &parameters);
187
188 /**
189 * @brief Update the time-dependent contributions of the element in a sparse
190 * system
191 *
192 * @param system System to update contributions at
193 * @param parameters Parameters of the model
194 */
195 void update_time(SparseSystem &system, std::vector<double> &parameters);
196
197 /**
198 * @brief Update the solution-dependent contributions of the element in a
199 * sparse system
200 *
201 * @param system System to update contributions at
202 * @param parameters Parameters of the model
203 * @param y Current solution
204 * @param dy Current derivate of the solution
205 */
206 void update_solution(SparseSystem &system, std::vector<double> &parameters,
207 const Eigen::Matrix<double, Eigen::Dynamic, 1> &y,
208 const Eigen::Matrix<double, Eigen::Dynamic, 1> &dy);
209
210 /**
211 * @brief Modify the solution after solving it
212 *
213 * @param y Current solution
214 */
215 void post_solve(Eigen::Matrix<double, Eigen::Dynamic, 1> &y);
216
217 /**
218 * @brief Number of triplets of element
219 *
220 * Number of triplets that the element contributes to the global system
221 * (relevant for sparse memory reservation)
222 */
224
225 private:
226 // Below variables change every timestep and are then combined with
227 // expressions that are updated with solution
228 double AA; // Atrial activation function
229 double Elv; // LV elastance
230 double Erv; // RV elastance
231 double psi_ra, psi_la, psi_ra_derivative,
232 psi_la_derivative; // Expressions for atrial activation
233 double valves[16];
234
235 /**
236 * @brief Update the atrial activation and LV/RV elastance functions which
237 * depend on time
238 *
239 * @param parameters Parameters of the model
240 */
241 void get_activation_and_elastance_functions(std::vector<double> &parameters);
242
243 /**
244 * @brief Compute sub-expressions that are part of atrial elastance and
245 * depends on atrial volume from the solution vector
246 *
247 * @param parameters Parameters of the model
248 * @param y Current solution
249 */
250 void get_psi_ra_la(std::vector<double> &parameters,
251 const Eigen::Matrix<double, Eigen::Dynamic, 1> &y);
252
253 /**
254 * @brief Valve positions for each heart chamber
255 *
256 * @param y Current solution
257 */
258 void get_valve_positions(const Eigen::Matrix<double, Eigen::Dynamic, 1> &y);
259};
260
261#endif // SVZERODSOLVER_MODEL_CLOSEDLOOPHEARTPULMONARY_HPP_
model::Block source file
BlockType
The types of blocks supported by the solver.
Definition BlockType.h:42
BlockClass
The classes/categories of blocks supported. Some classes require special handling (e....
Definition BlockType.h:64
SparseSystem source file.
Base class for 0D model components.
Definition Block.h:101
const int id
Global ID of the block.
Definition Block.h:103
const Model * model
The model to which the block belongs.
Definition Block.h:104
Heart and pulmonary circulation model.
Definition ClosedLoopHeartPulmonary.h:107
ClosedLoopHeartPulmonary(int id, Model *model)
Construct a new ClosedLoopHeartPulmonary object.
Definition ClosedLoopHeartPulmonary.h:115
void setup_dofs(DOFHandler &dofhandler)
Set up the degrees of freedom (DOF) of the block.
Definition ClosedLoopHeartPulmonary.cpp:35
void update_solution(SparseSystem &system, std::vector< double > &parameters, const Eigen::Matrix< double, Eigen::Dynamic, 1 > &y, const Eigen::Matrix< double, Eigen::Dynamic, 1 > &dy)
Update the solution-dependent contributions of the element in a sparse system.
Definition ClosedLoopHeartPulmonary.cpp:129
void post_solve(Eigen::Matrix< double, Eigen::Dynamic, 1 > &y)
Modify the solution after solving it.
Definition ClosedLoopHeartPulmonary.cpp:310
ParamId
Local IDs of the parameters.
Definition ClosedLoopHeartPulmonary.h:137
@ KXV_RA
Right atrium pressure-volume relationship (?)
Definition ClosedLoopHeartPulmonary.h:158
@ RRA_V
Right atrium outflow resistance.
Definition ClosedLoopHeartPulmonary.h:145
@ KXP_RA
Right atrium pressure-volume relationship (?)
Definition ClosedLoopHeartPulmonary.h:157
@ VRV_U
Right ventricle unstressed volume.
Definition ClosedLoopHeartPulmonary.h:152
@ LRV_A
Right ventricle inductance.
Definition ClosedLoopHeartPulmonary.h:146
@ RLV_AO
Left ventricle outflow resistance.
Definition ClosedLoopHeartPulmonary.h:151
@ EMAX_RA
Right atrium elastance (?)
Definition ClosedLoopHeartPulmonary.h:161
@ TPWAVE
Fraction of cardiac cycle (P-wave)
Definition ClosedLoopHeartPulmonary.h:139
@ IML
Scaling for intramyocardial pressure (left coronaries)
Definition ClosedLoopHeartPulmonary.h:142
@ IMR
Scaling for intramyocardial pressure (right coronaries)
Definition ClosedLoopHeartPulmonary.h:143
@ VASO_LA
Left atrium rest volume (?)
Definition ClosedLoopHeartPulmonary.h:164
@ LLV_A
Left ventricle inductance.
Definition ClosedLoopHeartPulmonary.h:150
@ ELV_S
Scaling for left ventricle elastance.
Definition ClosedLoopHeartPulmonary.h:141
@ CP
Pulmonary capacitance.
Definition ClosedLoopHeartPulmonary.h:155
@ VLV_U
Left ventricle unstressed volume.
Definition ClosedLoopHeartPulmonary.h:153
@ RLA_V
Left atrium outflow resistance.
Definition ClosedLoopHeartPulmonary.h:149
@ EMAX_LA
Left atrium elastance (?)
Definition ClosedLoopHeartPulmonary.h:162
@ TSA
Fractions of cardiac cycle (not sure)
Definition ClosedLoopHeartPulmonary.h:138
@ KXP_LA
Left atrium pressure-volume relationship (?)
Definition ClosedLoopHeartPulmonary.h:159
@ CPA
Aortic capacitance.
Definition ClosedLoopHeartPulmonary.h:156
@ RPD
Pulmonary resistance.
Definition ClosedLoopHeartPulmonary.h:154
@ LLA_V
Left atrium inductance.
Definition ClosedLoopHeartPulmonary.h:148
@ ERV_S
Scaling for right ventricle elastance.
Definition ClosedLoopHeartPulmonary.h:140
@ VASO_RA
Right atrium rest volume (?)
Definition ClosedLoopHeartPulmonary.h:163
@ KXV_LA
Left atrium pressure-volume relationship (?)
Definition ClosedLoopHeartPulmonary.h:160
@ RRV_A
Right ventricle outflow resistance.
Definition ClosedLoopHeartPulmonary.h:147
@ LRA_V
Right atrium inductance.
Definition ClosedLoopHeartPulmonary.h:144
void update_time(SparseSystem &system, std::vector< double > &parameters)
Update the time-dependent contributions of the element in a sparse system.
Definition ClosedLoopHeartPulmonary.cpp:106
TripletsContributions num_triplets
Number of triplets of element.
Definition ClosedLoopHeartPulmonary.h:223
void update_constant(SparseSystem &system, std::vector< double > &parameters)
Update the constant contributions of the element in a sparse system.
Definition ClosedLoopHeartPulmonary.cpp:41
Degree-of-freedom handler.
Definition DOFHandler.h:48
Model of 0D elements.
Definition Model.h:75
Sparse system.
Definition SparseSystem.h:57
Handles the properties of input parameters.
Definition Parameter.h:127
The number of triplets that the element contributes to the global system.
Definition Block.h:52