Design of an in vitro mock circulatory loop to reproduce patient-specific vascular conditions: towards precision medicine

2/4/041004/975461/Design-of-an-In-Vitro-Mock-Circulatory-Loop-to DOI: https://doi.org/10.1115/1.4044488 *VOR (version of record)


Design of an in vitro mock circulatory loop to reproduce patient-specific vascular conditions: towards precision medicine
Patient-specific haemodynamic studies have attracted considerable attention in recent years due to their potential to improve diagnosis and optimise clinical treatment of cardiovascular diseases. Personalised computational models have been extensively investigated as a tool to improve clinical outcomes and are often validated against in vitro experimental data. Replicating patient-specific conditions in vitro is thus becoming increasingly important in cardiovascular research; experimental platforms can not only allow validation of in silico approaches but can also enable physical testing of various intervention scenarios and medical devices. Current experimental approaches suffer from shortcomings regarding personalisation and biomimicry. To address some of these limitations we have designed and developed a novel in vitro platform for the study of complex patient-specific vascular pathologies. This is achieved by using novel tunable 3element Windkessel vasculature simulators and a computer controlled pulsatile pump, coupled with mathematical models and computer routines to calibrate the parameters according to the available clinical datasets. In particular, the vessel inlet flow rate waveform and the afterload resistances and compliances are tuned in order to obtain target systolic and diastolic pressures, and cardiac output distribution. Pulse frequency (40-70 bpm), cardiac output (2-5 l/min), resistance (0.03-10.6 mmHg s/ml) and compliance (0.07-1 ml/mmHg) values have been tested and the overall reliability of the platform components as well as its computer routines to reproduce controlled physiological conditions demonstrated.

Introduction
In the past decades there has been an overall shift towards precision medicine in the bioengineering and clinical communities; an attempt to move away from the 'one-size fits all' approach for the management and treatment of patients. Personalised medicine is the tailoring of treatment to the individual characteristics of a specific subgroup of patients, or individual, with the goal of adapting drugs, medical devices and surgical procedures so that the best available option can be provided [1,2].
Patient-specific studies have become commonplace in different bioengineering areas, including the cardiovascular one. Vascular pathologies can be extremely complex; key parameters such as morphological features, flow patterns, pressures, velocities and shear rates are extremely patientspecific and, therefore, diagnosis management and treatment can vary greatly amongst individuals [3,4]. Aortic dissection, limb ischaemia and arteriovenous malformations, are all examples of complex clinically challenging conditions that could benefit from the development of precision vascular medicine tools and methodologies. Such methodologies can provide cost-effective treatment and design tailored medical devices, and have the potential to revolutionise healthcare [5][6][7].
Personalised computational models (CFD) have been extensively investigated as a tool to improve clinical outcome. Previous patient-specific in silico studies have studied a specific pathology [8,9], simulated the effect of a medical device or surgical procedure [10][11][12] and analysed the impact of personalised features on the haemodynamics [13,14]. The recent clinical adoption of computational tools in coro-nary disease [15] is a testament of the potential of such approaches.
In vitro studies have also evolved towards more complex, and personalised, setups capable of simulating haemodynamic physiological and pathological conditions, providing a valuable tool to develop algorithms to treat vascular pathologies, design medical implants, and validate the results of CFD studies. Various mock circulatory loops have been designed and employed since the seventies in order to simulate the human systemic circulation and applied to a variety of applications; for instance, to study mechanical or biological heart valves [16,17], to perform LVAD tests [18,19], to simulate surgical procedures [20] or to study the haemodynamics in different sites, such as in the ascending aorta [21], curved arteries [22], or stenotic vessels [23]. Pulsatile pump systems are typically employed to provide physiological flow rates [24][25][26] and afterload simulators are often used to reproduce the downstream vasculature, employing tunable resistances and compliances [27] or Windkessel models [28][29][30].
The importance of patient-specific features has been largely demonstrated, however in vitro simplified assumptions are often made. Significant efforts have been made to accurately reproduce the geometry of the vessel of interest from clinical images [31][32][33], and to mimic the mechanical properties of the tissue (i.e. elasticity of the vessel wall). However, the heamodynamics of a specific vascular pathology is dependent on the flow rates and pressures, as well as highly influenced by the boundary conditions, in a manner similar to patient-specific CFD modelling; hence, accounting for the morphology of a vessel alone cannot guarantee personalisation. Outlet boundary conditions are key to delivering physiological flows and oversimplification may limit reproducing the physiological downstream vasculature. Nonetheless, even though this is widely acknowledged, resistances or two-element Windkessel models are still employed in recent works, with their values tuned based on literature [21,30].
While in silico studies have demonstrated the use of dynamic boundary conditions tuned by in vivo data [8,34], to the best of our knowledge, equivalent in vitro apparatuses with both patient-specific tuned inflow and afterloads have not been reported, despite atempts to generate complex experimental mock loops [35][36][37].
In this context, a computer-controlled mock circulatory loop has been designed to emulate aortic vascular pathologies informed by in vivo patient-specific data. It comprises a pulsatile pump system and three element Windkessel afterloads and is coupled with mathematical models and computer routines translating available clinical data into tuning parameters for the physical components in order to reproduce personalised conditions.

Material and Methods
The ability to replicate a wide range of haemodynamic conditions and configurations was considered essential in order to simulate patient-specific features. This translated into two key design aspects: first, the components had to be tunable according to the available clinical data; second, the components had to be modular so that they could be assembled in several configurations according to the test section of interest (e.g. number of branches in the vessel). Therefore, three design constraints were considered: 1. The mock loop should be able to reproduce the correct inflow into the vessel. For aortic pathologies, this entails primarily the systolic phase of the cycle; 2. An afterload simulator is required in order to achieve physiological working conditions in the test section (i.e. physiological pressures and flows) and it should be tunable in order to reproduce different vessels and patients; 3. The mock loop should be coupled to computer routines to calibrate the input parameters to the physical models based on available clinical data. Figure 1a shows a schematic of the components typically found in mock circulatory loops reported in the literature. In this work, unlike the majority of published studies, a pulsatile pump system combined with a three element Windkessel simulator was used to reproduce in vitro vascular physiological conditions. Moreover, as represented in Fig.  1b, the designed devices can be tuned and assembled in different configurations to adapt to the case study under investigation.

Pulsatile pump system
Several types of pulsatile pumping systems can be found in reported mock circulatory loops, such as ventricular assist devices [4,21], collapsible ventricles [38], combined use of centrifugal pumps and piston pumps [25,39] and, lastly, the use of computer driven piston pumps only [22,40,41]. In this work, the latter option was chosen due to the necessity to reproduce the inlet flow curve of different patients as closely as possible.
The designed pump system is illustrated in Fig. 2, which shows its key elements, namely a mechanical actuator (motor, ball screw system and piston), a ventricular chamber and connection/support components. The mechanical actuator was designed based on previous work [42]. It consists of a brushless servo motor (AKM22E, Kollmorgen, VA, USA) electronically controlled by a driver unit (Single axis AKD drive, Kollmorgen, VA, USA). The rotating motor moves a ball screw (SH12.7-12.7R, SKF, Goteborg, Sweden) which transforms the rotational motion into the axial movement of a piston (CPA Series Cylinders, Waircom-MBS, Granby, QC, Canada). Compared to other options, the ball screw system was chosen to minimise friction and maximise efficiency. Figure 3 shows a schematic of the procedure adopted to dimension and design the pump system. Considering the design of the system, illustrated in Fig. 2, three components had to be initially dimensioned: (i) the motor, (ii) the piston, and (iii) the ball screw shaft and nut.
First, to support the design and dimension the whole mock circulatory loop, a simplified computational lumped parameter model was utilised. Namely, a 0D analogue of the in vitro circulation was simulated using 20-Sim (Controllab, The Netherlands), and solved using the backward Euler method. At this initial stage, physiological flow waveforms (i.e. input to the model) and resistance and compliance values for the afterload, were obtained from literature. The model provided an estimate for the pressure in the ventricular chamber required to specify the power of the rotational motor.
Subsequently, the stroke and diameter required for the piston were chosen according to the maximum Stroke Volume (SV) desired, which is equivalent to the maximum volume that the piston is able to move with a single forward movement of the gasket. The piston chosen for this work has been experimentally evaluated in a previous study [43].  Fig. 3. Flow diagram of the procedure followed to design the three main components of the pulsatile pump system: the piston, the ball screw system and the rotational motor The specifications for the ball screw system and motor were obtained as follows. The total force acting on the piston rod was obtained from: where F piston is the pressure force acting on the piston, F inertia is the inertia force to the acceleration of the piston, and F f riction is the friction force. The first term is calculated as: where P chamber is the pressure acting in the ventricular chamber estimated from the 0D analogue model, and A the cross sectional area of the chosen piston. The second force term can be obtained as: wherev is the acceleration of the piston, m piston its mass, anḋ Q the derivative in time of the flow rate. The last force term was experimentally evaluated in [42] to be 115 N, measuring the pressure in the cylinder at which the piston starts to move. Using the equations above, the required torque τ and angular speed ω were obtained from the force acting on the piston rod F rod and from the linear speed of the piston, v piston = Q A , respectively: where p screw and η are the screw pitch and efficiency, respectively. An efficiency value of 0.9 [42] was assumed for the mechanical transmission. The ball screw system (i.e. p screw ) was dimensioned based on the specifications of the motor as according to equations 4 and 5. Finally, once the maximum expected torque τ and angular speed ω were obtained, i.e. for the most demanding operating condition (i.e. highest total force), a suitable brushless motor was selected with a characteristic curve above the maximum operating condition requirements. Various support plates and connectors were used to ensure that the pumping system was stable, aligned and free of vibrations as shown in Fig. 2.
The mechanical actuator is connected to the ventricle unit, comprising a rigid chamber and two valve housings (i.e. inlet and outlet valves). The unit was designed to mimic the function of the left ventricle with the help of two one-way passive mechanical heart valves (tilting disk, Carbomedics): when the piston moves backwards, the inlet valve (i.e. mitral valve) opens and the fluid enters the chamber from the reservoir simulating the ventricle filling phase (diastole); when the piston moves forward, the inlet valve closes, the outlet valve opens (i.e. aortic valve) and the fluid is ejected into the aortic phantom, simulating the systolic phase.
The pump system is computer controlled and it includes a feedback mechanism on both velocity and position of the piston gasket. Safety procedures are implemented including over-current limiters, maximum allowed values of velocity, and position switch limiters. A mechanical stop button was also implemented as additional safety measure.

Afterload impedance simulator
In this study, a novel afterload simulator was designed making use of the 3-element Windkessel (WK3) model, to account for the downstream vasculature in vitro. The controllability of its parameters, ease of manufacturing and low encumbrance were considered the main requirements. Figure 4 shows the configuration of the impedance simulator, composed of five main elements: namely, the characteristic resistance (B), the main air chamber (C), the peripheral resistance (D) and the two additional air chambers (G and F). The fluid flows from the characteristic resistance, through the main air chamber and lastly into the peripheral resistance. Figure 4c shows the WK3 as an electrical analogue. The flow (Q) and the pressure (P) over this element are related by where R 1 and R 2 represent the proximal and distal resistances, respectively, and C is the compliance of the distal vasculature.

Compliance
Compliance was simulated using closed air chambers with a controlled volume. The initial volume of air V 0 is pressurised with an initial pressure P 0 . In the time interval required for the compliance to reach the steady state condition, the volume of fluid in the chamber increases until the operating pressure (P e ), corresponding to the operating volume V e , is reached. From this moment, a continuous oscillation of the fluid free surface is observed, due to the cyclic compression and expansion or air. As standard practice [42] , this process can be considered as a series of thermodynamic processes: an isothermal process followed by an adiabatic one. The first phase, which describes the slow filling phase of the chamber until V e and P e are reached, is represented by the isothermal process (i.e. P e V e = P 0 V 0 ), where the unknown variable of interest is the volume of air initially present in the chamber (V 0 ). The second phase, when oscillations around the operating pressure are observed, is approximated as an adiabatic process (i.e. PV k = P e V k e ). P and V indicate generic values of pressure and volume in the working range under consideration and k is the gas constant, which for bi-atomic ideal gas is equal to 1.4. The generic volume V is then obtained as: Knowing that: differentiation was performed in order to find the expression for the compliance of the system under consideration: Substituting Eq. 10 into Eq. 9 and re-arranging terms the following expression was obtained: The initial volume of air, to be used to control the value of compliance, was finally obtained: A computer routine was developed in Matlab (MathWorks, USA) to calculate the minimum and maximum initial volumes corresponding, respectively, to the minimum and maximum values of compliance of the considered range. To reproduce the range of interest, the compliance unit consisted of three different air chambers connected among them. The main central one, dimensioned with the minimum air volume needed, and the other two dimensioned so that the sum of the three volumes reproduces the maximum value of compliance. The compliance can therefore be adjusted by varying the total initial air volume in the chambers; this was done by partially filling the auxiliary chambers with water.

Resistances
The hydraulic resistances were designed so that a linear relationship between Q and ∆P is achieved; analytical calculations aided the design before manufacturing. Two design options were investigated to generate the required resistances: (i) bundles of parallel tubings, and (ii) generic porous material. A cylindrical geometry was considered in both cases.
Parallel tubings have commonly been employed in the literature to reproduce the resistance in Windkessel in vitro models [29,43], and hence considered first. Permeability (k) and resistance (R) parameters can be readily estimated combining Darcy's and Poiseuille's equations.
To facilitate patient-specific adjustment of resistance values experimentally and simplify the manufacturing process, a porous material generated resistance approach was also investigated. 3D printing was explored to manufacture the porous material with a geometry whose permeability can be varied and characterised analytically allowing the resistance to be computed using Darcy's equation. A computer routine was developed in Matlab ( Fig. 5) to evaluate the performance of the two approaches by generating all possible combinations of parameters to achieve a target range of resistance values and selecting the ones that allow both the minimum and maximum resistance values to be reached by varying only one variable at a time. The minimum and maximum target values were selected based on reported values from the literature and were set to 0.098 mmHg s/ml and 1.510 mmHg s/ml for the characteristic resistance and to 1.665 mmHg s/ml and 14.59 mmHg s/ml for the peripheral one, respectively. The total length of the resistance was chosen as the variable to be adjusted to modify the values because, among all the parameters, it was the easiest to control experimentally. The process described in Fig. 5 was applied to both resistance-emulating approaches considered and the relationship between resistance R and total length L was derived. The generic porous material option was found to reproduce both the characteristic and peripheral resistances in the impedance simulator and hence implemented. Two different geometries were designed and 3D printed for the two re-Variation of L tot of the porous material > R < R n 1 n 2 n 3 Fig. 6. Schematic illustrating the process used to tune resistance values by varying the total length of the porous material. The length is directly proportional to the resistance and its value is controlled through the number of 3D printed porous units assembled together sistances in the Windkessel physical model. A modular design approach was followed whereby several units were 3D printed via Selective Laser Sintering (SLS) to be assembled to the desired length (Fig. 6) for a given resistance value. SLS was chosen because it allows printing small structures without the need of support material.

Model personalisation
In order to deliver patient-specific, physiological flows the mock circulatory loop components described above were coupled with mathematical models and computer routines to calibrate the parameters to be given as input according to the specific vascular pathology and patient under investigation. This process was termed model personalisation and makes use of available clinical data from the patient. In this section the clinical dataset necessary to calibrate the parameters is described, and details regarding the personalisation methodologies developed for both the pump system and the Widkessel models are presented.

Clinical Dataset
The experimental platform was developed on the basis of non invasive, routinely acquired clinical data to be used to inform key input parameters in order to emulate vascular pathologies and create patient-specific haemodynamic flows.
The clinical data typically required to achieve this are as follows. First, volumetric medical images are required in order to reconstruct the vessel geometry and create realistic phantoms of the area of interest; these images can be acquired through CT or MRI. Information regarding the inflow conditions into the aorta can be obtained from 2D or 4D PC MRI or, when unavailable, from Doppler Ultrasound, in order to calculate the Cardiac Output (CO) and heart frequency, and subsequently generate a patient-specific flow wave starting from a typical aortic one. Lastly, two types of information are required in order to calibrate the outflow boundary conditions (i.e. physical Windkessel models) in this work: (i) the brachiocephalic pressure that can be used to obtain the target systolic and diastolic pressure values using the procedure secribed in [44], and (ii) the available flow rate waves  Finally, when the elasticity of the vessel is also taken into account, in vivo information of the vessel-area variation is also required (e.g. from 2D-cine-MRI).

Personalised flow wave generator
A mathematical routine to create the analytical waveform that reproduces the desired aortic flow rate was developed (Fig. 7).
The first part, in vivo parameter extraction, is aimed at obtaining the systolic time (T s ), the cycle duration (T c ) and the value of the systolic peak flow (Q peak ) from the in vivo flow wave. These parameters are used as inputs to the second part, the mathematical flow generator that creates a systolic waveform formula for a physiological and a simplified case, respectively. A simplified waveform is described by whereas for a physiological waveform the analytical expression introduced by [45] is used: Further adjustments can be performed to better match the slope of the target flow wave by changing the number of terms of the Fourier transform (Eq. 14). Once the correct systolic waveform has been generated in the time interval 0 < t < T s, the stroke volume is calculated through mathematical integration and the diastolic wave is subsequently created. If only the systolic phase is of interest, a simplified trapezium waveform can be used to minimise flow acceleration whereas for a physiological flow the analytical form proposed by [46] can be used: In both cases, the curve generated attempts to reproduce the diastolic phase with a smooth transition from the systolic one (i.e. same derivative) and is cyclically optimised so that the diastolic volume (i.e. integral of the diastolic waveform) matches the stroke volume previously calculated. Once this condition is verified, the full cardiac cycle curve is generated.
The third and last step in the personalisation procedure of Fig. 7, driver code and input parameter generator, serves as the interface between the obtained analytical curve and the code to be provided to the driver in order to produce the correct motion of the piston. Velocity, acceleration and displacement are calculated from the flow waveform generated in the previous step considering the bore of the piston. The maximum and minimum values of these parameters and the normalised displacement curve are extracted to be used as input to the driver.

Windkessel model parameter calibration
A physical Windkessel model is used in each outlet of the vessel phantom under investigation (e.g. four for an aortic phantom) and its parameters R 1 ,R 2 and C need to be individually defined. A tuning procedure developed in our previous work [8] to calibrate the Windkessel boundary conditions for a patient-specific CFD model of a diseased aorta was modified and coupled with the physical models as illustrated in Fig. 8. The goal of the tuning procedure was to obtain the target systolic (P sys ) and diastolic (P dia ) pressures, and the correct cardiac output distribution, set using available PC MRI data and literature data, at the outlets of the test section.
The process starts with extracting the geometry from clinical images (Step 1). The 3D CAD model of the patient vessel can be then divided into several segments, and each one can be modelled as a 0D circuit comprising one resistance (R) and one inductance (L) in order to account for the inertial and viscous effects. L is estimated as: where r is the vessel radius, l is the segment length and ρ is the blood density (i.e. ρ = 1060Kg/m 3 ). The resistance R is calculated taking into account both the geometry of the  vessel and the geometry of the output connectors (i.e. conic connectors). The first contribution is calculated as: where r is the vessel radius, l is the segment length and µ is the blood density (i.e. µ = 3.4 · 10 −3 Pas). In order to estimate the second contribution, the minor losses for each outlet connector have to be calculated as follows: where v is the mean velocity, calculated as Q/A, and k l is the minor losses coefficient defined as a value dependent on the angle and on the diameters' ratio. The corresponding pressure drop is obtained through Eq. 19: Finally, the second contribution of the resistance is obtained as the ratio between the pressure drop and the mean flow rate for each outlet.
In step 2 of the tuning procedure, a three-element 0D analogue (i.e. R sys,1 ,C sys , R sys,2 ) of the vascular system of interest is used to determine to total compliance of the system C sys . The model is informed using the flow rate wave as an input and the three parameters are varied iteratively until the target pressure values are obtained. A second 0D model reproducing the vascular system of interest is then used to tune R 1 ,C and R 2 for each Windkessel model (Step 3). The vessel of interest is represented through 0D segments of resistance and inertances, previously calculated in Step 1. The resistances and inertances due to hydraulic connectors and tubing present in the experimental setup are also taken into account.
The compliance of every branch C i is calculated by distributing the total compliance of the system (C sys ) among each W K3 i proportionally to the mean flow of the correspondent branch ( Q i ). The ratio R 1 R tot (where R tot = R 1 + R 2 ) is set equal to 5.6 [8]. R i tot is calculated as the ratio between P i and Q i , where P i is the mean pressure at the inlet of W K3 i obtained with the 0D model. R i tot is adjusted until the target flow distribution among the outlets is reached.
An additional stage (Step 4) is required when the phantom to be used in the in vitro experiments accounts for the vessel compliance. As described in [8], in this scenario, a multi-scale model of the 3D geometry and 0D peripheral circulation is needed. The parameters of the Windkessel models are then retuned considering that the compliance of the vessel has to be subtracted from C sys .
Once the values of the parameters of the Windkessel models are obtained, a Matlab routine was used to translate the resistance and compliance values into information for the physical models (Step 5). The initial volume of air to be used in the different compliance chambers is calculated from the compliance values and the number of 3D printed units to be used for each resistance is computed according to the total length needed for the porous material.

Functional assessment
The mock loop was tested to verify its overall mechanical and hydraulic behaviour. Pressure and flow rate waves were acquired in real time in different experiments using pressure transducers (Omega Engineering, UK) and an ultrasound flow meter (Sonotec, Halle, Germany) respectively. The signals were digitised at 200 Hz sampling frequency, recorded with an AD converter, and handled by a purposedesigned LabView virtual instrument. Figure 9 illustrates the locations where the data was acquired. A smoothing spline method was performed on the pressure curves (smoothing parameter p = 0.9999) to attenuate high frequencies.
First, the volumetric behaviour of the pump system, i.e. its ability to provide the desired flow regardless of the afterload or the imposed speed, was tested by performing experiments using the configuration shown in Fig. 9. A first test was performed without any afterload, imposing the same SV at different pulse frequencies ( f ) (i.e. 50, 60 and 70 bpm) and measuring the flow rate downstream of the aortic valve. In a second test, the behaviour of the system was tested at where t 1 and t 2 are the beginning and end time points of diastole, respectively. Experimental tests were performed to verify the efficiency of both the protocol to create a displacement curve from in vivo data, and the mechanical system movement. An exemplar target in vivo aortic inlet flow rate [8] was compared to the analytical curve constructed by the computer algorithm both in terms of shape and SV to verify the adaptability of the routine. The correspondent experimental flow rate acquired with the flow meter was also compared to the target one. Then, the reliability of the mechanical actuator was verified comparing the target velocity and displacement curves imposed to the driver to the ones obtained by the feedback of the motor, at 40, 60 and 70 bpm.
The behaviour of the Windkessel model, and its controllability, was analysed. Experiments were performed both by varying the inlet conditions ( f , 40-70 bpm; CO, 2-5 l/min)

Results and Discussions
The assembled pump system is shown in Fig. 10. The results obtained through the design calculations were in accordance with the ones reported in [42]. Considering an inlet Q=20 l/min, the obtained pressure in the ventricular chamber was 250 mmHg (33330 Pa) and the force acting on the piston (A = 50 mm) was estimated to be 65.41 N. Including the inertia and friction forces, the total estimated force was 183 N. Finally, increasing the pump system demand and considering an inlet Q=30 l/min, the correspondent total force was 202 N. In order to calculate the torque required by the motor, a force of 200 N was therefore considered and a screw pitch of 12.7 mm was chosen. Finally, according to equations 4 and 5, a torque τ=0.5 Nm and an angular speed ω=1420 rpm were found.
The system is able to reproduce a wide range of physiological and pathological SV. To accurately calculate the CO and guarantee agreement with the analytical calculations, the piston volume was experimentally evaluated and corrected due to a mismatch with the expected value. Indeed, because of the irregular geometry of the piston rod, the inner volume cannot be calculated ad a cylinder.
A series of experiments was conducted to evaluate the performance of the pulsatile pump system. Target displacement and velocity curves obtained by the mathematical flow wave generator were successfully compared to the mechanical actuator's feedback and a maximum error of 3.2% was obtained at 70 bpm in the displacement peak value because of the inertia of the physical system. No difference was found in the pulse frequency value for the tested range. At different pulse frequencies, a negligible error (0.78%) was found between the minimum and maximum SV calculated from the curves, demonstrating the volumetric behaviour of the system. Moreover, a very good overlap was found between the flow rate waves acquired in the absence and presence of an afterload, both upstream the mitral valve and downstream the aortic valve. The passive behaviour of the  valves was tested and no signs of wear or misalignment were found when inspected. As expected, the AVBF value increased proportionally to the CO and to the diastole duration. It was calculated as percentage of SV as about 1.5% for the static (SBF) and 5% for the dynamic one (DBF) at 70 bpm and CO = 4 l/min. The personalised flow generator algorithm, coupled with the pump system, was found to be able to successfully identify the correct diastolic equation to reproduce the same systolic SV, and to calculate the analytical cycle equation.  The computer routine was tested in the range 40-80 bpm and 1-6 L/min. As an example, Fig. 11a shows the comparison between a target flow rate and three different analytical curves obtained with the computer routine. The first analytical curve (Analytical 1) was reproduced with the simplified function. To achieve a closer match between the generated and target curves, the equation for the systolic waveform was modified, increasing the number of terms of the Fourier transform (Eq. 14), creating the analytical curves 2 and 3. These better reproduce the target systolic flow wave demonstrating the adaptability of the flow generator routine. A smooth slope transition between the systolic and diastolic phase can also be noted in the figure. It can be noted that in this case the best match is not reached with the commonly adopted physiological expression (Analytical 2, Eq. 14), demonstrating that to reproduce patient-specific conditions further adjustments have to be made. Furthermore, unlike other pump systems described in the literature that employ a fixed relation between the duration of the systole and the total cycle, the developed flow curve generator adapts not only to the desired pulse frequency, SV and CO, but also to different input duration ratios, increasing the personalisation of the output curve.
It should be noted that while the target of the systolic phase is the correct match of the shape and values of the curve, the diastolic phase (grey area in Fig. 11) is here approximated with a trapezoid calculated only to reproduce the correct SV and thus allowing the piston to return to its initial position. The negative flow rate values observed in the figure in this phase do not represent a back flow through the aortic valve, but an inflow from the mitral valve, allowing the filling phase of the simulated ventricle. The physiological shape of the resulting experimental flow rate curve is indeed achieved thanks to the passive action of the two one-way valves. Figure 11b shows a comparison between the analytical flow rate curve calculated by the algorithm and the correspondent one acquired in vitro. It can be noted that during the diastolic phase the closure of the aortic valve enables the piston to return to its initial position while obtaining a physiological small back flow in the aorta. The higher difference in the SV values between the target flow rate curve and the acquired one increases proportionally with the CO and is mainly due to the AVBF (about 6% at 70bpm and CO=4 l/min).
Similarly to most piston-based pulsatile pumps, the employment of passive valves does not allow continuous controlled positive flow output because the piston has to move backwards at the end of each stroke. This limitation could be circumvented by implementing automatic computer control of valves in relation to the piston position.
The fine personalisation of the systolic phase of the inlet flow rate curve (acceleration of the systolic phase, SV, peak flow etc) represents a significant advancement of in vitro flow generators, essential to reproduce patient-specific conditions, as demonstrated by the use of personalised curves as inlet conditions of in silico models.
The final design of the Windkessel simulator is shown in Fig. 12. The simulator is reasonably compact measuring 310 mm (length) x 134 mm (height) x 66 mm (width). Leak proof operation was ensured by testing the models under air pressure of up to 300 mmHg. The modular 3D porous units used to control the resistance values R 1 and R 2 by varying the total material length can also be seen in the figure. This novel approach was found not only to be very efficient in changing the resistance values but also highly versatile as various geometries can be designed and 3D printed to emulate resistances in different target areas of the human circulation. After the calibration procedure, fine parameter adjustment was performed to compensate for the resistances and inductances caused by the physical components of the system (i.e. hose connectors, tubings).
The developed Windkessel model addresses some of the limitations of previously reported systems (based on 2element Windkessel models) by simulating also the characteristic impedance of the proximal aorta that allows, in particular, to capture the behaviour at high frequencies. The choice of using a 3-element circuit in order to reproduce the afterload is supported by [28], where it was identified as the best compromise between accuracy and complexity of the system. Moreover, the compact design, contained dimensions and ease of controllability of the parameters of the model developed in this work, allow the personalisation of in vitro apparatus, overcoming some of the limitations of previously reported devices. Indeed, several Windkessel models can be manufactured and tuned to represent different boundary conditions (corresponding to different branches) of the vessel under investigation. Varying the resistance and compliance values in the analysed range allowed to obtain the target systolic and diastolic pressure values with a maximum difference of about 5 mmHg.
The overall behaviour of the platform and its ability to reproduce physiological and pathological conditions of the systemic circulation were tested. The platform was found to be fairly robust and able to reproduce pressure and flow waveforms of the circulation under different experimental conditions (i.e. f, SV, Afterload). Figure 13 shows an example of in vitro aortic and ventricular pressures, and aortic flow waveforms. The obtained curves were comparable with the known physiological ones. The oscillations observed in the experimental unfiltered data are due to the mechanical action of the valves and are present at the beginning of the systolic and diastolic phases. Experimentally, the difference of the pressure values between the ventricular and the aortic pressures during systole can be attributed to the pressure drop caused by both the mechanical valve and the conical connector which is directly proportional to Q 2 . Such difference is indeed minimised at low CO when, indeed, the two curves present the same shape and values. The additional resistance caused by the physical connector does not represent a limitation when investigating the haemodynamics of downstream vessels. However, when the aortic valve dynamic is of inter-est, it can be minimised by varying the design or simulating the aortic root.
The controllability of all the parameters of the platform and the ability to combine different devices according to the test section of interest allows the reproduction of a wide range of conditions. The obtained experimental behaviour of the resistances and compliances in relation to the flow rates and pressures corresponds to the known theoretical one. For instance, when low resistances and compliances are present, and the inertance is the main observed effect, the pressure wave follows the derivative of the inlet flow rate curve; an increase of the resistance values follows an increase of the average pressure while higher compliances result in a smaller pulse pressure and slower diastolic pressure decay.
Depending on the complexity of the test section under investigation (e.g. vessel geometry and mechanical properties, number of branches considered), the time required to produce the physical representation from patient-specific data varies from a minimum of 4 days to about 2-3 weeks.

Conclusions
In vitro experiments are fundamental for the development of cardiovascular devices, haemodynamic investigations of pathologies, surgical intervention simulation and validation of CFD models. Their reliability and efficacy depends on the ability to both perform highly controlled tests and reproduce physiological conditions.
In this work, a novel approach for in vitro personalised models is presented. A mock circulatory loop coupled with mathematical models and computer routines, and informed by in vivo data to reproduce patient-specific conditions was designed and developed. The apparatus was able to reproduce physiological flow rates and pressures under several experimental conditions and the observed behaviour of the system under different configurations of parameters (resistances, compliances, CO, etc) was in agreement with the known theoretical one. The computer algorithms that support the personalisation of the devices proved to be efficient, both in terms of ease of applicability and quality of the results; and allowed not only to reproduce physiological conditions but also to adapt in order to match a specific one. The platform is currently being employed to develop a case study of a complex vascular pathology, aortic dissection, and to test third party prototypes under physiological flow rates and pressures. Future works also include the validation of in silico patient-specific models with complex boundary conditions.
The controllability of the devices and the personalisation approach represent a significant advancement in experimental mock circulatory loops, bringing in vitro setups closer to patient-specific in silico simulations, moving towards personalised medicine.
work is gratefully acknowledged. The authors are grateful to Professor Pat Lawford for providing the mechanical heart valves used in this work.