 Research
 Open Access
 Published:
Analytical method for calculation of deviations from intended dosages during multiinfusion
BioMedical Engineering OnLine volume 16, Article number: 18 (2017)
Abstract
Background
In this paper, a new method is presented that combines mechanical compliance effects with Poiseuille flow and pushout effects (“dead volume”) in one single mathematical framework for calculating dosing errors in multiinfusion setups. In contrast to existing numerical methods, our method produces explicit expressions that illustrate the mathematical dependencies of the dosing errors on hardware parameters and pump flow rate settings.
Methods
Our new approach uses the Ztransform to model the contents of the catheter, and after implementation in Mathematica (Wolfram), explicit expressions are produced automatically. Consistency of the resulting analytical expressions has been examined for limiting cases, and three types of invitro measurements have been performed to obtain a first experimental test of the validity of the theoretical results.
Results
The relative contribution of various factors affecting the dosing errors, such as the Poiseuille flow profile, resistance and internal volume of the catheter, mechanical compliance of the syringes and the various pump flow rate settings, can now be discerned clearly in the structure of the expressions generated by our method. The invitro experiments showed a standard deviation between theory and experiment of 14% for the delay time in the catheter, and of 13% for the time duration of the dosing error bolus.
Conclusions
Our method provides insight and predictability in a large range of possible situations involving many variables and dependencies, which is potentially very useful for e.g. the development of a fast, bedside tool (“calculator”) that provides the clinician with a precise prediction of dosing errors and delay times interactively for many scenario’s. The interactive nature of such a device has now been made feasible by the fact that, using our method, explicit expressions are available for these situations, as opposed to conventional timeconsuming numerical simulations.
Background
If a patient needs intravenous administration of medication, it is important for the clinician to understand the basic pharmacokinetics of these medications. However, in recent years, ample evidence has been found [1,2,3,4], including a number of clinical cases and in vivo studies [5,6,7,8,9], indicating that physical effects related to the infusion hardware are equally important to understand [10, 11]. Especially in critical care, on the ICU and the OR, where multiple medications are typically delivered through one single lumen of a thin catheter, these physical effects may cause ambiguous and counterintuitive discrepancies between the intended dose and the dose that has been actually delivered. This is particularly true if the dosing rate is adjusted ad hoc [12]; for example, fastacting and critical inotropic drugs are often titrated, based on the mean arterial blood pressure (MABP).
There are three major factors that can produce significant deviations from the intended medication dosing rate scheme: (i) the length of the catheter causes a delay in the administration of the medication into the patient, and therefore, a mixture, corresponding to previous medication dosing rates, may still be present inside the catheter, i.e. the contents inside the catheter constitutes a “memory” in which the effects of previous medication dosing rates may be stored. This has been called “dead volume effect” in the literature [7, 13]. (ii) The syringes that are used in clinical practice are far from ideal, because these syringes have a significant mechanical compliance primarily due to the compressibility of the rubber plunger inside the syringe [14]. Therefore, whenever the flow rate setting of one of the infusion pumps in a multiinfusion setup is altered, pressure changes within the entire system cause a change in the deformation of compressible and expandable parts (i.e., the other syringes) within the system, and hence may cause a deviation from the intended flow rates. (iii) The low velocities of the fluids in catheters ensure that the flow inside these lines will be laminar (low Reynolds number), and hence exhibit a Poiseuille flow profile, in which the fluid particles near the central longitudinal axis of the catheter travel faster than fluid particles near the wall of the catheter, thus giving rise to a “mixing effect” of its own.
Several infusion simulation studies have shown the importance of these problems by demonstrating the influence of the physical effects of dead volume, compliance and the Poiseuille effect on drug delivery in isolation: There are a number of studies that describe the “dead volume effect” in isolation from the syringe compressibility (“compliance”) effect [7]. In simple cases, calculation of the “dead volume effect” is straightforward: If the actual volumetric flow rate \(u_{cath}\) inside a single lumen catheter is assumed to be constant, and the internal volume \(V_{cath}\) of the catheter (i.e., the internal volume as measured starting from the mixing point, at which the medications from all syringes come together, up to the distal catheter tip inside the vasculature of a patient) is known, then calculation of the “delay time” \(t_{delay}\), i.e. the time needed for a droplet of medication to travel through the dead volume before it reaches the blood stream of the patient, is simply \(t_{delay} = V_{cath}/u_{cath}\).
In clinical practice, however, the situation is often more complex: the flow rate may vary during the delay time due to changes in the pump flow rate settings, and the actual flow rate also typically differs from the pump flow rate setting value temporarily, due to the effects of mechanical compliance as mentioned above. These compliance effects have been studied in isolation as well [2, 3, 15]. Its basic mechanism is a “capacitor” effect, which has been quantified using the electric analog of a system containing capacitors, and the calculations of the dosing errors have been performed using the Laplace transform.
The third major factor, the mixing effect due to the Poiseuille profile [16], can be described as a convolution, as will be explained in more detail in “Appendix”. It causes a spreading out of the dosing error in time, in which the first arrival of the dosing error occurs sooner then it would have without Poiseuille mixing effect.
As a result, due to the combination of the effects (i), (ii), and (iii) mentioned above, nontrivial deviations from the intended medication dosages scheme may occur, as will be explained below.
In an earlier paper, we have shown that the “dead volume effect” on one hand, and syringe compliance (“capacitor”) effect on the other hand, produce opposite deviations from the pump flow rate settings in the actual drug output concentrations, making the net result hard to predict and often counterintuitive [4].
In this paper, however, we focus on the mathematical method of calculating dosing errors in multiinfusion setups in complex situations. The aim of the new method described in this paper is to obtain analytical expressions for the deviations from the intended dosages during multiinfusion that result from the combination of the effects (i) and (ii) described above. These analytical expressions contain dependencies on parameters that represent the physical variables in a multiinfusion setup, such as: intended flow rates of the various pumps, compliances of the syringes, resistances of the various tubes, height differences within the multiinfusion setup, etc.
As opposed to conventional numerical simulations, our objective therefore is to obtain explicit analytical expressions for these deviations, in which physical characteristics of the multiinfusion setup are represented explicitly as variables.
We see the need for a tool for clinicians and medical physicist that provides understanding of the role that various physical parameters (i.e., characteristics of the multiinfusion hardware) play in the emergence of a dosing error. In order to make these roles explicit, we will use our new mathematical approach, presented in the “Method: analytical model, and invitro setup” section of this paper, to explicate these roles in the form of direct mathematical relationships between physical hardware parameters and dosing errors.
For our mathematical model we used the strong symbolic calculation capabilities of the Mathematica package (Mathematica 10, Wolfram^{®} Inc., USA) to process Laplacetransforms and Ztransforms analytically. This will be explained in “Method: analytical model, and invitro setup” section. Furthermore, we have performed invitro measurements in order to verify our mathematical results experimentally.
Method: analytical model, and invitro setup
To start with, we analyse the total chain of causality, from a change in the pump flow rate setting value, up to the moment that a dosing error enters the bloodstream of the patient. As has been mentioned in the introduction, it is essential to recognize that the lumen of the catheter (i.e. the internal volume \(\mathcal {C}\) inside the catheter, starting from the mixing point \(\mathcal {M}\), at which the medications from all syringes come together, up to the distal tip \(\mathcal {P}\) of the infusion line inside the vasculature of a patient) constitutes a memory in which the effects of previous changes in pump flow rate settings are stored. In the explanation of our method, we focus on a singlelumen catheter. This however can be easily extended to a multilumen catheter. In order to represent this memory in the mathematical model, the internal volume of the catheter \(\mathcal {C}\) between the points \(\mathcal {M}\) and \(\mathcal {P}\) is divided into tiny voxels \(A_k\), in which the index k runs from \(k=0\) at point \(\mathcal {P}\) to \(k=N\) at point \(\mathcal {M}\) (see Fig. 1), in which N is a very large number. Let L denote the length of the catheter. Hence, the length of a single voxel equals \(\gamma \), with \(\gamma = L/N\). In our mathematical model, we use a general setup containing any number of infusion pumps. For simplicity, however, we start with three pumps, in which each of these pumps contains a solution of a different medication, in which we used colours (“Red”, “Green”, and “Blue”, or R, G, and B, respectively) to denote the three different solutions. It is important to note that the R, G, and B denote the three different solutions as stored in their respective syringes, not the medications themselves. Furthermore, in this paper, the dosing errors will be expressed as volumes of the R, G, and B solutions, instead of dosages of the medications themselves. In Fig. 1a, the situation directly after \(t=0\) is depicted, i.e. the point in time at which a change in pump flow rate setting value of one of the pumps takes place. As a result, the first small voxel at k=N near the mixing point \(\mathcal {M}\) is now being filled with a droplet \(\xi \) featuring a new mixing ratio between the solutions R, G, and B, resulting from the new pump flow rate setting values at \(t=0\). The rest of the voxels inside the line \(\mathcal {C}\), however, contain a fluid mixture that still has the old mixing ratio corresponding to the steadystate situation before the change in pump flow rate setting value. In Fig. 1b, the situation at \(t=t_{delay}\) is depicted, i.e. the point in time at which this specific droplet \(\xi \), that has been situated inside voxel \(A_{k=N}\) at \(t=0\), has now reached the distal tip of the catheter at point \(\mathcal {P}\) at \(t=t_{delay}\). Hence, at \(t=t_{delay}\), the entire line \(\mathcal {C}\) has now been filled with the new mixing ratio. The contents of voxel \(A_k\), for any value of k, is described by the 3Dvector \(\bar{a}_k\), in which
in which
The ratio \(a_{k}^{(R)}\) represents the volume fraction of fluid (R) inside voxel k, i.e. the volume of solution (R) inside voxel k divided by the total volume inside voxel k. The voxel for which \(k=0\) is the very last voxel at point \(\mathcal {P}\) inside the catheter, i.e. the distal tip of the catheter. This specific voxel at \(k=0\) releases its contents \(\bar{a}_{k=0}\) directly into the blood stream of the patient.
Now consider the situation in which the flow rate setting of the “Green” pump is changed suddenly at \(t=0\), whereas the flow rate settings of the “Red” and “Blue” pumps are never altered. If the “old” total flow rate inside the catheter \(\mathcal {C}\) before \(t=0\) has been constant and equal to \(u_{old}\) for all \(t<0\), and the “new” total flow rate is also constant and equals \(u_{new}\) for all \(t>0\), then the partial flow rate \(u^{(R)}_{patient}\) of the “Red” fluid leaving the catheter at point \(\mathcal {P}\) (and thus entering the patient) at \(t=0\) equals:
in which the “old” mixing ratio refers to the mixing ratio produced by the pump rate settings of \(t<0\).
As is explained in Fig. 2, the presence of the “old” mixing ratio in Eq. (3) can give rise to effects that may be unexpected to physicians, in which a temporary increase in the partial flow rate of (R) entering the blood stream occurs, although the setting value of the flow rate of the pump of (R) has never been changed. We will refer to this unwanted temporary increase in the \(u^{(R)}_{patient}\) as the “pushout” effect, which has sometimes been called “dead volume” effect in the literature [7]. In Fig. 2a, the situation before \(t=0\) is depicted. In Fig. 2b, a temporary, unplanned, and undesirable increase in \( u^{(R)}_{patient}\) occurs (see figure legend for explanation). After \(t=t_{delay}\), the value of \( u^{(R)}_{patient}\) returns to its original value (Fig. 2c).
Furthermore, due to the mechanical compliance of the syringes, however (as mentioned at point (ii) in the introduction), the actual flow rate \(u_{cath}(t)\) in the catheter does not follow the changes in pump flow rate setting values immediately. As a result, the actual flow rate \(u_{cath}(t)\) is not exactly equal to the sum of the pump flow rate setting values \(u_{pump}^{(R)}(t)\), \(u_{pump}^{(G)}(t)\), and \(u_{pump}^{(B)}(t)\), due to the transient mechanical compliance effects that cause extra fluid to be stored in, or released from, the syringes directly after changes in pump flow rate setting values:
Therefore, the general expression from which the “delay time” \(t_{delay}\) needs to be established is:
in which \(u_{cath}\) depends on the compliance characteristics of the syringes as well.
New method for incorporating the memory effect of the catheter into an analytical model
In this paper, we will follow the following strategy to calculate dosing errors analytically (see Fig. 3):
First, standard techniques [17] (involving the Laplace Transform) are used to calculate pressures and total flow rates (see Fig. 3a), in which the differences between the “Red”, “Blue” and “Green” solutions are disregarded; i.e. in this “colorblind” calculation, only the total flow rate inside the catheter \(u_{cath}(t)\) is calculated. The mechanical compliance of the catheter is very small with respect to the mechanical compliance of the syringes [18], and therefore the mechanical compliance of the catheter is neglected in our model. As a result, the flow rate entering the tube at time t and the flow rate leaving the tube (entering the blood stream) at the same time t are both equal to \(u_{cath}(t)\), disregarding the different constitutions that the fluids entering, and leaving, respectively, may have in terms of the partial fluids R, G, and B.
Secondly, our new analytical method is introduced (see Fig. 3b), which enables incorporation of the “memory effect” of the catheter into the model. In order to make the analytical approach possible throughout the entire calculation up to the end of point (d) in the figure, a formulation in the Zdomain (using the Ztransform, which is a discrete variant of the Laplace transform) is introduced. This analytical method uses the results form the standard (Laplacebased) method [see (a)] as an input. This input has the form of general expressions for the various total flow rates as function of time. The Ztransform formulation in our method [see (b)] therefore does not replace the Laplacebased method from (a), but is used after it.
Third, the Zdomain formulation from (b) enables an easy, analytical, incorporation of the Poiseuille mixing effect (c) into the model, because the convolutionlike nature of the Poiseuille mixing effect is reduced to a mere multiplication in the Zdomain. The mathematical details of the Poiseuille mixing effect are described in “Poiseuille mixing effect” section.
Finally (d), we derive expressions for the volume, and the first and second moment, of the dosing error as function of hardware parameters. These moments can be calculated directly within the Zdomain. From the first and second moment, key characteristics concerning timing and duration of the dosing error are derived.
Throughout the entire process, we have used the symbolic calculation capabilities of the Mathematica package (Mathematica 10, Wolfram^{®} Inc., USA) to process Laplacetransforms and Ztransforms analytically, as can be seen in “Results” section.
Starting point: Laplace transform and Kirchhoff’s laws
In this paper, we focus on the multiinfusion of aqueous solutions, in which all solutions have a viscosity close to the viscosity of water. Therefore, in our mathematical model, all fluids are assumed to have approximately the same viscosity. For the calculation of the pressures and total flow rates, we use an electric circuit model that is analogous to the infusion setup, in which current sources (the infusion pumps), resistances (infusion lines and catheter), and capacitors (mechanical compliance of the syringes) are present, and in which Kirchhoff’s laws are applied on the voltages (pressures) and currents (flow rates) in the Laplace domain. The output \(U_{cath}(s)\) of the calculation is the Laplace transform of the total flow rate \(u_{cath}(t)\) inside the catheter. This \(u_{cath}(t)\) equals the total flow rate entering the tube at time t as well as the total flow rate leaving the tube (entering the blood stream) at the same time t. Let \(\Gamma \) denote the complete set of these changes in pump flow rate setting values, together with the physical hardware parameters of the setup, which are present in the form of explicit parameters in the analytical expressions for \(U_{cath}(s)\). Using the Mathematica package (Mathematica 10, Wolfram^{®} Inc, USA), we were able to perform an inverse Laplace transform in order to retrieve an analytical expression for \(u_{cath}(t, \Gamma )\) from the \(U_{cath}(s, \Gamma )\), in which the above mentioned hardware parameters \(\Gamma \) are present as explicit variables. This result is rendered in Eq. (33) in “Results” section.
Key parameters of dosing errors
The dosing errors that are produced after changes in pump flow rate setting values in a multiinfusion setup are of a temporary nature, i.e. “dead volume” and “mechanical compliance” (points (i) and (ii), respectively, as mentioned in the introduction) give rise to only a temporary deviation of the actual dose rates (entering the blood stream) from the intended ones (see e.g. Fig. 2). As a result, each of these dosing errors has the shape of a “bolus”.
Let \(\beta ^{(R)}_{patient}(t)\) denote such a bolusshapes dosing error, i.e. let \(\beta _{patient}^{(R)}(t)\) denote the deviation from the intended partial flow rate [in this case of solution (R)], as a function of time, in which
Defining \(t=0\) as the moment that a change in pump flow rate setting value takes place, and defining \(\beta _{patient}^{(R)}(t)\) as the dosing error due to (only) this specific change in pump flow rate setting value at \(t=0\), then we have by definition:
in agreement with the bolusshaped nature of the dosing error \(\beta _{patient}^{(R)}(t)\).
There are two basic mechanisms that contribute to the dosing errors, and the effect of each of these two mechanisms is restricted to a particular time interval: During the time interval \(0< t < t_{delay}\), the dosing error entering the blood stream is caused by the pushout effect, as explained in Fig. 2. During the time interval \(t > t_{delay}\), however, the dosing error entering the blood stream is caused by the effects of mechanical compliance, because the “pushout effect” is by definition restricted to the time interval \(0< t < t_{delay}\). The effects of mechanical compliance are particularly strong directly after changing the pump flow rate settings, and thus the effects of mechanical compliance (in the form of deviating mixtures) are entering the catheter at point \(\mathcal {M}\) at \(t=0\), and hence entering the patient at point \(\mathcal {P}\) directly after \(t = t_{delay}\).
In the following, we will use \(t^{POIS}_{delay}\) instead of \(t_{delay}\), because, as will be demonstrated below, the Poiseuille profile of the flow affects the \(t_{delay}\), resulting in a modified delay time that we will refer to as \(t^{POIS}_{delay}\). Furthermore, we will consider a situation with two pumps (“Red” and “Green”) only, and will change the flow rate setting of the “Green” pump only (at \(t=0\)), leaving the flow rate setting of the “Red” pump (denoted as \(u_{pump}^{(R)}\)) constant. This does not affect the applicability of our method in more general cases. The flow rate setting of the “Green” pump changes from \(u_{old}^{(G)}\) to \(u_{final}^{(G)}\) at \(t=0\).
The \(\beta _{patient}^{(R)}(t)\) for \(0< t < t^{POIS}_{delay}\) is due to the pushout effect, and hence its calculation is relatively easy. The \(\beta _{patient}^{(R)}(t)\) for \(t > t^{POIS}_{delay}\) however is due to the mechanical compliance, and hence its calculation is much more complex. Therefore, we split the \(\beta _{patient}^{(R)}(t)\) into two parts according to the two corresponding distinct time intervals:
in which \(\tau = t  t^{POIS}_{delay}\) and \(\psi _{patient}^{(R)}(\tau )\) represents the effects of mechanical compliance for \(t > t^{POIS}_{delay}\), and in which the dimensionless \(\psi _{patient}^{(R)}(\tau )\) is bolusshaped of its own:
In the following, we concentrate on finding analytical expressions for (the moments of) \(\psi _{patient}^{(R)}(\tau )\). To this end, the discretized version \(\psi _k^{(R)patient}\) of \(\psi _{patient}^{(R)}(\tau )\) is defined as:
Contents of the catheter without poiseuille mixing effect
In this subsection, we derive a generic expression for the ZTransform of the contents \(\{a_k^{(R)}\}\) of the catheter, without incorporating the Poiseuille mixing effect. The Poiseuille mixing effect will be incorporated in the model later on.
First, let \(\lambda _{tot}(t)\) denote the distance that the fluid inside the catheter has traveled as a result of the total flow rate \(u_{cath}\), i.e.:
Furthermore, let \(u_{\mathcal {M}}^{(R)}\) denote the actual partial flow rate of the Red fluid entering the catheter at the mixing point \(\mathcal {M}\).
For some voxel k inside the catheter, the value of the ratio \(a_k^{(R)}\) reflects the \(u_{\mathcal {M}}^{(R)}(t^{*})\) divided by the total flow rate \(u_{cath}(t^{*})\), at the time \(t^{*}\) that the contents of that particular voxel entered the catheter at \(\mathcal {M}\), i.e.:
in which \(\gamma = L/N\). In order to facilitate a Ztransform, we rewrite Eq. (12) in the following form:
The value of the Dirac delta \(\delta (\gamma k  \lambda _{tot}(t))\) is zero everywhere, except around \(t = t^{*}\); hence evaluation of the integral in Eq. (13) shows the equivalence of Eqs. (12) and (13):
in which \(\varepsilon \) and \(\eta \) are both very close to zero.
Using the displacement rule from Ztransform theory [17], Eq. (13) now yields the following expression A(z) in the Zdomain:
Here it needs to be mentioned that, since the Ztransform is merely a discrete form of the Laplace transform, it would be, theoretically, possible to replace this Ztransform by yet another Laplace transform. However, for clarity, we have chosen to use the Ztransform, because it connects to the “voxelbased” description of the contents inside the catheter in an intuitive way, and because it creates a clear distinction from the standard Laplace method used to calculate the pressures and flows that enter the mixing point \(\mathcal {M}\).
Similarly, we may define the basic input for the calculation of the dosing error of fluid (R) as the difference between \( u_{\mathcal {M}}^{(R)}\) (i.e. the actual partial flow rate of the Red fluid entering the catheter at the mixing point \(\mathcal {M}\)), and the intended \(u_{pump}^{(R)}\) (i.e. the intended partial flow rate of the Red fluid entering the catheter at the mixing point \(\mathcal {M}\)):
yielding
This approach now enables easy incorporation of the Poiseuille mixing effect into our model using a simple multiplication in the Zdomain, as is explained below.
Poiseuille mixing effect
The basic input for calculation of the Poiseuille mixing effect is the input from the Red line at \(\mathcal {M}\) in the form of
which corresponds to \(A^{diff}(z)\) in the zdomain [see Eq. (17)]. The mixing effect due to the Poiseuille profile causes a spreading out of the dosing error in time, in which the first arrival of the dosing error occurs sooner then it would without Poiseuille mixing effect. In a Poiseuile flow, the velocity along the centerline of the tube equals two times the average velocity. This implies that, after a change in pump flow rate settings at \(t=0\), the tip of the parabolic flow profile reaches point \(\mathcal {P}\) at the end of the catheter already at \(t = t_{delay}/2\). Hence, we define \(t^{POIS}_{delay}\) as: \(t^{POIS}_{delay} = t_{delay}/2\). Because the Poiseuille mixing effect is a phenomenon that is created within the catheter during the very passage of the liquid through the catheter from point \(\mathcal {M}\) to point \(\mathcal {P}\), an accurate description of the full extent of this effect can only be given in the form of the contents of the very last voxel at point \(\mathcal {P}\), just before this contents is released into bloodstream of the patient. This is illustrated in Fig. 4, in which the top diagram (Fig. 4a) shows the Poiseuille profile along the length of the catheter from point \(\mathcal {M}\) to point \(\mathcal {P}\), in which the tip of the red, innermost, parabola has just finished its journey through the catheter, and is just to be released into the blood stream. Each parabola in Fig. 4a constitutes a boundary between two volumes of liquid that once were subsequent, undistorted, “regular” voxels when they entered the catheter at \(\mathcal {M}\). These initially flat boundaries between voxels near \(\mathcal {M}\) are being stretched and distorted into the parabolashaped boundaries during their travel to \(\mathcal {P}\).
In “Appendix”, it is derived that the Poiseuille flow profile causes a mixing effect that features a simple linear relation between concentration and distance along the line from \(\mathcal {M}\) to \(\mathcal {P}\) (see Fig. 4b), which is consistent with earlier work by e.g. Taylor [19] and Hutton and Thornberry [20].
Using this simple linear relationship, the constitution of the final voxel of the catheter at \(\mathcal {P}\) can be defined in terms of the original, undistorted, voxels that once entered the catheter at \(\mathcal {M}\), as is illustrated in Fig. 4c. For instance, the contents of the “voxel \(*\)”, indicated by the symbol \(*\), which is voxel nr i counting from point \(\mathcal {P}\), and which is limited by two downward sloping lines and the thick blue horizontal line segment in Fig. 4c, contributes to the final voxel at \(\mathcal {P}\) with weight factor \(w_i\) (indicated by the small vertical thick blue line segment). The dimensionless weight factor \(w_i\) equals the voxel length \(\gamma \) times the first derivative of the expression \({{x}\over {L + x}}\), which is depicted in Fig. 4b, and which is established on the basis of geometrical reasoning using the triangles involved:
which yields
Furthermore, the contents of each “voxel” i, with i running from 0 to infinity, equals that of voxel \(b^{(R)}_{ki}\), in which \(b^{(R)}\) refers to the green curve in Fig. 4d, which is a “stretchedout” version of the original voxels \(a^{(R)diff}_j\) that entered the catheter at \(\mathcal {M}\), i.e.:
which reflects the fact that the velocity along the centerline of the tube equals two times the average velocity. As a result, the summation of all contributions \(b^{(R)}_{ki}\) with weight factors \(w_i\) yields a (discrete) convolution, which is visible in the first term of the following equation, in which the fraction of the Red solution inside the last voxel at point \(\mathcal {P}\) is denoted by the dimensionless scalar \(\psi _k^{(R)patient}\), and in which k refers to the position of the beginning of the green curve b (which shifts in time):
The second term in Eq. (22) does not depend on the compliances and resistances in the multiinfusion setup, but represents the remnant fluid from the old mixture due to the poiseuille flow profile.
When calculating the first and second moment of \(\psi _k^{(R)patient}\) (or \(t_{central}\) and the variance \(\sigma ^2\) from Fig. 5, respectively), these moments are simply the sum of the moments of the two distinct terms in Eq. (22). The moments of the second term can be calculated directly in the spatial domain (using an upper limit to the amount of fluid). For the first term, however, the Ztransform is needed, as will be explained below. Therefore, we split the \( \psi _k^{(R)patient}\) according to:
in which \( \psi _k^{(R)POIS} = \sum _{j=0}^k b_{kj}^{(R)} w_{j} \; \text{ and } \; \psi _k^{(R)remnant} = {{L}\over {(L + \gamma k)}}\frac{ u_{pump}^{(R)} (u_{final}^{(G)}u_ {old}^{(G)})}{ \left( u_{final}^{(G)}+u_{pump}^{(R)}\right) \left( u_ {old}^{(G)}+u_{pump}^{(R)}\right) } \) and from now on will concentrate on the calculation of the first term, i.e. \(\psi _k^{(R)POIS}\), only.
Performing the Ztransform on \(\psi _k^{(R)POIS}\) yields:
in which B(z) is the Ztransform of \(b_j\), and W(z) is the Ztransform of \(w_j\), which yields:
Since \(L/{\gamma } = N\) is a large number, this constitutes a rapidly converging series for \(z \approx 1\). Furthermore, since \(b_{2j}^{(R)} = a_j^{(R)diff}\) (Eq. 21), we have \(B(z)=A^{diff}(z^2)\), in which \(A^{diff}(z)\) is the Ztransform of the original \(a_k^{(R)diff}\), in the form in which it entered the catheter at \(\mathcal {M}\), before being distorted due to the Poiseuille mixing effect [see Eq. (17)].
Therefore, in summary, we now have
in which \(\Psi ^{(R)POIS}(z)\) is the Ztransform of \(\psi _k^{(R)POIS}\). This expression for \(\Psi ^{(R)POIS}(z)\), i.e. Eq. (26), constitutes the central equation in our method, and enables derivation of explicit expressions for the first and second moment of \(\psi _k^{(R)POIS}\) (or \(t_{central}\) and \(\sigma \) from Fig. 5, respectively), without the need of performing an inverse ZTransform.
Calculation of moments of the dosing error distribution after T\(_{delay}\)
We will now use the “theorems of moments” from Ztransform theory to provide expressions for key characteristics of the \(\psi _{patient}^{(R)}(\tau )\) that enters the patient.
The zero’th moment yields:
and hence the total volume of the dosing error \(\psi _{patient}^{(R)}(\tau )\) equals:
The first moment yields:
and hence
Furthermore, using the second moment, the following general expression for \(\sigma \) is obtained, in which the typical “duration” or “width” of a dosing error is characterized as \(2\sigma \) (see Fig. 5):
In “Results” section, Eqs. (27)–(31) will be used to produce analytical expressions for various situations. Furthermore, also in “Results” section, Eqs. (27)–(31) will be evaluated for situations without the Poiseuille mixing effect, in which case the \(\Psi ^{(R)POIS}(z)\) in these equations is replaced by merely A(z).
Spectrometric invitro setup for feasibility tests
A schematic of the invitro setup is rendered in Fig. 6. Two Perfusor B. Braun^{®} syringe pumps (B. Braun, Melsungen AG, Germany) and 50ml syringes (Melsungen AG, Germany) were used in the experimental setup. The syringes contained Tartrazine (TT) and Indigo Carmine (IC) solution, solved in distilled water. Each syringe was connected to an infusion line (d = 1 mm) of 200 cm (Cair LGL^{®}, France) and subsequently combined using the using a 3needlefree Yconnector (20038E7D, Cardinal Health^{®}, Switzerland). The Yconnector outflow was connected to a flowcell (Z flowcell w/SMA 905, 10mm pathlength, FIAlab^{®}, Seattle, WA, USA), this flowcell was also connected to a visual light spectrum (250–2500 nm) DT1000 light source (Ocean Optics^{®}, Dunedin, FL, USA). This allowed a spectrometer QE65000 (Ocean Optics^{®}, Dunedin, FL, USA), also connected to the same flowcell, to continuously measure an absorption spectrum of the dye mixture. Concentrations of each dye were acquired from this absorption spectrum [4]. After the flow cell, the cumulative flow rate was measured using three M12p flowmeters (Bronkhorst^{®}, Ruurlo, The Netherlands). Sample time of all the measurements was 1 s.
Results
General result from laplace transform as a starting point
Performing the first step in the scheme rendered in Fig. 3, a general expression for \(u_{\mathcal {M}}^{(R)diff}\) has been derived in the Laplace domain, using a known technique, i.e. Laplace Transform in combination with Kirchhoff laws, on the basis of the electric analog depicted in Fig. 7:
in which \(\breve{u}_{\mathcal {M}}^{(R)diff}(s)\) is the Laplace transform of \(u_{\mathcal {M}}^{(R)diff}(t)\), and \(u_{downstep}^{(G)}\) equals the change in pump flow rate setting of the Green pump at \(t=0\), and \(C_1\) and \(C_2\) are the mechanical compliances of the “green” and “red” syringe, respectively, and \(R_1\) and \(R_2\) the resistances of the “green” and “red” feeding lines, respectively. Meanwhile, the pump flow rate setting of the Red pump, \(u_{pump}^{(R)}\) remains unchanged all the time. Using the inverse Laplace transform function of Mathematica, the following general expression was obtained in the time domain:
in which
Using Eq. (32), the total volume of the dosing error [i.e., the timeintegral of Eq. (33)] is easily calculated using the general rules for integration and for calculation of limits from Laplace transform theory. The resulting expression \( Q^{(R)}_{dosingerror}(\Gamma )\), i.e. the volume of the dosing error, depends on the set of hardware parameters \(\Gamma \). In this specific case, the parameter set \(\Gamma \) comprises of: \(\Gamma = \{Resistances, Compliances, Flow Rates\}\).
The subtraction of two exponentially decaying functions in Eq. (33) yields a characteristic bulbshaped graph of \(u_{\mathcal {M}}^{(R)diff}(\tau )\) as function of the time \(\tau \); see Fig. 8. In this figure, Eq. (33) has been evaluated for a number of parameter settings of \(\Gamma \).
From the general equation, more specific and practical expressions can be derived by substituting standard values into the equation, leaving the parameters of interest as variables. The standard values for \(R_{cath}\), \(R_1\), \(R_2\), \(C_1\), \(C_2\), \(u_{pump}^{(R)}\), \(u_{old}^{(G)}\) and \(u^{(G)}_{downstep}\) are rendered in Table 1.
Typical example: resulting dosing error during syringe exchange
We now consider the specific case in which the “Green” syringe is exchanged, during which the green pump is stopped and green line has been clamped (with a Kocher), but the red line is not clamped, and the settings of the red pump remain unchanged. Let \(T_{restart}\) denote the time between clamping the green line and the reopening the green line, i.e. \(T_{restart}\) denotes the time duration of the entire procedure of exchanging the green syringe.
During the syringe exchange time interval, the red line still oozes red fluid in undiluted form into the first voxels of the catheter directly after the mixing point \(\mathcal {M}\). We now consider two scenario’s within the general scenario of the syringe exchange with the green line being clamped off during the entire exchange time interval. In the first scenario, the red line is not clamped off, but the red pump is switched off (pump flow rate setting value is zero; still connected to the red line) simultaneously at the beginning of the syringe exchange procedure, and restarted to its original pump flow rate setting value at the very end of the syringe exchange procedure. In the second scenario, no actions are performed in relation to the red pump or red line at all. The first scenario entails that accumulation of undiluted red fluid directly beyond the mixing point is caused only by the compliance of the red syringe, in combination with the drop in pressure caused by the clamping of the green line. This is indicated by the text “compliance” under the second brace in Eq. 36. In the second scenario, however, the regular pumping action of the red pump (at flow rate \(u_{pump}^{(R)}\)) adds to the accumulation of undiluted red fluid beyond the mixing point as well (indicated by the text “red pump on” under the first underbrace in Eq. 36). As a result, in the second scenario, the total dosing error \(Q^{(R)}_{dosingerror}(T_{exchange})\) is the sum of the “compliance” and the “red pump on” effects. See Fig. 9. After \(t_{delay}\), this total dosing error is entering the blood stream of the patient within a very short time interval, due to the fact that the green pump has a high flow rate.
In Eq. (36), the second term (“compliance”) has been calculated on the basis of integrating the right side of Eq. (33) over time.
Analytical results for reducing the flow rate of the green pump, without Poiseuille mixing effect
We now consider the more general case in which the pump flow rate setting of green pump, i.e. of the fast pump, is lowered with an amount \(u^{(G)}_{downstep}\), without clamping any line. The volume \(Q^{(R)}_{dosingerror}\) does not depend on the specific distribution of the dosing error over time, and hence yields the same result as in Eq. (35).
Here we present the general result of the calculation of the first moment of A(z) without the Poiseuille flow effect:
Applying this result to the specific case described in Eq. (32), i.e. lowering the pump flow rate setting of the green pump, and substituting Eq. (33) into Eq. (37), yields the following result:
In many clinical situations we have \(C_1 = C_2 = C\) and \(R_1 = R_2 = R\), in which case Eq. (38) reduces to:
Applying Eq. (31) from the “Method: analytical model, and invitro setup” section to the present situation, without the Poiseuille flow effect reads:
Evaluation of Eq. (40) yields a very long and unwieldy expression. However, as has been noted before, in many clinical situations we have \(C_1 = C_2 = C\) and \(R_1 = R_2 = R\), and, furthermore, for most catheters we have \(R \ll R_{cath}\). Applying these assumptions, we obtain the following result:
This result will be compared with the findings from the invitro experiments in “Results of invitro experiments, and comparison with theoretical predictions” section.
If, in Eq. (41), the \(u^{(G)}_{downstep}\) would be very small with respect to \(u_{final}\), then this expression for \(\sigma \) approaches the familiar \(\sigma \approx 2 C R_{cath}\).
Analytical results for reducing the flow rate of the green pump, incorporating the Poiseuille mixing effect
Evaluation of Eq. (31) from “Method: analytical model, and invitro setup” section to the situation of reducing the pump flow rate setting of the green pump (still without any clamping) once more, but now substituting \(\Psi ^{(R)POIS}(z)\) [instead of A(z)] in order to incorporate the Poiseuille flow effect, yields a very long and unwieldy expression. Therefore, again, we apply the assumptions \(C_1 = C_2 = C\) and \(R_1 = R_2 = R\) and \(R \ll R_{cath}\), which yields the following result:
As can be seen in Eq. (42), the strength of the contribution of the Poiseuille mixing effect to the width \(\sigma ^{POIS}\) depends on the internal volume of the catheter \(V_{cath}\) with respect to the time \(C R_{cath}\) in combination with the flow rates \((u^{(G)}_{downstep})\) and \(u_{final}\); or, more precisely, on the ratio between \(3 (V_{cath})^2\) and \(C^2 R_{cath}^2 ((u^{(G)}_{downstep})^2+12 u^{(G)}_{downstep} u_{final}+48 (u_{final})^2 )\).
If the catheter has a large internal volume but the time \(C R_{cath}\) is short and the flow rates are low, then the \(\sigma ^{POIS}\) in Eq. (42) reduces to: \(\sigma ^{POIS} \approx V_{cath} / u_{final}\).
If, however, in Eq. (42), the \(u^{(G)}_{downstep}\) would be very small with respect to \(u_{final}\), and the time \(C R_{cath}\) would be very large and the \(V_{cath}\) would be small, then this expression for \(\sigma \) approaches \(\sigma ^{POIS} \approx 4 C R_{cath}\), which is two times larger than the \(\sigma \approx 2 C R_{cath}\) that was calculated before when omitting the Poiseuille mixing effect. This factor two is a result from the fact that, in a Poiseuille flow, the velocity at the centerline of the catheter equals two times the average velocity.
Results of invitro experiments, and comparison with theoretical predictions
Three types of invitro experiments have been performed, in order to compare theoretical results with actual measurements. These three types of experiments are:

(i)
measurement of the flow rate of the Green fluid as function of time, immediately after a change in pump flow rate setting value of the green pump, corresponding to \(u^{(G)}_{downstep} = \) 6 ml/h; see Fig. 10.

(ii)
measurement of \(t_{delay}\) for a set of three different syringes and three resistances, corresponding to various values of C and \(R_{cath}\); see Fig. 11. The experimental measurement of the \(t_{delay}\) showed a standard deviation between theory and experiment of 16.2 s, which is 14% of the total range of variations in \(t_{delay}\).

(iii)
measurement of \(\Delta t_{central}\) as function of various values of the relative resistance of the catheter with respect to its standard value R0, in which \(\Delta t_{central}\) was defined as the measured \(t_{central}\) minus the standard value of \(t_{central}\) for \(R_{cath}=R0\). Three different values of \(R_{cath}/R_0\) have been used; and the measurement \(\Delta t_{central}\) of was repeated three times for each value of \(R_{cath}/R_0\). See Fig. 12. The experiment measuring the \(\Delta t_{central}\) showed a reasonable agreement between theory and measurements, in which standard deviation between theory and experiment was 8.4 s, which is 13% of the total range of variations in \(\Delta t_{central}\).
Discussion
In this paper, explicit expressions have been derived for the total volume (Q), central time point (\(t_{central}\)) and “width” or “duration” (\(2\sigma \)) of dosing errors as function of hardware parameters such as mechanical compliance of syringes, resistance of the catheter, length of the catheter, for two general cases of a change in pump flow rate setting value in a multiinfusion setup, as well as for a typical case of syringe changeover. Consistency of the resulting analytical expressions has been examined for limiting cases, and, more importantly, three types of invitro measurements have been performed to obtain a first experimental test of the validity of the theoretical results derived in this paper. The experimental measurement of the \(t_{delay}\) (Fig. 11) showed a standard deviation between theory and experiment of 16.2 s, which is 14% of the total range of variations in \(t_{delay}\). The experiment measuring the \(\Delta t_{central}\) (Fig. 12) showed a reasonable agreement between theory and measurements, in which standard deviation between theory and experiment was 8.4 s, which is 13% of the total range of variations in \(\Delta t_{central}\).
General findings from the expressions derived
In many of the resulting expressions as presented in “Results” section, the relative contribution of various factors affecting the dosing errors, such as the poiseuille mixing effect, can be discerned clearly in the structure of these expressions. The characteristic “time constant” \(R_{cath} C\) is clearly recognisable in the expressions. Furthermore, the ratio between the \(u^{(G)}_{downstep}\) (the size of the change in pump flow rate setting value) and the \(u_{final}\) (the stabilized final flow rate) appears as an important factor in determining the nature of the dosing error bolus, particulary how the dosing error will spread out in time (the “width” or “duration” (\(2\sigma \)) of the dosing error bolus). The contribution of the Poiseuille mixing effect is visible in the resulting expressions. Of particular clinical importance is the fact that the Poiseuille flow profile, once fully developed, causes a significant reduction of the time that is needed for a newly administered medication to reach the patient. This may come as an unexpected effect for the clinician, as may also the fact that the value of \(2\sigma \) (the “spread out”) is increased at the same time. The results in this paper may help to determine the magnitude of this “spreadout” effect as function of the hardware parameters (most notably, the characteristic \(R_{cath} C\) time), and e.g. the length and resistance of the catheter.
Limitations of the method, and possible extensions
A number of limitations can be identified in our method in its present form; most apparent is the assumption that when a time interval of duration \(t^{POIS}_{delay}\) has lapsed after having changed a pump flow rate setting value, the actual flow rate has already stabilized and reached the value \(u_{final}\). This needs not to be true in a general case. However, using the same reasoning as presented in this paper, our method can be extended to include nonstable flow rates at \(t=t^{POIS}_{delay}\) as well. Another limitation may arise by the fact that in our method we did not examine other elements (other than just syringes, infusion lines, catheters, and pumps) that may be present in an infusion setup, such as a nonreturn valve (preventing flow back from the catheter into a line towards a pump), or filters. The nonreturn valves may be incorporated into the method by restricting flow rates in the catheter to positive values or zero, whereas filters may be modeled on basis of a resistance and a mechanical compliance, as indicated in the literature (e.g. [21]). Another complicating factor may be the use of high viscosity fluids in the infusion setup. In our method, as it stands now, all fluids in the system were assumed to be of the same viscosity. Incorporation of deviating (high) viscosities into our method would entail an extension of the laminar flow profile used in this paper, because differences in viscosity of mixing liquids may produce destabilization of laminar flow. An easy, but useful, extension of our method is the incorporation of the possibility that the height at which the syringes are positioned in a setup near the patient, is altered during the period that a patient receives medication using the multiinfusion setup. Since a change in height may be modeled as the addition of an extra pressure source within the model during the process of infusion, this results in an extra bolus of dosing error, and, hence, in the addition of some extra terms and factors in the equations in “Results” section. Finally, the parabolic flow profile of a laminar flow does not come into existence instantaneously at the mixing point; it has been calculated [20] how long, c.q. what distance, it takes for the flow profile to approach the parabolic shape. All of these complicating factors need to be incorporated into the model to make it more realistic. As far as we can anticipate now, we do not expect these extra factors to be incompatible with our general approach outlined in this paper; however, further research is needed.
Potential use of the results in clinical practice
Healthcare professionals working with infusion technology in critical care have expressed the desire for a realtime tool that visualizes the multiinfusion drug therapy, e.g. continuously calculates predictions to indicate when the drugs will be entering the blood stream of the patient and in what dose. Such a tool may also visualize the causal consequences of an intervention (i.e., change in a pump flow rate setting value), before a clinician decides to proceed with such an intervention. In order to develop such a tool in the future, a fast and generic model will be necessary, combining all the relevant physical effects. The desirability of such a visualization tool, and the mathematical modeling that is a prerequisite for the development of such a tool in the future, is what prompted us to go beyond the stateoftheart and to develop the fully analytical method described in this paper. We envision three types of developments in which the results from this paper may be useful: (i) an interactive tool, running synchronized with the multiinfusion system on a smartphone device or on a bedside display, e.g., next to the vital signs display monitors. Such a tool could then be used in several cases, such as inotropic titration, or to visualize the effects of a syringe changeover, or even the consequence of changing the height of a pump during infusion. (ii) Another application of the method presented in this paper could be actual computer control of an infusion system. It has been shown that a computercontrolled pump with “knowledge” about the dead volume and the mixing effect within the dead volume can be useful in preventing overshoot [16]. Moreover, a control system with a feedback approach has also been attempted, where the mean arterial blood pressure was used to control the administration of a fastacting vasodilator [22, 23]. In both cases, however, increasingly complex situations and infusion setups were encountered where only the incorporation of all the interdependent physical effects, as described in this paper, would provide a sufficiently accurate prediction in order to make computer control feasible. (iii) The model can potentially also be used as an analysis and design tool, prior to investing in potential new infusion hardware. It is known that flow characteristics are influenced by valves [24, 25], syringes [26], infusion lines and catheters [27–29,28,], and filters [21]. By using the method from this paper, benefits of these components can be compared against potential tradeoffs.
Conclusions
We have developed a new method that combines mechanical compressibility (compliance) effects with poiseuille flow and pushout effects (“dead volume”) in one single mathematical framework for calculating dosing errors in multiinfusion setups.
In contrast to existing numerical methods, our method produces explicit expressions that indicate the mathematical dependencies of the dosing errors on hardware parameters and pump flow rate settings.
The results from the invitro experiments show a reasonable to good agreement between measurements and theoretical results.
The relative contribution of various factors affecting the dosing errors, such as the poiseuille mixing effect, resistance and internal volume of the catheter, mechanical compliance of the syringes and the various pump flow rate settings, can now be discerned clearly in the structure of the expressions generated by our method.
This enables insight and predictability in a large range of possible situations involving many variables and dependencies, which is potentially very useful for e.g. the development of a fast, bedside tool “calculator” that provides the clinician with a precise prediction of dosing errors and delay times interactively for many scenario’s. The interactive nature of such a device has now been made feasible by the fact that, using our method, explicit expressions are available for these situations, as opposed to conventional timeconsuming numerical simulations. Other potential applications of our method involve analysis and design tools for new infusion hardware, and interactive devices connected to the infusion hardware that counteract impending dosing errors using predictive calculations.
Abbreviations
ICU: intensive care unit; OR: operation room; MABP: mean arterial blood pressure.
List of symbols
 \(V_{cath}\) :

internal volume of the catheter
 \(t_{delay}\) :

time needed to travel through the catheter
 \(u_{cath}\) :

volumetric flow rate inside the catheter
 \(\mathcal {M}\) :

mixing point where fluid enters the catheter
 \(\mathcal {P}\) :

tip of catheter where fluid enters the blood stream
 L :

length of the catheter
 N :

large dimensionless number, representing the number of voxels inside the catheter
 \(\gamma \) :

L/N, i.e. length of a single voxel inside the catheter
 \(t=0\) :

point in time at which a change in pump flow rate setting value takes place
 \(a_k^{(R)}\) :

fraction of the volume of voxel k that contains “red” solution
 \(u_{patient}^{(R)}(t)\) :

partial flow rate of the “red” solution that enters the blood stream at time t
 \(u_{pump}^{(R)}\), \(u_{pump}^{(G)}\) :

pump flow rate setting of the “red”, and “green”, pump, respectively
 \(\Gamma \) :

complete set of flow rate settings and physical hardware parameters
 \(\beta _{patient}^{(R)}(t)\) :

dosing error of “red” solution, i.e. deviation from intended partial flow rate, as function of t
 \(t_{delay}^{POIS} = {{1}\over {2}}t{delay} \) :

time needed for the tip of the Poiseuille flow to travel through the catheter along the central longitudinal axis of the catheter
 \(u_{old}^{(G)}\) :

pump flow rate setting of the “green” pump before \(t=0\)
 \(u_{final}^{(G)}\) :

pump flow rate setting of the “green” pump after \(t=0\)
 \(\tau \) :

\(t  t_{delay}^{POIS}\)
 \(\beta _{pushout}^{(R)}(t)\) :

\(\beta _{patient}^{(R)}(t)\) for \(0< t < t_{delay}^{POIS}\)
 \(\psi _{patient}^{(R)}(\tau )\) :

\(\beta _{patient}^{(R)}(t) / (u_{final}^{(G)} + u_{pump}^{(R)} )\) for \(t_{delay}^{POIS} < t\)
 \(\lambda _{tot}(t)\) :

distance that the fluid would have travelled without the Poiseille mixing effect
 \(u_{\mathcal {M}}^{(R)}(t)\) :

the actual partial flow rate of “red” fluid entering the catheter at the mixing point \(\mathcal {M}\) at time t
 \(u_{\mathcal {M}}^{(R)diff}(t)\) :

\(u_{\mathcal {M}}^{(R)}(t)  u_{pump}^{(R)}\)
 \(w_i\) :

weight factor for contribution of original voxel nr i to the voxel at the tip \(\mathcal {P}\)
 \(b_{2j}\) :

\(a_j^{(R)diff}\)
 \(\psi _k^{(R)patient}\) :

discrete form of \(\psi _{patient}^{(R)}(\tau )\)
 \(\psi _k^{(R)POIS}\) :

part of \(\psi _k^{(R)patient}\) that is calculated in the Zdomain
 \(\Psi ^{(R)POIS}(z)\) :

Ztransform of \(\psi _k^{(R)POIS}\)
 Q :

total volume of the \(\psi _k^{(R)POIS}\)part of the dosing error
 \(t_{central}\) :

see Fig. 5c
 \(\sigma \) :

see Fig. 5c
 \(R_{cath}\) :

resistance of the catheter
 \(R_1\) :

resistance of the “green” feeding line
 \(R_2\) :

resistance of the “red” feeding line
 \(C_1\) :

mechanical compliance of the “green” syringe
 \(C_2\) :

mechanical compliance of the “red” syringe
 \(u_{downstep}^{(G)}\) :

change in pump flow rate setting of the “green” pump at \(t=0\)
 \(\vartheta _{first}\), \(\vartheta _{second}\), a, b, c :

see Eq. (34)
References
 1.
Foinard A, Decaudin B, Barthelemy C, Debaene B, Odou P. Impact of physical incompatibility on drug mass flow rates: example of furosemidemidazolam incompatibility. Ann Intensive Care. 2012;2:28.
 2.
Kern H, Kuring A, Redlich U, Dpfmer UR, Sims NM, Spies CD, Kox WJ, Dopfmer UR, Sims NM, Spies CD, Kox WJ. Downward movement of syringe pumps reduces syringe output. Br J Anaesth. 2001;86(6):828–31.
 3.
Weiss M, Fischer J, Neff T, Baenziger O. The effects of syringe plunger design on drug delivery during vertical displacement of syringe pumps. Anaesthesia. 2000;55:1094–8.
 4.
Timmerman AM, Snijder RA, Lucas P, Lagerweij MC, Radermacher JH, Konings MK. How physical infusion system parameters cause clinically relevant dose deviations after setpoint changes. Biomed Tech (Berl). 2015;60(4):365–76.
 5.
Sherwin CM, Medlicott NJ, Reith DM, Broadbent RS. Intravenous drug delivery in neonates: lessons learnt. Arch Dis Child. 2014;99:590–4.
 6.
Peterfreund RA, Philip JH. Critical parameters in drug delivery by intravenous infusion. Expert Opin Drug Deliv. 2013;10:1095–108.
 7.
Lovich MA, Kinnealley ME, Sims NM, Peterfreund RA. The delivery of drugs to patients by continuous intravenous infusion: modeling predicts potential dose fluctuations depending on flow rates and infusion system dead volume. Anesth Analg. 2006;102:1147–53.
 8.
Stowe CD, Storgion SA, Lee KR, Phelps SJ. Hemodynamic response to intentionally altered flow continuity of dobutamine and dopamine by an infusion pump in infants. Pharmacotherapy. 1996;16:1018–23.
 9.
Lovich MA, Wakim MG, Wei A, Parker MJ, Maslov MY, Pezone MJ, et al. Drug infusion system manifold deadvolume impacts the delivery response time to changes in infused medication doses in vitro and also in vivo in anesthetized swine. Anesth Analg. 2013;117:1313–8.
 10.
Snijder RA, Konings MK, Lucas P, Egberts TC, Timmerman AD. Flow variability and its physical causes in infusion technology: a systematic review of in vitro measurement and modeling studies. Biomed Eng Biomed Tech. 2015;60:277–300.
 11.
Lovich MA, Doles J, Peterfreund RA. The impact of carrier flow rate and infusion set deadvolume on the dynamics of intravenous drug delivery. Anesth Analg. 2005;100:1048–55.
 12.
Tsao AC, Lovich MA, Parker MJ, Zheng H, Peterfreund RA. Delivery interaction between coinfused medications: an in vitro modeling study of microinfusion. Paediatr Anaesth. 2013;23:33–9.
 13.
Decaudin B, Dewulf S, Lannoy D, Simon N, Secq A, Barthelemy C, et al. Impact of multiaccess infusion devices on in vitro drug delivery during multiinfusion therapy. Anesth Analg. 2009;109:1147–55.
 14.
Neff T, Fischer J, Fehr S, Baenziger O, Weiss M. Startup delays of infusion syringe pumps. Paediatr Anaesth. 2001;11:561–5.
 15.
Murphy RS, Wilcox SJ. The link between intravenous multiple pump flow errors and infusion system mechanical compliance. Anesth Analg. 2010;110:1297–302.
 16.
Parker MJ, Lovich MA, Tsao AC, Wei AE, Wakim MG, Maslov MY, et al. Computer control of drug delivery by continuous intravenous infusion. Anesthesiology. 2015;122:647–58.
 17.
James G. Advanced modern engineering mathematics. Reading: AddisonWesley; 1993.
 18.
Batista E, Almeida N, Furtado A, Filipe E, Sousa L, Martins R, Lucas P, Petter HT, Snijder R, Timmerman A. Assessment of drug delivery devices. Biomed Tech (Berl). 2015;60(4):347–57.
 19.
Taylor G. Dispersion of soluble matter in solvent flowing slowly through a tube. Proc R Soc Lond Ser A Math Phys Sci. 1953;219(1137):186–203.
 20.
Hutton P, Thornberry EA. Factors affecting delivery of drug through extension tubing. Br J Anaesth. 1986;58:1141–8.
 21.
Brotschi B, Grass B, Weiss M, Doell C, Bernet V. Inline filter included into the syringe infusion pump assembly reduces flow irregularities. Intensive Care Med. 2012;38:518–22.
 22.
Sheppard LC. Computer control of the infusion of vasoactive drugs. Ann Biomed Eng. 1980;8:431–44.
 23.
Ying H, Sheppard LC. Realtime expertsystembased fuzzy control of mean arterial pressure in pig with sodium nitroprusside infusion. Med Prog Technol. 1990;16:69–76.
 24.
van der Eijk AC, van der Plas AJ, van der Palen CJ, Dankelman J, Smit BJ. In vitro measurement of flow rate variability in neonatal IV therapy with and without the use of check valves. J Neonatal Perinat. 2014;7:55–64.
 25.
Hall JM, Roberts FL. An investigation into the reduction in flow rate of intravenous fluid by antireflux valves. Anaesthesia. 2005;60:797–800.
 26.
Weiss M, Hug MI, Neff T, Fischer J. Syringe size and flow rate affect drug delivery from syringe pumps. Can J Anaesth. 2000;47:1031–5.
 27.
Weiss M, Banziger O, Neff T, Fanconi S. Influence of infusion line compliance on drug delivery rate during acute line loop formation. Intensive Care Med. 2000;26:776–9.
 28.
Lovich MA, Peterfreund GL, Sims NM, Peterfreund RA. Central venous catheter infusions: a laboratory model shows large differences in drug delivery dynamics related to catheter dead volume. Crit Care Med. 2007;35:2792–8.
 29.
Ma H, Lovich MA, Peterfreund RA. Quantitative analysis of continuous intravenous infusions in pediatric anesthesia: safety implications of dead volume, flow rates, and fluid delivery. Paediatr Anaesth. 2011;21:78–86.
 30.
Potter M, Wiggert D. Fluid mechanics. Schaum’s outline series. New York: McGrawHill; 2008.
Authors' contributions
MKK has invented the new mathematical approach, performed the calculations, and largely wrote the manuscript, including all mathematical derivations. RAS has revealed the importance of Poiseuille flow characteristics in the evolution of dosing errors, and has initiated the incorporation of this effect into the model. Furthermore, RAS has placed the calculations within meaningful clinical context, and contributed significantly to the “Background” and “Discussion” sections and the invitro experiments. JHR has performed the invitro experiments, and has contributed to early stages of the modelling. AMT has initiated the study described in this manuscript, identified its targets, defined its clinical purpose, designed the invitro experiments, and coauthored the manuscript. All authors read and approved the final manuscript.
Acknowledgements
This research was partly funded by the EMRP project Metrology for drug delivery. The EMRP is jointly funded by the EMRP participating countries within EURAMET and the European Union.
Competing interests
The authors declare that they have no competing interests.
Data availability statement
Please contact author for data requests.
Ethics approval and consent to participate
Not applicable: No human or animal tissues involved.
Funding
This research was partly funded by the EMRP project Metrology for drug delivery. The EMRP has not been involved in the design of this particular study, and has not been involved in writing the manuscript, or in the collection, analysis, interpretation of the data.
Author information
Affiliations
Corresponding author
Appendix: Volume distribution in a Poiseuille flow
Appendix: Volume distribution in a Poiseuille flow
Laminar flow produces a Poiseuille flow profile, in which the velocity u depends on the distance r to the central longitudinal axis of the catheter [30]:
in which \(\eta \) is the viscosity of the fluid, \(\Delta p / L\) is the pressure drop over a length L of catheter, and \(2r_0\) is the diameter of the catheter.
Furthermore, the average velocity \(u_{avg}\) equals
and the maximum velocity (located along the central longitudinal axis of the catheter) equals
Evidently, if thermal diffusion is left out of the model and the fluid is considered to be completely homogeneous at the mixing point \(\mathcal {M}\), then the length \(\lambda _{max}(t)\) of the protrusion of the top of the Poiseuille profile into the catheter along the central longitudinal axis, is calculated easily by combining Eqs. (45) and (44):
Let x denote the position measured along the central axis along the catheter, starting withy \(x=0\) at \(\mathcal {M}\), and let \(r_{boundary}(x,t)\) denote the distance with respect to the central longitudinal axis at position x at time t, at which the boundary (parabola) is situated. Multiplying both sides of Eq. (43) with t, yields the relation:
in which x is a function of \(r_{boundary}\) and represents the distance that the boundary (see Fig. 13) has travelled in the longitudinal direction (along the length of the catheter), for a given radial position \(r_{boundary}\). The relation between x and r in Eq. (47) is depicted in Fig. 13:
Conversely, solving \(r_{boundary}\) from Eq. (47) for all values \(x \in [0, t \,\, u_{max}]\), we have:
which is depicted in Fig. 14.
In order to calculate the amount of the fluid R that is present inside a thin, diskshaped, crosssectional volume S of the catheter at position x, expressed as a fraction \(f_R(x,t)\) of the total fluid in the thin volume S, in which the diskshaped volume S is perpendicular to the central longitudinal axis of the catheter at point x, we need to integrate over the crosssectional area of the catheter at point x, from \(r=0\) to \(r=r_{boundary}(x,t)\), and divide this by the total crosssectional area \(\pi r_0^2\):
An interesting feature of this surprisingly simple result is that, within the region of the Poiseuille profile, there is a linear relationship between x and \(f_R(x,t)\).
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
About this article
Cite this article
Konings, M.K., Snijder, R.A., Radermacher, J.H. et al. Analytical method for calculation of deviations from intended dosages during multiinfusion. BioMed Eng OnLine 16, 18 (2017). https://doi.org/10.1186/s1293801603094
Received:
Accepted:
Published:
Keywords
 Infusion
 Catheter
 Dosing error
 Mathematical model
 Safety
 Poiseuille flow