Growth Behavior of Human Adipose Tissue-Derived Stromal/Stem Cells at Small Scale: Numerical and Experimental Investigations

Human adipose tissue-derived stromal/stem cells (hASCs) are a valuable source of cells for clinical applications, especially in the field of regenerative medicine. Therefore, it comes as no surprise that the interest in hASCs has greatly increased over the last decade. However, in order to use hASCs in clinically relevant numbers, in vitro expansion is required. Single-use stirred bioreactors in combination with microcarriers (MCs) have shown themselves to be suitable systems for this task. However, hASCs tend to be less robust, and thus, more shear sensitive than conventional production cell lines for therapeutic antibodies and vaccines (e.g., Chinese Hamster Ovary cells CHO, Baby Hamster Kidney cells BHK), for which these bioreactors were originally designed. Hence, the goal of this study was to investigate the influence of different shear stress levels on the growth of humane telomerase reversed transcriptase immortalized hASCs (hTERT-ASC) and aggregate formation in stirred single-use systems at the mL scale: the 125 mL (=SP100) and the 500 mL (=SP300) disposable Corning® spinner flask. Computational fluid dynamics (CFD) simulations based on an Euler–Euler and Euler–Lagrange approach were performed to predict the hydrodynamic stresses (0.06–0.87 Pa), the residence times (0.4–7.3 s), and the circulation times (1.6–16.6 s) of the MCs in different shear zones for different impeller speeds and the suspension criteria (Ns1u, Ns1). The numerical findings were linked to experimental data from cultivations studies to develop, for the first time, an unstructured, segregated mathematical growth model for hTERT-ASCs. While the 125 mL spinner flask with 100 mL working volume (SP100) provided up to 1.68 × 105 hTERT-ASC/cm2 (=0.63 × 106 living hTERT-ASCs/mL, EF 56) within eight days, the peak living cell density of the 500 mL spinner flask with 300 mL working volume (SP300) was 2.46 × 105 hTERT-ASC/cm2 (=0.88 × 106 hTERT-ASCs/mL, EF 81) and was achieved on day eight. Optimal cultivation conditions were found for Ns1u < N < Ns1, which corresponded to specific power inputs of 0.3–1.1 W/m3. The established growth model delivered reliable predictions for cell growth on the MCs with an accuracy of 76–96% for both investigated spinner flask types.

Disposable Corning ® spinner flasks (Corning, Corning, NY, USA), which are commercially available in different sizes (125 mL and 500 mL, see Figure 1), were used for all investigations in this study. These rigid culture containers are made of polycarbonate and were equipped with two angled side ports. The slightly opened lids provided gas exchange (CO2/O2) for surface aeration in a standard cell culture incubator during the cultivations. The main physical dimensions and ratios of the spinner flasks are summarized in Table 1. The maximum working volumes that were used for the cultivation studies were 100 mL (= SP100) and 300 mL (= SP300). Both spinner flasks were equipped, as standard, with paddle-like impellers consisting of a blade and a magnetic bar. As can be seen in Table 1, the used spinner flasks have comparable geometrical ratios.

CFD
For all CFD simulations, the fluid flow and the MC distribution were modeled using the Fluent 16.2 finite volume solver (ANSYS Inc., Canonsburg, PA, USA). The numerical technique was based on the subdivision of the fluid domain into a finite number of control volumes and the discretization of the time-averaged mass and momentum equations. This approach provided the algebraic equations, which were solved iteratively [19]. All simulations were run in parallel and solved on a These rigid culture containers are made of polycarbonate and were equipped with two angled side ports. The slightly opened lids provided gas exchange (CO 2 /O 2 ) for surface aeration in a standard cell culture incubator during the cultivations. The main physical dimensions and ratios of the spinner flasks are summarized in Table 1. The maximum working volumes that were used for the cultivation studies were 100 mL (=SP100) and 300 mL (=SP300). Both spinner flasks were equipped, as standard, with paddle-like impellers consisting of a blade and a magnetic bar. As can be seen in Table 1, the used spinner flasks have comparable geometrical ratios.   the time-averaged mass and momentum equations. This approach provided the algebraic equations, which were solved iteratively [19]. All simulations were run in parallel and solved on a computational cluster (up to 16 Intel ® Xeno ® E5-2630 v4 CPU's @ 2.2 GHz, 64 GB RAM) in order to speed up the computational turnaround time.

Euler-Euler (EE) approach:
Multiphase simulations were carried out using the EE Reynolds-averaged Navier-Stokes (RANS) approach, which considered water (=culture medium) as the continuous and the MCs as the dispersed phase. The continuity equation for the q th phase was written as shown in Equation (1), where → u q , ρ q , α q (q = L, the liquid phase and q = P, the particle phase) represented the phase velocity vector, the phase density and the phase volume fraction, respectively. The phases were assumed to share space in proportion to their volume fractions (see Equation (2)), whereas the maximum volume fraction of the dispersed phase was restricted to 0.63 (=maximum packing limit) due to the spherical shape of the MCs.
The conservation equation for the momentum of the qth phase was based on an extended Navier-Stokes equation (see Equation (3)). The momentum exchange resulted in the coupling of the two phases.
In addition to all of the forces (i.e., viscous stresses, overall pressure gradient, gravitational force, interphase momentum forces) acting on a fluid element of the qth phase in the fluid domain, the drag force was considered to be the most important interphase force. In general, drag force results from the relative velocity between the two phases (see Equation (4)).
The drag coefficient C D was modeled in all of the EE simulations using the Symlal & O'Brien sub-model [34]. The Symlal and O'Brien sub-model was used due to numerical stability issues in the granular model. However, the computed C D values did not differ significantly from the expected Re P range compared to values derived from the standard correlation given by Schiller and Neumann [35], which was used for the Euler-Lagrange simulations (see Euler-Lagrange approach). The term u r,p in Equation (5) represents the terminal velocity correlation for the solid phase and is, therefore, dimensionless (see also References [34,36]). C D = (0.63 + 4.8 Re P /u r,P ) 2 (5)

Euler-Lagrange (EL) approach:
The EL model is suitable for dispersed flows, where the particles are non-homogeneously distributed. This model approach combines the description of the continuous phase (=culture medium) with a segregated description of the dispersed phase (=MCs) in the Lagargian frame. Each particle in the flow was characterized by its location ( → χ P ), velocity ( → u P ), and other mechanical and thermodynamic variables, while the fluid phase was treated as a continuum by solving the Navier-Stokes equations. The computational effort for the EL model was higher than the EE model because of the separate treatment of each particle. Therefore, the number of particles in each simulation was set between 116,000 and 125,000 (phase fraction = 0.1-0.2%). The velocities of the dispersed phase were obtained by integrating the force balance on each individual particle. Consequently, this approach equates the particle inertia with forces acting on the particles and can be written as shown in Equation (6) (e.g., for the x-direction in Cartesian coordinates).
The term → F χ denotes additional forces in the particle force balance, such as the Coriolis force, centrifugal force, virtual mass force, Saffman lift force, Basset force, Magnus force, and pressure gradient-dependent forces. As mentioned for the EE model, the drag force was considered to be the important interphase force. Therefore, the drag force → F D was also introduced in the EL simulations. The drag coefficient C D was calculated using the standard correlation given by Schiller and Neumann [35] (see Equation (7)).

C D =
Re P ≤ 1000 : 24 Re P ·(1 + 0.15·Re 0.687 P ) Re P > 1000 : 0.44 (7) Numerical details and boundary conditions: The EE and EL simulations were carried out for maximum working volumes of either 100 mL (=SP100) or 300 mL (=SP300). In both cases, a hybrid mesh consisting of unstructured tetrahedral elements and a prism layer at the vessel walls was used. The meshes were generated using the ANSYS Meshing Tool, which was implemented in ANSYS Workbench 16.2 (ANSYS Inc., Canonsburg, PA, USA). A previous grid sensitivity study had confirmed that grids used with 710,000 (SP100) and 2,100,000 (SP300) fluid elements result in grid-independent results (data not shown). For both spinner flasks, the vessel and impeller walls were treated as non-slip boundaries with standard wall functions. The impeller rotation was implemented using the sliding mesh approach and the turbulence was modelled using the k-ω SST or k-ε realizable turbulence models, depending on the expected Re range. All multiphase simulations were performed transiently (dt = time corresponding to 1-3 • of impeller motion) and simulation convergence within a time step was assumed when the residuals dropped below 10 −5 . The simulated time corresponded to at least three times the mixing time (Θ 95% ) of the individual spinner flask (SP100 Θ 95% = 7-17 s/25-120 rpm; SP300 Θ 95% = 6-17 s/20-100 rpm) at the defined impeller speed. The phase-coupled SIMPLE algorithm (Semi-implicit Method for Pressure Linked Equations) [36] was used for pressure-velocity coupling in all cases.

Segregated Growth Model (SGM)
Based on the findings from the CFD simulations (e.g., P/V, τ nt ) and their combination with the cell culture data, an unstructured, segregated, simplistic growth model was developed. The general concept of the growth model and the influencing factors are shown in Figure 2. Because hMSC growth is anchorage-dependent, the cell aggregation and the formation of spheroids in the suspension was ignored. It was assumed, that cells in suspension do not contribute to the increase in overall cell numbers, with cell growth restricted to the MC surface. However, cells in suspension do affect the glucose balance, since they still have a maintenance metabolism. To define the starting conditions, it was assumed that initial cell attachment took place during the cell attachment phase. The cell attachment rate was described by the constant k at (see Table 2). Once the cells had attached to the MC surface, it was assumed that they immediately start to proliferate.
The specific cell growth rate (µ) was calculated based on Monod-type kinetics, with the consumed substrate (Glc), the produced metabolites (Lac, Amn), and the available growth surface (X max ) considered as influencing factors (see Equation (8)). A comparable approach was already successfully However, the detachment constant −k det is strongly affected by hydrodynamic forces, and therefore, variable for different specific power inputs. As mentioned before, cell growth in the suspension was negligible and therefore changes in cell concentration will only be affected by cell attachment to or detachment from the MC surface (see Equation (10)).
Bioengineering 2018, 5, 106 6 of 30 The specific cell growth rate (μ) was calculated based on Monod-type kinetics, with the consumed substrate (Glc), the produced metabolites (Lac, Amn), and the available growth surface (Xmax) considered as influencing factors (see Equation (8)). A comparable approach was already successfully used by Möhler et al. [37] and Bock et al. [38] to simulate anchorage-dependent Madin-Darby Canine Kidney (MDCK) cell growth in a MC-based culture.
The concentration of the cells on the MC surface increased through mitotic cell division and the attachment of cells from the suspension (see Equation (9)). However, this increase was reduced by cell detachment from the MC surface and was accounted for by the detachment constant −kdet.
However, the detachment constant −kdet is strongly affected by hydrodynamic forces, and therefore, variable for different specific power inputs. As mentioned before, cell growth in the suspension was negligible and therefore changes in cell concentration will only be affected by cell attachment to or detachment from the MC surface (see Equation (10)). Contrary to the growth restriction based on the specific growth rate, glucose consumption was only limited by the glucose concentration itself (see Equation (11)). Consequently, glucose consumption was the result of the glucose uptake by the mitotic cells and the maintenance metabolism of the mitotic and non-mitotic cells (X v ).
In order to avoid negative glucose concentrations, a step response was implemented (see Equation (12)).
Due to stability issues in the culture medium, ultra-glutamine (alanyl-L-glutamine) was used in all of the cultivation studies. However, glutamine consumption was not considered in the growth model, since our previous cultivation studies have shown that glutamine (Gln) is not a limiting factor in the medium, especially under the investigated conditions (data not shown). The production of the two main metabolites (Lac, Amn) was accounted for by a growth associated assumption, i.e., they were only produced during cell growth when sufficient substrate concentration and free growth surface were available (see Equations (13) and (14)). This meant that any increase in concentration due to maintenance metabolism was ignored.
All calculations and simulations were performed with MATLAB R2018a (MathWorks Inc., Natick, MA, USA). The model equations were solved using the ode15s solver.

Suspension Studies
The suspension criteria (N s1u , N s1 ) for the ProNectin ® F-COATED MCs (ρ P 1026 ± 4 kg/m 3 , d P 169 ± 43 µm, A P 360 cm 2 /g) were determined in both spinner flasks. The methodology for the determining of the suspension criteria was in accordance with Kaiser et al. [25]. In brief, the N s1 (=N js or just suspended) suspension criterion was defined as the impeller speed required to just fully suspend the MCs in the spinner flasks [40]. N s1u described the suspension state at which some of the MCs were still in contact with the reactor bottom, but none of them were at rest [41]. The suspension experiments were carried out at different MC concentrations (2.5-20 g/L) and with a specially developed cell culture medium from Lonza (w/5% FBS). While the impeller was in motion, the suspension state was recorded by two digital cameras (from the side and below) and the recordings were subsequently evaluated by visual observation. The optical accessibility to the spinner flask bottom was improved by a mirror which was placed below the flask.

Particle Image Velocimetry (PIV)
Stereoscopic PIV measurements were carried out using a FlowMaster PIV system (LaVision, Göttingen, Germany) in order to verify the CFD simulation results. The illumination of the field of investigation was performed by a double-pulsed Nd:YAG laser (Litron Laser Ltd., Rugby, UK), which generated a laser light sheet at a wavelength of λ 532 nm. In both cases, the laser light sheet was vertically aligned through the spinner flasks and allowed radial, axial, and tangential fluid velocity components to be measured (see Figure 3a). The fluid flow field was measured at different impeller positions by phase-resolved measurements using a photoelectric barrier (see Figure 3b). Two Imager Pro X4 CCD cameras (LaVision, Göttingen, Germany) with a resolution of 2048 × 2048 pixels were used for image capturing. A Scheimpflug set-up with backward/forward scattering was used to align the cameras, which were spatially calibrated using a two-level calibration plate (106-10, LaVision, Germany). In order to reduce the effect of refraction/diffraction, which was induced by the laser light beam hitting the cylindrical spinner flask surface, the spinner flasks were placed in a water-filled cubical box. The DaVis 8.3 (v 8.3, LaVision, Göttingen, Germany) software was used for image acquisition and image processing, with up to 1,500 double-frame images per camera considered for analysis by cross correlation. Cross correlation was performed by means of 32 × 32 pixel-interrogation windows and an overlap of 50%. Fluid flow visualization was achieved by adding rhodamine coated polymethylmethacrylat beads (PMMA, 1190 kg m −3 , 20-50 µm, λ EM,max 584 nm) to the spinner flasks. Finally, laser light reflections were reduced by means of corresponding long pass optical fluorescence edge filters (100% transmission >545 nm) mounted on each CCD camera.

Microcarrier Measurement by Shadow Imaging (Shadowgraphy)
To verify the predicted MC characteristics from the multiphase simulations, MC distribution and velocities were measured experimentally by means of shadow imaging techniques (see Figure  4a). For this purpose, the ParticleMaster shadowgraphy system (LaVision, Göttingen, Germany) in conjunction with a double-pulsed Nd:YAG laser and a high-efficiency light diffusor (λ 550-600 nm) were used. The light beam was aligned parallel to the impeller shaft and directly opposite the CCD camera (Imager Pro X 4M, LaVision, Göttingen, Germany), which was equipped with a telephoto

Microcarrier Measurement by Shadow Imaging (Shadowgraphy)
To verify the predicted MC characteristics from the multiphase simulations, MC distribution and velocities were measured experimentally by means of shadow imaging techniques (see Figure 4a). For this purpose, the ParticleMaster shadowgraphy system (LaVision, Göttingen, Germany) in conjunction with a double-pulsed Nd:YAG laser and a high-efficiency light diffusor (λ 550-600 nm) were used. The light beam was aligned parallel to the impeller shaft and directly opposite the CCD camera (Imager Pro X 4M, LaVision, Göttingen, Germany), which was equipped with a telephoto lens (Nikon, Tokyo, Japan). Similar to the PIV measurements, the spinner flasks were placed in a water-filled cubical box to reduce the effects of refraction/diffraction. DaVis 8.3 software (v 8.3, LaVision, Göttingen, Germany) and ParticleMaster toolbox were used for image acquisition and analysis. Particle recognition was performed based on an image segmentation algorithm with subsequent analysis of the pixel intensity profile [42]. To calculate the MC velocities, individual particles were tracked over two images and the velocities were calculated based on the corresponding pixel shift in the x-y directions (v xy ). Since only one CCD camera was used for the investigations, the direct measurement and calculation of → w was not possible. However, the effects of the particle shift in the z-direction were considered by means of depth-of-field calibration. Depending on the size of the spinner flask, 1,500 double-frame images were captured from up to 4 positions in the vertical direction. Spatial scanning of the fluid domain in the vertical direction was performed by a traverse system. All measurements were carried out with cell-free and cell-loaded MCs at different impeller speeds (N s1u /N s1 ) and impeller positions (see Figure 4b). The cell-loaded measurements were immediately performed at the end of 7 days of cultivating hTERT-ASCs. Synchronization of the camera and the laser to the impeller motion was performed by means of a trigger signal obtained from a photoelectric barrier in order to perform phase-resolved measurements.
Bioengineering 2018, 5, 106 9 of 30 lens (Nikon, Tokyo, Japan). Similar to the PIV measurements, the spinner flasks were placed in a water-filled cubical box to reduce the effects of refraction/diffraction. DaVis 8.3 software (v 8.3, LaVision, Göttingen, Germany) and ParticleMaster toolbox were used for image acquisition and analysis. Particle recognition was performed based on an image segmentation algorithm with subsequent analysis of the pixel intensity profile [42]. To calculate the MC velocities, individual particles were tracked over two images and the velocities were calculated based on the corresponding pixel shift in the x-y directions (vxy). Since only one CCD camera was used for the investigations, the direct measurement and calculation of was not possible. However, the effects of the particle shift in the z-direction were considered by means of depth-of-field calibration. Depending on the size of the spinner flask, 1,500 double-frame images were captured from up to 4 positions in the vertical direction. Spatial scanning of the fluid domain in the vertical direction was performed by a traverse system. All measurements were carried out with cell-free and cell-loaded MCs at different impeller speeds (Ns1u/Ns1) and impeller positions (see Figure 4b). The cell-loaded measurements were immediately performed at the end of 7 days of cultivating hTERT-ASCs. Synchronization of the camera and the laser to the impeller motion was performed by means of a trigger signal obtained from a photoelectric barrier in order to perform phase-resolved measurements.

Cells, Microcarriers, and Medium
Human telomerase reversed transcriptase immortalized hASCs (hTERT-ASCs) from the American Type Culture Collection (SCRC-4000 TM , ATCC, Manassas, VA, USA) were used for all cultivation experiments. Prior to the cultivation studies, the hTERT-ASCs (p = 23, PDL 33) were adapted to the serum-reduced cell culture medium (5% FBS) from Lonza (data not shown) and a cryovial-based working cell bank was established (storage in liquid nitrogen). For the inoculation of each spinner flask, an initial average cell density of 3,000 cells/cm 2 (corresponding to 10,800 cells/mL) and a MC concentration of 10 g/L (ProNectin ® F-COATED, Pall SoloHill, New York, NY, USA) was used. The MC concentration of 10 g/L was selected based on previous investigations of Schirmaier et al. [16]. The required number of MCs was prepared and sterilized as recommended by the vendor one day before inoculation.

Analytics
Off-line samples were taken daily to measure the substrate (Glc, Gln) and metabolite (Lac, Amn) concentrations with a BioProfile 100Plus (Nova Biomedical, Waltham, MA, USA) and/or a CedexBio (Roche Diagnostics, Risch-Rotkreuz, Switzerland). After the cells had detached from the MC surface (15 min with TrypLE Select; Gibco by Life Technologies, Carlsbad, CA, USA), the hTERT-ASC cell number was measured (n = 3) using a NucleoCounter ® NC-200 TM (ChemoMetec, Allerød, Denmark). The measured cell specific values were used to calculate the growth-related parameters (maximum specific growth rate, µ max , minimum doubling time, td, expansion factor, EF). In addition to the cell measurements, 2 mL of the MC-cell suspension was fixed immediately after sampling with a 3% paraformaldehyde solution for 4',6-diamidin-2-phenyliondol (in short DAPI) staining (see Figure 5a). The fixed MC-cell suspension was also used to analyze the MC-cell aggregates (see Figure 5b). For this purpose, the MCs were scanned prior to staining with a document scanner (Epson Perfection 1650) with an image resolution of 4,200 dpi (=5 µm/pixel). The images were then processed and analyzed with ImageJ (particle analysis toolbox) and Matlab in order to obtain the MC-cell aggregate size distribution (see Figure 5c Flow cytometric measurements were carried out with cells from the inoculum and with microcarrier-free, purified cell samples from the different spinner flasks. For flow cytometric analysis, the cells were stained with fluorochrome-conjugated anti-human CD14, CD20, CD34, CD45, CD73, CD90, and CD105 (according to the recommendations of the International Society for Cellular Therapy ISCT and the International Federation for Adipose Therapeutics and Science IFATS) antibodies (Miltenyi Biotec, Germany) and measured with a MACSQuant device. Flow cytometric data were analyzed with Flowlogic (Inivai Inc., Mentone Vicoria, Australia) and presented as mean ±SD. A one-way ANOVA and a Holm-Sidak test (multiple comparisons versus control group) were used to compare numeric data amongst the different experimental groups. p-values less than 0.05 were accepted as statistically significant. Flow cytometric measurements were carried out with cells from the inoculum and with microcarrier-free, purified cell samples from the different spinner flasks. For flow cytometric analysis, the cells were stained with fluorochrome-conjugated anti-human CD14, CD20, CD34, CD45, CD73, CD90, and CD105 (according to the recommendations of the International Society for Cellular Therapy ISCT and the International Federation for Adipose Therapeutics and Science IFATS) antibodies (Miltenyi Biotec, Germany) and measured with a MACSQuant device. Flow cytometric data were analyzed with Flowlogic (Inivai Inc., Mentone Vicoria, Australia) and presented as mean ±SD. A one-way ANOVA and a Holm-Sidak test (multiple comparisons versus control group) were used to compare numeric data amongst the different experimental groups. p-values less than 0.05 were accepted as statistically significant.

Spinner Flask Cultivations
For each condition, three spinner flasks (n = 3) of both spinner flask types were inoculated with hTERT-ASCs and cultivated for 10 days at 37 • C, 5% CO 2 , and 80% humidity. All spinner flasks were inoculated with the same inoculum (P24, PDL 35), which was prepared using cryopreserved hTERT-ASCs (one passage post thawing in T 75 flasks; 5000 cells/cm 2 ). In addition, the cells were also expanded in a T75-flask in parallel as a static control. Prior to inoculation, the MC suspension was equilibrated for 12 h. Post inoculation, no agitation was performed for 4 h in order to support cell attachment on the MC surface. After the cell attachment phase, the impeller speed was set according to the individual experimental conditions (SP100 = 25 rpm, 49 rpm N s1u , 60 rpm N s1 , 120 rpm; SP300 = 20 rpm, 41 rpm N s1u , 52 rpm N s1 , 100 rpm). On days 4 and 8, partial medium exchanges (50%) were performed. For this purpose, the impellers were switched off and the MCs were allowed to settle. Fifty percent of the working volume of the spinner flasks was replaced with fresh preheated medium, and the impellers were restarted. No MC feeds were performed during the cultivations.

Suspension Studies
Different MC suspension states are represented in Figure 6a for the Corning ® SP100 with a MC concentration of 10 g/L. As can be seen from the images, different suspension states can be detected: (1) transport of the MCs to the vessel center and the formation of a clear outer zone, (2) swirling up of the MCs from the center of the vessel bottom and further reduction of the MC solid fraction at the reactor bottom, and (3) maintaining the MCs in suspension (Ns1u < N < Ns1). Comparable suspension states were also observed for the SP300 spinner flask. This was not surprising, since the two spinner flasks have comparable geometrical ratios and are equipped with identical impellers: a large paddle-like impeller (d/D 0.58-0.68). Due to the shape of the paddle-like impeller and the absence of probes, which probably act as baffles, a mainly tangentially oriented fluid flow was induced. The secondary flow resulted in the transport of the MCs to the vessel center (see Section 3.2). From this area, the MCs were swirled up as the impeller speed was further increased.
The impeller speeds required to fulfill the two suspension criteria for different MC concentrations (2.5-20 g/L) were in the range of 30 rpm to 81 rpm (Ns1u) and 35 rpm to 90 rpm (Ns1) for the SP100 and 22 rpm to 65 rpm (Ns1u) and 29 rpm to 75 rpm (Ns1) for the SP300. Therefore, it can be concluded that Ns1 = (1.1-1.3)·Ns1u. Within the investigated MC concentration range, a linear correlation was found between the impeller speeds required to achieve the suspension criteria for the corresponding MC concentration. Thus, an average increase of 27 rpm per g MCs for the SP100 and an average increase of 29 rpm per g MCs for the SP300 are required to maintain the suspension state defined by Ns1u and Ns1. The impeller speeds required to ensure the suspension criteria in the SP100 were on average 25.25 ± 6.5% higher than for the SP300. However, when comparing the corresponding tip speeds (u tip = π n d), similar values for the SP100 and SP300 were calculated (see Figure 6b). This can be explained by the comparable geometrical ratios of the two systems and the resulting fluid flow conditions, which are described in Section 3.2. The maximum deviations of the measured data from the predicted values were between 8% and 10%, depending on the linear regression analysis for Ns1u (u tip,Ns1u = 0.0067·c MC + 0.0379) and Ns1 (u tip,Ns1 = 0.0069·c MC + 0.0581). However, these statistical correlations are only valid within the tested MC concentration range (2.5-20 g/L). The accuracy of the two regression models is also shown in Figure 6c, where the experimentally measured values were plotted against the predicted model values. It can be clearly seen that the values only slightly scatter around the diagonal center line, and the maximum deviations were between ±10% and ±12%. Nevertheless, the results demonstrate the linear relationship between the tip speed and the MC concentration as well as the applicability of the regression models for estimating the suspension criteria within the tested MC concentration range for the two spinner flasks. In fact, the dependency of the suspension criteria on the geometrical dimensions is undeniable. The determined impeller speeds and regression models served as a basis for the CFD simulations (see Sections 3.2-3.5).
correlation was found between the impeller speeds required to achieve the suspension criteria for the corresponding MC concentration. Thus, an average increase of 27 rpm per g MCs for the SP100 and an average increase of 29 rpm per g MCs for the SP300 are required to maintain the suspension state defined by Ns1u and Ns1. The impeller speeds required to ensure the suspension criteria in the SP100 were on average 25.25 ± 6.5% higher than for the SP300. However, when comparing the corresponding tip speeds (utip = n d), similar values for the SP100 and SP300 were calculated (see Figure 6b). This can be explained by the comparable geometrical ratios of the two systems and the resulting fluid flow conditions, which are described in Section 3.2. The maximum deviations of the measured data from the predicted values were between 8% and 10%, depending on the linear regression analysis for Ns1u (utip,Ns1u = 0.0067•cMC + 0.0379) and Ns1 (utip,Ns1 = 0.0069•cMC + 0.0581). However, these statistical correlations are only valid within the tested MC concentration range (2.5-20 g/L). The accuracy of the two regression models is also shown in Figure 6c, where the experimentally measured values were plotted against the predicted model values. It can be clearly seen that the values only slightly scatter around the diagonal center line, and the maximum deviations were between ±10% and ±12%. Nevertheless, the results demonstrate the linear relationship between the tip speed and the MC concentration as well as the applicability of the regression models for estimating the suspension criteria within the tested MC concentration range for the two spinner flasks. In fact, the dependency of the suspension criteria on the geometrical dimensions is undeniable. The determined impeller speeds and regression models served as a basis for the CFD simulations (see Section 3.2-3.5).

Single-Phase Fluid Flow Pattern
The main objective of the single-phase CFD simulations was to characterize and compare the fluid flow under the same conditions investigated in the cultivation studies (see Section 3.5). As shown in Figure 7, the fluid flow profiles in the two spinner flask types were similar, due to their comparable geometrical ratios. In both cases, the highest fluid velocities occurred at the edges of the impeller blades and corresponded quite well to the theoretical tip speed. An area with relatively weak fluid velocities was generated directly below the impeller (r/R ± 0.3) in both systems. Thus, this area represents a critical zone for MC sedimentation. The observed MC transport form the outer part of the vessel to the vessel center was mainly driven by the induced secondary flow, previously discussed in Section 3.1. Similar findings were also reported by Berry et al. [31], Liovic et al. [28], and Venkat et al. [43] in other types of small scale spinner flasks.

Single-Phase Fluid Flow Pattern
The main objective of the single-phase CFD simulations was to characterize and compare the fluid flow under the same conditions investigated in the cultivation studies (see Section 3.5). As shown in Figure 7, the fluid flow profiles in the two spinner flask types were similar, due to their comparable geometrical ratios. In both cases, the highest fluid velocities occurred at the edges of the impeller blades and corresponded quite well to the theoretical tip speed. An area with relatively weak fluid velocities was generated directly below the impeller (r/R ± 0.3) in both systems. Thus, this area represents a critical zone for MC sedimentation. The observed MC transport form the outer part of the vessel to the vessel center was mainly driven by the induced secondary flow, previously discussed in Section 3.1. Similar findings were also reported by Berry et al. [31], Liovic et al. [28], and Venkat et al. [43] in other types of small scale spinner flasks.
A more quantitative comparison of the volume-weighted fractions of the individual velocity components (x, y, z-directions) is shown in Figure 8a. The profiles of the velocity components are very similar for both systems. The results indicate that the fluid flow in both systems was mainly tangentially oriented. As defined by the boundary conditions, maximum tangential fluid velocity of 0.99 u tip (SP100) and 1.0 u tip (SP300) were expected. It can be concluded that the velocity distribution is relatively homogeneous without large fluctuations in the volume-weighted velocity profile. The axial part of the fluid flow was not particularly pronounced and the maximum values were between 0.41 u tip (SP100) and 0.59 u tip (SP300). The highest axial fluid velocities occurred directly below the impeller blade and were crucial for MC dispersion. The comparison of the dimensionless velocity magnitude (u/u tip ) along the dimensionless radial coordinates (r/R) at c/2 (=4 mm) indicated a greater velocity increase and a higher maximum fluid velocity in the SP300, which can be ascribed to the wider impeller blade (≈20 mm) (Figure 8b). However, the spatial dimensions of the critical fluid zone have similar courses in both spinner flasks.
(b) (c) Figure 6. Microcarrier suspension dynamics. (a) Photographic pictures of the MC distribution during the suspension studies (e.g., SP100, 10 g/L). (b) Graphical representation of the required tip speeds to achieve the suspension criteria at different MC concentrations for the SP100 and SP300. (c) Experimental vs. predicted Ns1u and Ns1 (=shown as tip speed) for the Corning spinner flasks (merged data for SP100 and SP300).

Single-Phase Fluid Flow Pattern
The main objective of the single-phase CFD simulations was to characterize and compare the fluid flow under the same conditions investigated in the cultivation studies (see Section 3.5). As shown in Figure 7, the fluid flow profiles in the two spinner flask types were similar, due to their comparable geometrical ratios. In both cases, the highest fluid velocities occurred at the edges of the impeller blades and corresponded quite well to the theoretical tip speed. An area with relatively weak fluid velocities was generated directly below the impeller (r/R ± 0.3) in both systems. Thus, this area represents a critical zone for MC sedimentation. The observed MC transport form the outer part of the vessel to the vessel center was mainly driven by the induced secondary flow, previously discussed in Section 3.1. Similar findings were also reported by Berry et al. [31], Liovic et al. [28], and Venkat et al. [43] in other types of small scale spinner flasks.  A more quantitative comparison of the volume-weighted fractions of the individual velocity components (x, y, z-directions) is shown in Figure 8a. The profiles of the velocity components are very similar for both systems. The results indicate that the fluid flow in both systems was mainly tangentially oriented. As defined by the boundary conditions, maximum tangential fluid velocity of 0.99 utip (SP100) and 1.0 utip (SP300) were expected. It can be concluded that the velocity distribution is relatively homogeneous without large fluctuations in the volume-weighted velocity profile. The axial part of the fluid flow was not particularly pronounced and the maximum values were between 0.41 utip (SP100) and 0.59 utip (SP300). The highest axial fluid velocities occurred directly below the impeller blade and were crucial for MC dispersion. The comparison of the dimensionless velocity magnitude (u/utip) along the dimensionless radial coordinates (r/R) at c/2 (=4 mm) indicated a greater velocity increase and a higher maximum fluid velocity in the SP300, which can be ascribed to the wider impeller blade (≈20 mm) (Figure 8b). However, the spatial dimensions of the critical fluid zone have similar courses in both spinner flasks.

Fluid Flow Field Verification
Since a number of mathematical assumptions were used for the CFD modeling of the two spinner flasks, stereoscopic PIV measurements were performed to validate the CFD-predicted fluid flow pattern (Figure 9). The contour plots show the fluid flow vectors in the x and y-directions as well as the fluid velocity component at a spatial position of 80° along mid plane of the vessel. A comparison of the fluid velocities in the SP100 was only possible for dimensionless radial coordinates between 0.50 and 0.82 (r/R), due to the pronounced curve of the vessel surface. Nevertheless, the qualitative comparison of the fluid velocity vectors showed good agreement between the CFD model and the PIV measurements. The only differences were slightly underestimated fluid velocities (0.76 utip) that were determined near the impeller bar. These differences can be accounted for by measurement uncertainties based on optical phenomena and the restricted measurement accuracy

Fluid Flow Field Verification
Since a number of mathematical assumptions were used for the CFD modeling of the two spinner flasks, stereoscopic PIV measurements were performed to validate the CFD-predicted fluid flow pattern (Figure 9). The contour plots show the fluid flow vectors in the x and y-directions as well as the w (up to 8.7%). However, the CFD velocity profiles were well captured and the overall agreement of PIV and CFD was satisfactory, findings are consistent with those published by Kaiser et al. [25]. All three velocity components in the SP300 were well captured by the stereoscopic PIV-measurements. The greatest differences (7.9-15%) were found for CFD was satisfactory, findings are consistent with those published by Kaiser et al. [25]. All three velocity components in the SP300 were well captured by the stereoscopic PIV-measurements. The greatest differences (7.9-15%) were found for between 0.7 and 0.85 (r/R). Hence, it can be postulated that, regardless of the differences, the established CFD model provides reliable fluid flow predictions in both spinner flask types.

Microcarrier Distribution
After determining and validating the fluid flow, the MC distribution was simulated (EE) for a MC solid fraction of 0.1% and for different impeller speeds, in order to compare the two systems. Figure 10 shows an example of the volume-weighted frequency distribution of the dimensionless MC solid fractions ( / mean) in the two spinner flasks for Ns1u. As expected, the highest MC volume

Microcarrier Distribution
After determining and validating the fluid flow, the MC distribution was simulated (EE) for a MC solid fraction of 0.1% and for different impeller speeds, in order to compare the two systems. Figure 10 shows an example of the volume-weighted frequency distribution of the dimensionless MC solid fractions (α/α mean ) in the two spinner flasks for Ns1u. As expected, the highest MC volume fractions (up to 2.8·α mean for the SP100 and SP300) were, in both cases, found directly below the impeller in the weak mixing zone. This observation was again not surprising because of the definition of Ns1u. Even for Ns1, a higher MC volume fraction was found in this region, because the MCs must periodically pass through the region in order to swirl up. The spatial position of the CFD-predicted MC deposits agreed well with the observations made during the suspension studies. The CFD-derived volume-weighted frequency distribution of the dimensionless MC volume fractions showed comparable MC homogeneity for the two spinner flask types (α/α mean close to 1). Even though the conditions were comparable in both systems at the vessel bottom for Ns1u and Ns1, this does not necessarily result in a same MC distribution over the entire vessel volume. The fronting of the distributions clearly indicates zones with low MC volume fractions. These zones were mainly determined near the fluid surface and represent the sedimentation boundary. The similar conditions at the vessel bottom can mainly be explained by the same off-bottom clearance (c = 0.07-0.12), whereas the MC distribution over the entire vessel volume is mostly affected by the d/D ratio. The swirl up of the MCs can be negatively affected by the recirculating fluid flow below the impeller. However, this effect is partially compensated in the two spinner flask types by the low off-bottom clearance. Due to the higher d/D ratio and the slightly different impeller shape of the SP100, recirculating fluid flow were formed near the fluid surface and decreased overall MC homogeneity (see Section 3.2). It seems that this slightly reduced overall MC homogeneity might have an effect on the MC-cell-aggregate formation rate due to the higher probability of particle interactions. However, this effect is partially compensated in the two spinner flask types by the low off-bottom clearance. Due to the higher d/D ratio and the slightly different impeller shape of the SP100, recirculating fluid flow were formed near the fluid surface and decreased overall MC homogeneity (see Section 3.2). It seems that this slightly reduced overall MC homogeneity might have an effect on the MC-cell-aggregate formation rate due to the higher probability of particle interactions. Based on the findings from the EE simulations, additional EL simulations were performed in order to derive the spatial distribution of discrete MC particles and to additionally calculate the circulation time, residence time and the hydrodynamic stresses acting on the particles (see Section 3.5). For this purpose, the two spinner flask types were vertically divided into four zones (∆ h/HL ≈ 0.25) and particle recovery (= percentage of particles per zone) was tracked for each zone. In order to start with "stationary" conditions, the tracking was started after particle movement had been simulated for at least ≥2• Θ % . Afterwards, particle recovery in each zone was calculated and averaged over a period of 1•Θ % . Figure 11 shows the individual particle recoveries calculated for the four spinner segments. The results revealed that the highest probability of presence of MCs is in the lowest spinner segment. This qualitative observation corresponds well with the results of the suspension studies and the EE simulations. The number of particles in the lower part of the SP100 was 2.8-21.5% higher than in the SP300 for all simulations. However, a clear reduction in the number of particles was observed at Based on the findings from the EE simulations, additional EL simulations were performed in order to derive the spatial distribution of discrete MC particles and to additionally calculate the circulation time, residence time and the hydrodynamic stresses acting on the particles (see Section 3.5). For this purpose, the two spinner flask types were vertically divided into four zones (∆ h/H L ≈ 0.25) and particle recovery (=percentage of particles per zone) was tracked for each zone. In order to start with "stationary" conditions, the tracking was started after particle movement had been simulated for at least ≥2·Θ 95% . Afterwards, particle recovery in each zone was calculated and averaged over a period of 1·Θ 95% . Figure 11 shows the individual particle recoveries calculated for the four spinner segments. The results revealed that the highest probability of presence of MCs is in the lowest spinner segment. This qualitative observation corresponds well with the results of the suspension studies and the EE simulations. The number of particles in the lower part of the SP100 was 2.8-21.5% higher than in the SP300 for all simulations. However, a clear reduction in the number of particles was observed at higher impeller speeds in both cases, resulting in increased particle recoveries in the other zones (zones 2-4). Interestingly, the particle reduction in zone 1 and the increase in the particle number in zones 2-4 were more significant at lower impeller speeds, showing an asymptotic convergence to absolute homogeneous conditions. Furthermore, the results indicate that the hydrodynamic stresses at the lower part of the spinner flasks in particular have the most significant effect on the cells, because the MCs are most often in this zone. However, the effects of the hydrodynamic stresses in the different zones depended heavily on the particle circulation and residence times (see Section 3.5), demonstrating the dynamics and complexity of the systems. To verify the applicability of the EL model, the CFD-predicted MC distribution was representatively verified for zones 3 and 4 of the SP100 using shadow imaging measurements (see Figure 12). For this purpose, a discrete number of MCs were added to the spinner flask and the positions of the MCs were captured and subsequently compared with the CFD-derived data. The interrogation space (≈2500-3000 mm 3 ) for the particle counting and comparison was defined based on the depth-of-field of the optical lens used for the measurements. The comparison of the modeled and measured particle data in zones 3 and 4 showed good overall agreement, even though the CFDderived particle values deviated slightly from the measured values (see Figure 12). The observed deviations could be explained by the measurement accuracy of the shadowgraphy method itself or by an overprediction of the turbulence parameters in the transient fluid flow regime in the turbulence models. The reduction in particle recovery in zone 4 from Ns1u (utip = 0.106 m/s) to Ns1 (utip = 0.130 m/s) can be explained by the fluid flow transition and the generation of a downward oriented recirculating fluid flow structure in the zone. However, the overall agreement was satisfactory and demonstrates the applicability of the EL model for the prediction of particle distribution in the SP100 and SP300. The measured particle velocity components and were also in a comparable range to velocities calculated for particles in the Langarian framework. Thus, no significant differences were found between cell-free and cell-loaded MCs. Therefore, it may be hypothesized that the MCs followed the fluid flow in a predominantly slip-free manner. This statement is also supported by the low sedimentation velocities (1.46 ± 0.56 mm/s) of the MCs. To verify the applicability of the EL model, the CFD-predicted MC distribution was representatively verified for zones 3 and 4 of the SP100 using shadow imaging measurements (see Figure 12). For this purpose, a discrete number of MCs were added to the spinner flask and the positions of the MCs were captured and subsequently compared with the CFD-derived data. The interrogation space (≈2500-3000 mm 3 ) for the particle counting and comparison was defined based on the depth-of-field of the optical lens used for the measurements. The comparison of the modeled and measured particle data in zones 3 and 4 showed good overall agreement, even though the CFD-derived particle values deviated slightly from the measured values (see Figure 12). The observed deviations could be explained by the measurement accuracy of the shadowgraphy method itself or by an overprediction of the turbulence parameters in the transient fluid flow regime in the turbulence models. The reduction in particle recovery in zone 4 from Ns1u (u tip = 0.106 m/s) to Ns1 (u tip = 0.130 m/s) can be explained by the fluid flow transition and the generation of a downward oriented recirculating fluid flow structure in the zone. However, the overall agreement was satisfactory and demonstrates the applicability of the EL model for the prediction of particle distribution in the SP100 and SP300. The measured particle velocity components → u and → v were also in a comparable range to velocities calculated for particles in the Langarian framework. Thus, no significant differences were found between cell-free and cell-loaded MCs. Therefore, it may be hypothesized that the MCs followed the fluid flow in a predominantly slip-free manner. This statement is also supported by the low sedimentation velocities (1.46 ± 0.56 mm/s) of the MCs.

Circulation Times, Residence Times, and Hydrodynamic Stresses
The circulation times and residence times were calculated for each individual spinner segment based on the particle tracking data (from EL simulations) and were subsequently averaged over the four segments (= mean circulation and residence times). Fully sedimented MCs, especially at lower impeller speeds (≤ Ns1u), were not considered for the analysis. Figure 13a shows the relationship between the mean circulation times and the mean residence times. As expected, the circulation times (2.7-11.5 s) decreased proportionally to the residence times (0.74-4.94 s) as the impeller speed was increased. Interestingly, the proportionality constants for the SP100 = 0.54 and the SP300 = 0.49 were quite similar. This observation can be ascribed to the comparable fluid flow conditions resulting from the comparable geometrical ratios of the two systems. The calculated mean particle forces (see Table  3), which were calculated based on the particle force balance during the simulations, were inversely proportional to the circulation and residence times (see Figure 13a), and are indicated by the size of the circles. This finding is not unexpected, since the specific power input, which was derived from the CFD simulations, increased by approximately the 3rd potency in both spinner flask types. However, an experimental verification of this correlation is difficult due to the low acting torques in the two spinner flasks. Interestingly, the mean values of the particle forces did not change significantly between the lower impeller speeds (<Ns1u) and the two suspension criteria, even though the circulation and residence times decreased by up to 50%. Impeller speeds exceeding Ns1u and Ns1 resulted in a slight decrease of the circulation times, although the related particle forces increased logarithmically to the resulting specific power input (see Table 3). The product of the impeller speed and the circulation time resulted in values between 6.0-8.6 (SP100) and 4.5-5.7 (SP300) for impeller speeds exceeding Ns1u and Ns1. Thus, no constant values were achieved for the investigated conditions. Comparable observations as for the specific power input are also possible when considering the local normal and shear stresses, which were calculated according to Wollny et al. [44] (see Figure 13b). The volume-weighted mean values of the local normal (SP100: 1.15-1.51•10 -3 N/m 2 , SP300: 0.69-0.88•10 −3 N/m 2 ) and shear stresses (SP100: 4.96-6.62•10 −3 N/m 2 , SP300: 4.00-4.98•10 −3 N/m 2 ) were in a comparable range in both spinner flask types for impeller speeds between Ns1u and Ns1. Thus, similar conditions in terms of hydrodynamic stresses can be expected for cultivations in the resulting specific power input range of 0.3-1.1 W/m 3 , which was derived from the simulations. The defined specific power inputs were comparable to the data described by Schirmaier et al. [16], Grein et al. [45], Cierpka et al. [46], and Lawson et al. [17] for benchtop and pilot scale bioreactors, who

Circulation Times, Residence Times, and Hydrodynamic Stresses
The circulation times and residence times were calculated for each individual spinner segment based on the particle tracking data (from EL simulations) and were subsequently averaged over the four segments (=mean circulation and residence times). Fully sedimented MCs, especially at lower impeller speeds (≤Ns1u), were not considered for the analysis. Figure 13a shows the relationship between the mean circulation times and the mean residence times. As expected, the circulation times (2.7-11.5 s) decreased proportionally to the residence times (0.74-4.94 s) as the impeller speed was increased. Interestingly, the proportionality constants for the SP100 = 0.54 and the SP300 = 0.49 were quite similar. This observation can be ascribed to the comparable fluid flow conditions resulting from the comparable geometrical ratios of the two systems. The calculated mean particle forces (see Table 3), which were calculated based on the particle force balance during the simulations, were inversely proportional to the circulation and residence times (see Figure 13a), and are indicated by the size of the circles. This finding is not unexpected, since the specific power input, which was derived from the CFD simulations, increased by approximately the 3rd potency in both spinner flask types. However, an experimental verification of this correlation is difficult due to the low acting torques in the two spinner flasks. Interestingly, the mean values of the particle forces did not change significantly between the lower impeller speeds (<Ns1u) and the two suspension criteria, even though the circulation and residence times decreased by up to 50%. Impeller speeds exceeding Ns1u and Ns1 resulted in a slight decrease of the circulation times, although the related particle forces increased logarithmically to the resulting specific power input (see Table 3). The product of the impeller speed and the circulation time resulted in values between 6.0-8.6 (SP100) and 4.5-5.7 (SP300) for impeller speeds exceeding Ns1u and Ns1. Thus, no constant values were achieved for the investigated conditions. Comparable observations as for the specific power input are also possible when considering the local normal and shear stresses, which were calculated according to Wollny et al. [44] (see Figure 13b). The volume-weighted mean values of the local normal (SP100: 1.15-1.51·10 −3 N/m 2 , SP300: 0.69-0.88·10 −3 N/m 2 ) and shear stresses (SP100: 4.96-6.62·10 −3 N/m 2 , SP300: 4.00-4.98·10 −3 N/m 2 ) were in a comparable range in both spinner flask types for impeller speeds between Ns1u and Ns1. Thus, similar conditions in terms of hydrodynamic stresses can be expected for cultivations in the resulting specific power input range of 0.3-1.1 W/m 3 , which was derived from the simulations. The defined specific power inputs were comparable to the data described by Schirmaier et al. [16], Grein et al. [45], Cierpka et al. [46], and Lawson et al. [17] for benchtop and pilot scale bioreactors, who postulated that specific power inputs of up to 2.1 W/m 3 are suitable for the MC-based expansion of hMSCs. However, time-dependent hydrodynamic stresses might differ from the two spinner flask types investigated in this study, because the benchtop and pilot scale bioreactors were equipped with axial conveying impellers. It is worth mentioning that the mean values of the local shear stresses in the SP100 increased more than in the SP300 for impeller speeds exceeding Ns1. This is mainly caused by the larger d/D ratio, which results in a higher level of turbulence. It is evident that the local shear stresses, which are suspected to cause higher cell damage [47], have a dominant effect over the normal stresses. Moreover, different studies in laminar flow bioreactors have demonstrated that the shear stress can affect the cell morphology and the formation of the extracellular matrix (ECM) [48][49][50]. Yeatts et al. reported that continuous shear stresses of up to 0.15 dynes/cm 2 (=0.015 Pa) can cause higher expression levels of osteoblastic markers such as osteopontin and osteocalcin [51,52]. However, in most of the cases the hMSCs were exposed to these continuous shear stresses (up to 12 dynes/cm 2 ) over a long period (up to 28 days). Thus, the effect of the changing hydrodynamic conditions and the short exposure times in stirred bioreactors needs to be investigated in subsequent studies.
Bioengineering 2018, 5, 106 19 of 30 postulated that specific power inputs of up to 2.1 W/m 3 are suitable for the MC-based expansion of hMSCs. However, time-dependent hydrodynamic stresses might differ from the two spinner flask types investigated in this study, because the benchtop and pilot scale bioreactors were equipped with axial conveying impellers. It is worth mentioning that the mean values of the local shear stresses in the SP100 increased more than in the SP300 for impeller speeds exceeding Ns1. This is mainly caused by the larger d/D ratio, which results in a higher level of turbulence. It is evident that the local shear stresses, which are suspected to cause higher cell damage [47], have a dominant effect over the normal stresses. Moreover, different studies in laminar flow bioreactors have demonstrated that the shear stress can affect the cell morphology and the formation of the extracellular matrix (ECM) [48][49][50]. Yeatts et al. reported that continuous shear stresses of up to 0.15 dynes/cm 2 (=0.015 Pa) can cause higher expression levels of osteoblastic markers such as osteopontin and osteocalcin [51,52]. However, in most of the cases the hMSCs were exposed to these continuous shear stresses (up to 12 dynes/cm 2 ) over a long period (up to 28 days). Thus, the effect of the changing hydrodynamic conditions and the short exposure times in stirred bioreactors needs to be investigated in subsequent studies. Another popular method for evaluating hydrodynamic stress is based on Kolmogorov's theory of isotropic turbulence [53][54][55]. While cells in suspension are assumed to only be affected by turbulent eddies of comparable size, those growing on the surface of a MC appear to be more shear sensitive. This might be because they are attached to relatively large particles that are more prone to collisions that might damage the cells. Croughan et al. [56] found that cell damage became significant when the smallest turbulent eddies were approximately two-thirds of the size of a MC. However, to apply Kolmogorov's theory, the fluid flow must be very turbulent. Taking into account the fact that Re < 10 4 (see Table 3), the fluid flow is in the transition region of Reynolds numbers, between laminar and fully turbulent conditions. Thus, it would be more reasonable to describe it as moderately turbulent in both cases [31,43]. However, the calculated maximum dissipation rates were higher by a factor of one or two in the impeller swept volume than in the bulk, and therefore, agreed well with findings from the literature [28,57]. As expected, the smallest turbulent eddies were found for the highest tested impeller speeds, with values between 30 μm and 47 μm. In terms of the suspension criteria, the minimum values were predicted in the range of 60 μm and 76 μm, which is much lower than the proposed 2/3 MC size. In contrast, the volume-weighted mean values were slightly higher than the MC size, which demonstrated that only a small proportion of the turbulent eddies are comparable in size to the MCs. This lowers the risk that the MCs might come into contact with these detrimental eddies. However, this fact depends heavily on the resulting circulation times and residence times of Another popular method for evaluating hydrodynamic stress is based on Kolmogorov's theory of isotropic turbulence [53][54][55]. While cells in suspension are assumed to only be affected by turbulent eddies of comparable size, those growing on the surface of a MC appear to be more shear sensitive. This might be because they are attached to relatively large particles that are more prone to collisions that might damage the cells. Croughan et al. [56] found that cell damage became significant when the smallest turbulent eddies were approximately two-thirds of the size of a MC. However, to apply Kolmogorov's theory, the fluid flow must be very turbulent. Taking into account the fact that Re < 10 4 (see Table 3), the fluid flow is in the transition region of Reynolds numbers, between laminar and fully turbulent conditions. Thus, it would be more reasonable to describe it as moderately turbulent in both cases [31,43]. However, the calculated maximum dissipation rates were higher by a factor of one or two in the impeller swept volume than in the bulk, and therefore, agreed well with findings from the literature [28,57]. As expected, the smallest turbulent eddies were found for the highest tested impeller speeds, with values between 30 µm and 47 µm. In terms of the suspension criteria, the minimum values were predicted in the range of 60 µm and 76 µm, which is much lower than the proposed 2/3 MC size. In contrast, the volume-weighted mean values were slightly higher than the MC size, which demonstrated that only a small proportion of the turbulent eddies are comparable in size to the MCs. This lowers the risk that the MCs might come into contact with these detrimental eddies. However, this fact depends heavily on the resulting circulation times and residence times of the MCs. The mean volume-weighted values for the highest tested impeller speeds were in both cases much closer to the detrimental theoretical value of 141 µm. Even though such eddies occurred at the suspension criteria, the frequency with which the MCs were exposed to such eddies was much lower due to the lower circulation times and residence times (see Section 3.6). Furthermore, in both cases, the volume of λ < 141 µm increased from 0.03 to 52.72% and from 0.02 to 63.26% as the impeller speed increased from 25 to 120 rpm and 20 to 100 rpm, respectively. Table 3. Overview of the main biochemical engineering parameters for the SP100 and SP300. Parameters were obtained from the CFD simulations.

Cultivation Studies
Figure 14a-d shows the time-dependent profiles of living cell densities and MC-cell aggregates for the SP100 and SP300. It can be seen that the investigated hydrodynamic stresses have a significant effect on cell growth and MC-cell-aggregate formation. The highest living cell densities were achieved, of up to 1.68 ± 0.36 × 10 5 cells/cm 2 (=6.25 ± 0.35 × 10 5 cells/mL, EF 56.01) and 2.46 ± 0.16 × 10 5 cells/cm 2 (=8.77 ± 0.66 × 10 5 cells/mL, EF 81.14), in the SP100 and SP300 when working at the suspension criteria. The living peak cell densities in the SP300 were on average up to 40% higher than those in the SP100. Although the two spinner flask types had comparable geometrical ratios, the hydrodynamic stresses in the SP100 were higher at the suspension criteria (see Section 3.5). In fact, the time-dependent stresses were higher due to the lower circulation times, which increases the risk that the cells on the MCs are more frequently exposed to detrimental stresses. At the same time, the residence times, and therefore, the exposure times of the MCs to the time-dependent stresses were shorter. This observation is supported by the slightly lower level of homogeneity in the SP100, as shown in the CFD-simulations. However, in both cases, the peak cell densities were in the same range as cell densities measured in planar, static cultures at maximum confluency (2.9 ± 0.09 × 10 5 cells/cm 2 , data not shown), in which the cells were expanded in parallel. This result indicates that the cells cultivated at the suspension criteria are mainly restricted by the available growth surface. In contrast, significantly lower cell densities were achieved at lower and higher impeller speeds. A peak living cell density of 1.05 ± 0.06 × 10 5 cells/cm 2 (=4.49 ± 0.06 × 10 5 cells/mL, EF 35.05) and 1.36 ± 0.57 × 10 5 cells/cm 2 (=4.88 ± 0.57 × 10 5 cells/mL, EF 45.20) was determined for the SP100 and the SP300 at 25 rpm and 20 rpm, respectively. These peak cell densities are up to 84% lower than those at the suspension criteria. This observation may be caused by the higher amount of sedimented MCs and the increased MC-cell aggregate formation (see Figure 14a). Although the specific power input for the same tip speed in the SP300 was slightly lower than for the SP100, shorter circulation times and residence times occurred and resulted in reduced MC-cell aggregate formation due to the higher grade of homogeneity. The amount of aggregates with a size of >1.0 mm increased significantly after three days of cultivation at low impeller speeds, which impairs cell growth. We observed that the LDH activity values calculated relative to the values obtained for the suspension criteria on day seven increased by between 32% and 44% in the supernatant. Even, from day 7 to day 10, the LDH activity further increased by up to 60% in all cultivations with impeller speeds between 20 and 60 rpm, which is accompanied by stagnant cell growth after day seven. Comparable results were also found for the cell viability, which was measured for the cells in the supernatant by flow cytometry. The viability of the cells on the MCs was always >99%. This is not surprising as dead cells detach from the MC surface. Thus, the increase of dead cells in the supernatant depends on the cell detachment from the MC surface and the die off of cells in the supernatant. Interestingly, the MC-cell aggregate formation had a stronger influence on the number of dead cells in the supernatant than the hydrodynamic stresses. The percentage of dead cells in the supernatant increased to 58% at the end of the cultivations (day 10) for N < Ns1u. In contrast, the percentage of dead cells in the supernatant for N > Ns1 was only 30%. This means that day seven represents the optimal point for cell harvesting. In contrast, no significant MC-cell aggregate formation was observed for higher impeller speeds, due to the higher hydrodynamic stress. Although a peak cell density of between 0.6 × 10 5 cells/cm 2 and 1.25 × 10 5 cells/cm 2 was achieved, the amount of larger aggregates was low, meaning that higher hydrodynamic stresses affected MC-cell aggregation. This was expected and has also been reported by Takahashi et al. [58] and Ferrari et al. [59], who both recommend using impeller speed, among other parameters, to control aggregate size to some extent. activity values calculated relative to the values obtained for the suspension criteria on day seven increased by between 32% and 44% in the supernatant. Even, from day 7 to day 10, the LDH activity further increased by up to 60% in all cultivations with impeller speeds between 20 and 60 rpm, which is accompanied by stagnant cell growth after day seven. Comparable results were also found for the cell viability, which was measured for the cells in the supernatant by flow cytometry. The viability of the cells on the MCs was always >99%. This is not surprising as dead cells detach from the MC surface. Thus, the increase of dead cells in the supernatant depends on the cell detachment from the MC surface and the die off of cells in the supernatant. Interestingly, the MC-cell aggregate formation had a stronger influence on the number of dead cells in the supernatant than the hydrodynamic stresses.
The percentage of dead cells in the supernatant increased to 58% at the end of the cultivations (day 10) for N < Ns1u. In contrast, the percentage of dead cells in the supernatant for N > Ns1 was only 30%. This means that day seven represents the optimal point for cell harvesting. In contrast, no significant MC-cell aggregate formation was observed for higher impeller speeds, due to the higher hydrodynamic stress. Although a peak cell density of between 0.6 × 10 5 cells/cm 2 and 1.25 × 10 5 cells/cm 2 was achieved, the amount of larger aggregates was low, meaning that higher hydrodynamic stresses affected MC-cell aggregation. This was expected and has also been reported by Takahashi et al. [58] and Ferrari et al. [59], who both recommend using impeller speed, among other parameters, to control aggregate size to some extent.  Table 4 shows the main growth-dependent parameters including the cell specific glucose consumption rate, the lactate production rate and the ammonia production rate. By considering the calculated specific glucose consumption rates, it becomes clear that the lowest values were obtained at the suspension criteria in both cases. This is due to the efficient metabolization of glucose under these conditions. The calculated values for the hTERT-ASCs agreed well with those determined by  Table 4 shows the main growth-dependent parameters including the cell specific glucose consumption rate, the lactate production rate and the ammonia production rate. By considering the calculated specific glucose consumption rates, it becomes clear that the lowest values were obtained at the suspension criteria in both cases. This is due to the efficient metabolization of glucose under these conditions. The calculated values for the hTERT-ASCs agreed well with those determined by Rafiq et al. [60] for hMSCs in different culture media. The highest specific glucose consumption rates (20.79-35.00 pmol/cell/d) were found at the highest impeller speeds. The relationship between the specific glucose consumption rate and the specific power input can be expressed by a logarithmic function of the 3 rd order (f (−q gluc ) = −11.085 + [−1.311 × ln(P/V)] + [−3.529 × ln(P/V)] 2 + [−0.806 × ln(P/V)] 3 , R 2 = 0.963), whereas this correlation is only valid for the investigated P/V range. Similar correlations were also calculated for the cell specific lactate and ammonia production rates, where values of up to 193% and 170% higher than those in the spinner flasks at the suspension criteria were determined at the highest impeller speeds. These higher values indicate that the cells are more stressed at higher impeller speeds as a result of higher hydrodynamic loads. Moreover, such high metabolite production rates also result in a large accumulation of inhibitory metabolites during cultivation and may reduce cell yield [39,61,62]. The different obtained correlations were used as initial parameters for the growth modelling (see Section 3.6).  Figure 15a shows the relationship between the overall mean specific growth rate and the specific power input. The parabolic curve profile of the specific growth rate shows optimal cell growth between Ns1u and Ns1. For specific power inputs between 0.33 and 1.12 W/m 3 , maximum values for µ between 0.70 and 0.74 d −1 (=0.93-0.99 d) were achieved. Comparable growth rates for hTERT-MSCs were also described by Balducci et al. [63], Leber et al. [64], and Cierpka et al. [46]. Moreover, the maximum specific growth rates correlate quite well with the values (0.70 ± 0.02 d −1 ) from experiments in the planar and static culture systems (data not shown). This demonstrates the comparability of the two spinner flask types, even though slight deviations exist. Similar relationships to those for the specific growth rate and the specific power input were also found for other hydrodynamic stress parameters (i.e., l λ , LSS, F; see Table 3). The derived correlations represent the basis for growth modeling and future scale-up investigations.
Based on the measured MC-cell aggregates, a cell specific aggregation rate was derived from the data. For this purpose, the data of the class II aggregate (>1.0 mm) was correlated with cell growth over the time. The dependency of the MC-cell aggregation rate on specific power input is shown in Figure 15b. The determined MC-cell aggregation rates for the SP100 were higher than for the SP300 for all investigated conditions. However, this was not surprising because of the lower MC homogeneity, especially at lower impeller speeds (N < Ns1u and Ns1). for all investigated conditions. However, this was not surprising because of the lower MC homogeneity, especially at lower impeller speeds (N < Ns1u and Ns1).
(a) (b) Figure 15. Dependency of the specific growth rate (a) and the MC-cell aggregation rate (b) on the specific volumetric power input. The grey marked area indicates the range between Ns1u and Ns1. Table 5 shows the results of the flow cytometric measurements after cell harvesting for the four negative (CD14 − , CD20 − , CD34 − , and CD45 − ) and the three positive (CD73 + , CD90 + , and CD105 + ) markers. In order to visualize a change in the marker expression profile, the results were compared with those obtained from the cell inoculum. All positive markers were strongly expressed (>91%) and the negative markers exhibited a lack of expression (<2.7%). These results correspond with measurements for hTERT-ASCs by Balducci et al. [63], Yin et al. [65], and Wolbank et al. [66]. Moreover, the results are in good agreement with the minimal expression levels for hASCs recommended by the International Society of Cellular Therapy (ISCT) and the International Federation for Adipose Therapeutics and Science (IFATS) [67]. Statistically significant differences were found in all spinner flask cultivations for CD73, CD90, and CD105, when they were compared with the expression levels obtained from the cell inoculum. However, a correlation of the surface Figure 15. Dependency of the specific growth rate (a) and the MC-cell aggregation rate (b) on the specific volumetric power input. The grey marked area indicates the range between Ns1u and Ns1. Table 5 shows the results of the flow cytometric measurements after cell harvesting for the four negative (CD14 − , CD20 − , CD34 − , and CD45 − ) and the three positive (CD73 + , CD90 + , and CD105 + ) markers. In order to visualize a change in the marker expression profile, the results were compared with those obtained from the cell inoculum. All positive markers were strongly expressed (>91%) and the negative markers exhibited a lack of expression (<2.7%). These results correspond with measurements for hTERT-ASCs by Balducci et al. [63], Yin et al. [65], and Wolbank et al. [66]. Moreover, the results are in good agreement with the minimal expression levels for hASCs recommended by the International Society of Cellular Therapy (ISCT) and the International Federation for Adipose Therapeutics and Science (IFATS) [67]. Statistically significant differences were found in all spinner flask cultivations for CD73, CD90, and CD105, when they were compared with the expression levels obtained from the cell inoculum. However, a correlation of the surface marker expression levels with the hydrodynamic stress or MC-cell aggregate formation was not found.

Growth Modelling
To test the validity of the unstructured, segregated, simplistic growth model for the SP100 and SP300, independent cultivation runs (n = 3 per spinner flask type) were performed at the Ns1u criterion. To simulate the cell density, the substrate and the metabolites, the parameters determined in the previous cultivations were used (see Section 3.5). Figure 16 shows the measured values and the simulated time courses for the cell density (a-b), the substrate and the metabolites (c-d). The simulated time courses show good overall agreement with the experimentally measured values and demonstrate the applicability of the unstructured, segregated growth model. By using the determined growth parameters from the cultivation study, the cell growth, glucose consumption, lactate production, and ammonia production could be well approximated. The greatest deviations in cell density were in the range of 3-20% for the cells in suspension and 4-24% for the cells on the MCs. However, the accuracy of the measured cell densities decreased towards the end of the cultivations due to the formation of larger aggregates (see Section 3.5). The intensified aggregate formation and its effect on the growth surface was not considered in the model. However, the principal growth behaviors were well captured. Moreover, the cell densities measured in both systems were comparable to the previous cultivations (see Section 3.5), and therefore, demonstrated the reproducibility of the cultivation processes, especially when using the stable hTERT-ASCs cell line. Furthermore, the glucose, lactate and ammonia time courses, agreed well, even though the determined specific substrate consumption and metabolite production rates were prone to errors due to the accuracy of the measurement method. The lactate and ammonia time courses show that the growth associated assumption was valid and can be used for modelling because the production of these two substances was directly linked to the cell number increase.

Conclusions
In this study, extensive numerical and experimental investigations were carried out for two spinner flask types with comparable geometrical ratios. For this purpose, comprehensive multi-phase CFD simulations based on EE and EL models were performed and the derived hydrodynamic stress parameters were linked with growth-related parameters to create an unstructured, segregated, simplistic growth model. Providing the numerical models are valid, CFD modeling is one of the most effective techniques for characterization of flow fields. However, multi-phase CFD simulations, which take the MCs into account, have not been used very often for the prediction of MC-related stress parameters. Most of the previous CFD-based investigations only focused on single-phase flows, which might be driven by the higher computational requirements for multi-phase simulations and their economical relevance. However, it is widely accepted that hMSCs are exposed to different hydrodynamic stress levels during their lifespan in a stirred bioreactor. Thus, particle-related data, like the MC circulation and residence times in high shear zones, among others, are meaningful and will become more important in the future. They allow the bioreactor design to be improved and the optimum growth parameters to yield relevant amounts of therapeutic cells to be determined.
In order to test the influence of different hydrodynamic stress levels on cell growth and cell quality, CFD simulations were performed for different impeller speeds and for the suspension criteria (Ns1u, Ns1), which were determined experimentally. The suspension studies clearly showed the linear relationship between MC concentration and the tip speed required for Ns1u/Ns1. The resulting correlation from the multi-regression analysis can be used for further studies with Corning ® spinner flasks using either the same or other polystyrene-based MCs (with comparable properties). Due to the comparable geometrical ratios of the two spinner flask types used in this study, suspension

Conclusions
In this study, extensive numerical and experimental investigations were carried out for two spinner flask types with comparable geometrical ratios. For this purpose, comprehensive multi-phase CFD simulations based on EE and EL models were performed and the derived hydrodynamic stress parameters were linked with growth-related parameters to create an unstructured, segregated, simplistic growth model. Providing the numerical models are valid, CFD modeling is one of the most effective techniques for characterization of flow fields. However, multi-phase CFD simulations, which take the MCs into account, have not been used very often for the prediction of MC-related stress parameters. Most of the previous CFD-based investigations only focused on single-phase flows, which might be driven by the higher computational requirements for multi-phase simulations and their economical relevance. However, it is widely accepted that hMSCs are exposed to different hydrodynamic stress levels during their lifespan in a stirred bioreactor. Thus, particle-related data, like the MC circulation and residence times in high shear zones, among others, are meaningful and will become more important in the future. They allow the bioreactor design to be improved and the optimum growth parameters to yield relevant amounts of therapeutic cells to be determined.
In order to test the influence of different hydrodynamic stress levels on cell growth and cell quality, CFD simulations were performed for different impeller speeds and for the suspension criteria (Ns1u, Ns1), which were determined experimentally. The suspension studies clearly showed the linear relationship between MC concentration and the tip speed required for Ns1u/Ns1. The resulting correlation from the multi-regression analysis can be used for further studies with Corning ® spinner flasks using either the same or other polystyrene-based MCs (with comparable properties). Due to the comparable geometrical ratios of the two spinner flask types used in this study, suspension criteria were fulfilled at the same tip speed. However, the CFD simulations indicated that this does not necessarily result in the same overall MC homogeneity in the two systems, because the criteria only consider conditions at the bottom of the reactor. The slightly reduced overall MC homogeneity at Ns1u and Ns1 in the SP100 not only affects cell growth but also aggregate formation to some extent, even though the two systems had comparable geometrical dimensions. The observed deviations in cell growth and aggregate formation were mainly due to slight differences in the d/D and c/D ratios. Hence, the investigations demonstrated the dependency of the suspension criteria on the geometrical dimensions, findings that can serve as a basis for further scale-up studies. Nevertheless, optimal cultivation conditions were found for Ns1u < N < Ns1, which corresponded to specific power inputs between 0.3 and 1.1 W/m 3 . Although, all numerical methods, including CFD, use mathematical assumptions and empirical variables, the experimental PIV and shadowgraphy measurements demonstrated their accuracy and usefulness. In all cases, there was sufficient and satisfactory comparability between the experimental and the modeled data, underlining the meaningfulness of the derived hydrodynamic parameters (e.g., P/V, τ nt , F MC , etc.). Based on these parameters and their correlation with the growth-related data, stress dependent correlations (e.g., µ vs. P/V or τ nt , −q i vs. P/V or τ nt , etc.) were, for the first time, derived for hTERT-ASCs. The correlations were used to set-up an unstructured, segregated growth model with very well matching time courses compared to the experimental data for cell growth, glucose consumption, lactate production and ammonia production. Maximum deviations between the simulated and measured cell densities for cells growing on the MC surface were in the range of 4 to 24%. This means that the descriptiveness and predictive power of the model is satisfactory, especially when considering the accuracy of the experimentally measured values. The intensified aggregate formation was not considered in the growth model. Nevertheless, good agreement has already been achieved. Subsequent investigations are necessary to understand the aggregation dynamics under different conditions in more detail and to consider their effects in the growth model. In addition, the applicability of the derived stress-dependent correlations for hTERT-ASCs in other systems and their transferability to primary hMSCs should be studied in subsequent work. For primary hMSCs, physiological states such as the replicative senescence should be considered in the growth model. This might be done by incorporating a senescence-related inhibition factor or by using population modelling approaches (e.g., proliferating population vs. senescent cell population). However, further investigations are necessary to determine such senescence-related factors, which can then be incorporated in a mathematical model.