1. IntroductionBone tissue engineering (BTE) is an interdisciplinary field in which physicians, biologists and engineers collaborate to develop and produce artificial bone tissues [1]. One common strategy for the treatment of bone defects are in vitro applications of cultured autologous cells or stem cells on scaffolds followed by the implantation into specific tissues and tissue regions. Thereby, scaffolds play key roles, serving to mimic physiological in situ environments [2]. For this reason, they should have high interconnectivity, high porosity, good biocompatibility, as well as sufficient mechanical strength and reliability [3]. Besides metals and ceramics [4], plastic polymers such as Polycaprolactone (PCL) are in use as materials for the production of these three-dimensional (3D) structures [5].PCL is a 3D-printable semi crystalline synthetic polymer that can be reversibly deformed in certain temperature ranges [6,7,8]. Thereby, the melting point reaches from 55 to 63 °C [9]. PCL shows viscoelastic behavior [10] and a good biocompatibility [11]. Pure PCL is hydrophobic, which can negatively affect cell adhesion and cell proliferation [12,13]. Scaffolds can be manufactured using additive producing processes, which are also coined as 3D printing. Hereby, products are created by adding the material step by step. Thereby, the geometry of the scaffolds is modeled using computer-aided design (CAD) software [14]. Single layers can be produced by heated viscous materials that are pressed out of a needle by compressed air. The final construct geometry is created by stacking the PCL layers “layer by layer” [15]. A widely used method allowing the approximation of mechanical and structural behaviors via numerically solving partial differential equations in engineering and BTE is the finite element method (FEM). Hereby, in contrast to classical numerical approaches, the function is not defined over the entire functional area but is divided into a finite number of subareas. With a parametric FEM design, geometric and material properties can be varied [16,17]. In this study, we aimed to develop a parametric model for layered, 3D-printed PCL scaffolds, consisting of interconnected compartments. This kind of scaffold allows the development of 3D networks of cells. As architecture plays a crucial role for the successful development and usage of these layered scaffolds, a suitable and flexible numerical in silico twin can be an important tool [18,19,20,21]. Currently, there are no parametric numerical models available for such scaffolds. 2. Materials and Methods 2.1. Fabrication of Scaffolds

PCL scaffolds were produced with a 3D-Bioplotter-System (ET.B4.16.11.M092, Envisiontec, Gladbeck, Germany) using PCL (Mn = 45,000, Sigma-Aldrich, St. Louis, MA, USA) and 24 G stainless steel needle tops (Envisiontec, Gladbeck, Germany). 3D printings were performed with process temperatures of 80 °C, platform temperatures of 18 °C, pressures of 4 bar, pre- and post-flow times of 0.4 s, needle offsets of 0.2 mm and speed of 0.5 mm/s.

Collagen coating of scaffolds was performed as previously described [12,22]. In detail, as a first step, scaffolds were washed with ethanol (30, 50, 70, 80, 96 and 100%, for five seconds each). To improve cell–scaffold interaction, polar oxygen groups were created on PCL surfaces, using plasma devices (PiezoBrush PZ2, Relyon Plasma, Regensburg, Germany) [23,24]. Following this, collagen coating was performed in 24-well plates (Costar, Merck, Darmstadt, Germany) in collagen solutions (concentration: 3.82 mg/mL; Rat Tail Collagen Type I, Corning, Corning, NY, USA) for 72 h at 5 °C. Next, scaffolds were dried for 24 h at 37 °C. Crosslinking of collagen was performed using a 95% ethanol, 50 mM 1-ethyl-3-carbodiimid in water solution. 2.2. Cell SeedingMG-63 (MG63 ATCC CRL-1427, ATCC, Manassas, VI,), human osteoblasts, MLO-Y4 cells (MLO-Y4 EKC002, Kerafast, Boston, MA, USA) and MSC were used for the present work. MG-63 cells were purchased from ATCC and MLO-Y4 cells from Kerafast. In addition, MSC from human bone marrow and human osteoblasts from tibia resections (total knee arthroplasty) were used. Ethical approval was granted as ethics vote “FREEZE” by the Ethics Committee of the University Medical Center Freiburg. In a previous work, we observed different cell growth on collagen-I-coated PCL scaffolds for different cells [25].Scaffolds were seeded with MG63 cells (MG63 ATCC CRL-1427, ATCC), MLO-Y4 cells (MLO-Y4 EKC002, Kerafast), mesenchymal stem cells (MSC) and osteoblasts in 24-well plates. MSC and osteoblasts were gained during operative hip surgery from proximal femur tissue. Second passage of cells was used. Cells were cultured in 75 cm2 flasks (Cellstar, Greiner Bio-One, Frickenhausen, Germany) at 37 °C with 5% CO2 in specific nutrient solutions (see Supplementary Table S1). Osteoblasts, MG63 and MLO-Y4 cells were digested using 3 mL Trypsin-EDTA (Merck, Darmstadt, Germany), MSC using 5 mL Accutase (Sigma-Aldrich). Following this, cells were counted in a Neubauer counting chamber and diluted with nutrient solutions to the aimed cell concentrations (Osteoblasts, MG63 and MLO-Y4 cells: 50.000 cells/100 μL; MSC: 50.000 cells/100 μL). Subsequently 100 μL of the cell-containing dilution was dribbled on the scaffolds and incubated at 37 °C with 5% CO2 for 2 h. Following this, 1 mL of the cell-containing dilution was added to the wells, and specimens were further incubated at 37 °C with 5% CO2. Nutrient solutions were exchanged every 3 days. PCL scaffolds were stained with 1:10 diluted Giemsa staining and Live-Dead-Assays (Live/Dead Cell Staining Kit II, PromoCell, Heidelberg, Germany). Images were obtained using microscopes (Olympus BX51, KEYNECE Germany, Neu-Isenburg, Germany) with associated fluorescent lamps (X-Cities Series 1207, KEYNECE Germany, Neu-Isenburg, Germany). Collagen coating was assessed using electron microscopes (ESEM FEI Quanta 250 FEG; FEI, Hilsboro, TX, USA). 2.3. Characterization of ScaffoldsGeometrical properties of scaffolds were determined with a Keyence VK-X200 Series laser microscope (KEYENCE, Tokyo, Japan). The assessment approach is given in Figure 1.

Material properties of PCL were determined via tensile testing according to DIN EN-ISO-527-1 with a Z005-testing machine (Zwick/Roell, Zwick/Roell GmbH & Co. KG, Ulm, Germany) standardized test body shape: 5A, manufacturing temperature: 80 °C, test temperature: 23 °C, Testing speed: 10 mm/min). Further creep testing of PCL was performed with an MCR 702 MultiDrive system (Anton Paar GmbH, Graz, Austria, solid rectangular fixture 5, creep stress = 1 N/mm2, nitrogen atmosphere, manufacturing temperature: 80 °C, test temperature: 23 °C, test body shape: 35 mm × 14 mm × 2.5 mm).

Tensile testing of scaffolds was executed within 24 h after production with Zwick and Bose testing machines (BPS-LH0010.21.00, Bose, Friedrichsdorf, Germany) test temperature: 23 °C, Testing speed: 10 mm/min).

2.4. Finite Element MethodParametric FEM simulations were performed using ANSYS Mechanical APDL R1. The following assumptions for the numerical models were made: struts are idealized as beam elements; cross-sections are assumed to be ideally round and were modelled as solid cross-sections with variable diameters; material properties are constant and isotropic; a consistent distance between struts and rigid connections between struts are assumed. Table 1 gives the geometric parameters of the numerical models:FEM models of PCL scaffolds were implemented based on three-dimensional 2 nodes Beam188 Timoshenko elements, which can model normal and shear forces and bending moments. In contrast to Bernoulli elements, no shear rigidity of the cross-sections is assumed. Thus, stocky geometries can also be analyzed [26]. The size of elements were varied via the variable “element_size”. Parameterized geometries of numerical models were defined via key points and generated via lines. The basic geometric development is illustrated in Supplementary Figure S1. Elastic-plastic and viscoelastic material models based on a Prony-Row were implemented [27,28]. Material properties were implemented via the variables “emodulus_pcl”,”f_y”, “prony_value” and “prony_time”. Loads were applied as path-controlled node loads. The basic structure of a numerical model and the comparison of the scaffold and model node areas are demonstrated in Figure 2.For validation of the numerical model, the following force ratios between the results from the tensile testing of scaffolds and numerically determined values were defined (ε = strain, n = number of tensile tests). Equation (1) is defined as the ratio between the determined force during the tensile testing and the numerical calculated force at a certain strain. Equation (2) describes the average of Equation (1) for 1%, 2% and 3% strains. Ratio (3) is the ratio of the average of all tensile testing and the numerical determined force.

μi=FNumerical,ε=iFTensileTesting,ε=i

(1)

μsum=μ0.01+μ0.02+μ0.033

(2)

μ¯=FNumerical(∑i=1nFTensileTesting,i)/n

(3)

2.5. Statistical Analysis

Statistical analyses were performed using IRIS RheoHub 2018, Origin 2019, Excel 2017 and SPSS 26. All values in this paper are expressed as mean ± standard deviation.

3. ResultsThe results of the geometrical characterization of ten scaffolds printed on different days are shown as boxplots in Figure 3, demonstrating mean pore sizes of about 800 ± 52 μm in the xy plane and 197 ± 20 μm in the z direction. The mean thickness of the struts was 143 ± 26 μm and the nodes areas showed a mean diameter of 338 ± 43 μm.Overall, six PCL samples were tensile tested. The mean Young’s modulus of the PCL was 566 ± 45 N/mm2, the mean yield strength was determined to be 17.9 ± 1.0 N/mm2 and the mean strain at yield strength was 8.8 + 2.1%. Associated stress–strain diagrams are provided in Figure 4.The resulting Prony terms of the creep testing with a duration of 27 h are provided in Table 2.PCL scaffolds show an almost linear-elastic behavior up to tensile strains of 3%. From an elongation of 5% stress, the data show a plastic behavior. The mechanical behavior of PCL scaffolds is shown in Figure 5.Numerical calculations based on Timoshenko Beam188 elements (element sizes: a_xy/16) show valid results with a force ratio of μ¯sum = 1.0515 (μ¯0.01 = 1.0049, μ¯0.02 = 1.0412, μ¯0.03 = 1.1083). Thereby, the numerical results show stable results, independent of the implemented element types (rods vs. beams), grade of element functions (Beam188 vs. Beam189), size of elements and material laws (viscoelastic vs. elastic-plastic). The related data are provided in Supplementary Table S2.The parameters of numerical models can be varied. Figure 6 shows exemplary effects of strut thicknesses on numerical results.Studies with electronic microscopes show continuous collagen cover layers on all struts of PCL scaffolds. This is illustrated in Figure 7, showing, on the left side, a PCL scaffold before collagen coating and, on the right side, a PCL scaffold after coating. Without treatment with plasma devices, no collagen coating could be achieved due to less adherence on the PCL surface.MG63 cells, MLO-Y4 cells and MSC were marked via Giemsa staining on all PCL scaffold struts after three days of incubation, illustrated in Figure 8. MG63 cells exhibited the highest cell density on scaffolds. In preliminary tests, no cell growth was observed on specimens not coated with collagen.The growth of osteoblasts after three, seven and ten days is demonstrated in Figure 9. The Live/Dead (L/D) assays show vital osteoblasts in both strut and nodal areas of the PCL scaffolds. Further images of L/D assays of osteoblast-seeded scaffolds are provided in Supplementary Figure S3. 4. Discussion

The principal new finding of this study is the development of parameter-based numerical models for layered 3D-printed PCL scaffolds consisting of interconnected compartments. These in silico digital twins allow a fast, flexible and targeted scaffold design and production while conserving resources.

3D FEM models have been developed to describe the mechanical behavior of bone, cartilage and of scaffolds in BTE applications [29]. This is usually followed by creating 3D computer-aided design (CAD) models for 3D printing processes. The resulting datasets are commonly used as geometrical information for numerical FE simulations during preprocessing. Reversing that classical workflow, the use of parametric geometrical approaches allows us to widely vary these geometries to identify optimal structural shapes prior to fabrication. Such parametric geometrical numerical models are used in BTE applications, e.g., for describing human femur models [30], orthopedic implants [31] and the design of scaffolds [32]. Importantly, no parametric geometrical modeling studies for layered scaffolds based on PCL are available. So far, only a few FEM models for PCL have been reported. Those either focus on the influence of the unit cell, porosity and needle diameter [20] on the mechanical behavior, or the addition of calcium phosphates [33] or nanohydroxyapatite to improve biocompatibility [34]. In these studies, the compressive strength of the entire layered scaffold was investigated. Here, we fabricated layered PCL scaffolds with interconnected compartments and developed for the first time geometrical parametric numerical models that determined the mechanical behavior of such scaffolds. The rationale for developing layered PCL scaffolds with interconnected compartments is their flexible properties and applications, despite their simple manufacturing method. Compartments can be filled with growth factors, hydrogels, fluids or drugs and enable three dimensionally connected cell growths [18,19,20,21].FEM is an established numerical approach that is used to determine forces, stresses and stains [16,17]. Due to its compact geometry, the developed parameterized numerical models were implemented based on Beam188 Timoshenko elements. Comparative calculations with Strut-Elements and Euler–Bernoulli beams do not show considerable differences under the chosen boundary conditions. This is due to the fact that the length-to-width ratio is at a border between the use of Timoshenko and Euler–Bernoulli theory. Furthermore, the pure tensile load of elements under the present boundary conditions leads to the classical tensile strut behavior of the model. The developed numerical model allows for the flexibility study of different geometries in silico to specifically plan the load testing of PCL scaffolds. This is particularly relevant for studies with osteoblasts, as they react specifically to mechanical stimulations [35,36].Furthermore, we demonstrated the successful development and additive fabrication of layered scaffolds based on PCL [11,37,38], but with specific architectures that included interconnected compartments. In a previous study, we demonstrated that variations in the architectural design of tissue-engineered constructs, combined with choosing specific biomaterial properties, can strongly influence the biomechanical performance of constructs. For example, the failure load of lamellar constructs can be increased by altering the lamellar width, number and order, and the increase reached such an extent that the failure loads of the individual materials used for construction were greatly exceeded [39]. In the present study, the architectures had 3D interconnected pores and slim strut thicknesses of about 150 μm. The size of the pores and 3D connectivity play important roles to ensure sufficient nutrient transportation after in situ implantation [40,41]. Our data show pore sizes of about 800 μm and 200 μm. The recommended ranges of pore size for bone tissue regeneration are 150–900 μm [42,43,44,45,46]. Thus, we assumed the developed pore architecture to be suitable for BTE applications. As studies on PCL scaffolds have shown the evocation of osteogenesis in vivo after the seeding of osteoblasts [47] and the successful use in BTE applications [48,49], we expect comparable behavior of the scaffolds of the present study.PCL is a semi crystalline thermoplastic that is widely used in BTE applications [6,50] and shows sufficient biocompatibility [11]. The polymer induces only minor immune responses after implantation and is biodegradable [51,52]. Thereby, degradation starts in situ after 4 to 6 months. Full loss of mass can take more than 48 months [53,54]. Furthermore, PCL shows positive effects on bone healing processes and the growth of fetal bovine chondrocytes [55,56]. On mechanical testing, PCL shows viscoelastic behavior [10]. Thus, the material properties depend on the testing velocities and testing temperature, and all testing in this study was performed under the same conditions. Our mechanical data show plastic behavior of the polymer, with a molar mass of 45.000 and a production temperature of 80 °C, in the range of the comparable literature [9]. Tensile testing of scaffolds shows a linear-elastic elastic-plastic behavior, as shown in other studies [22]. The behavior of PCL scaffolds under compression loads was studied by Huber et al. [57]. Due to the small specimen size of the scaffolds, a perfect axial load application cannot be guaranteed; therefore, the occurrence of bending moments cannot be excluded. Fabricated scaffolds were used for the validation of the parametric numerical models. During testing of the real scaffolds, eccentric load application is not fully preventable. Thus, forces and stresses beyond the loading axis can arise. These effects are not represented in the numerical model and lead to a potential underestimation of the internal forces in the numerical model compared to the testing. Considering these effects, basic variations of the real scaffolds’ geometries and material behaviors of PCL, the calculated force ratios demonstrate a sufficient numerical approximation of displacements and forces for strains up to 3% by the developed parametric numerical models. Thus, in this range, forces and strains can be validly determined for specific loads with the present numeral model.Due to the hydrophobic properties of PCL [13], and to improve cell–scaffold interaction [58], scaffolds were treated with plasma generators and were coated with rat collagen type I. Various comparative studies showed a positive effect of plasma treatment and collagen coating on cell seeding and growth [59,60]. In the present study, collagen coating was visible in electron microscopy images. Studies with immunoassays on comparable scaffolds with the same coating technique show successful coatings [25]. We therefore assume a sufficient collagen coating on the scaffolds in this study.To improve the stability and durability, the crosslinking of collagen via EDC was performed [61]. Seeded cells (Osteoblasts, MSC, MYL-Y4, MG63) are widely used in BTE applications and research [62,63,64,65,66,67]. Thus, they can be considered suitable for this study. All cells show growth on PCL strands after three days. Thereby, MG63 cells show pronounced colonization. As these cells are derived from osteosarcoma cells, a more rapid cell division could be expected [68]. The results from the Live–Dead Assay analysis underline the successful seeding of cells, as vital osteoblasts can be demonstrated after three, seven and ten days of growth. In the current study, we aimed to show a qualitative representation of growth. In a further study with the same material and coating, but different scaffold geometry, we were able to demonstrate cell growth. The data are provided in Supplementary Figure S2.Our study is limited, as the determined material and mechanical properties of the PCL and scaffolds are only valid for the chosen molar weight of the PCL, process temperature and testing speed. Furthermore, eccentric load application can lead to inaccuracies in results. The polymers were processed in this study at the same day as fabrication. Embrittlement due to UV radiation exposure was not investigated and cannot be ruled out. For numerical modeling, various assumptions are made, such as the presence of an isotropic and constant material. These inhomogeneities of scaffolds are not represented by the numerical approximation. The developed models provide sufficient approximation in linear-elastic ranges up to 3% strains. For calculations at higher strains, further investigations and adapted non-linear material models are necessary. Furthermore, the collagen layer and seeded cells were not modeled. Due to their small dimensions compared to the scaffolds, it can be assumed that the effects of the collagen layers and cells on the mechanical behavior are negligible. Nevertheless, an influence cannot be ruled out based on this study. This study focused on the parametric numerical modeling and shows only a qualitative growth of bone cells on the fabricated scaffolds. For a detailed description and analysis of the cell growth, further quantitative data are needed. Furthermore, calcium staining or the detection of calcium or mineral deposits on cell-seeded scaffolds are needed to prove the final suitability for bone tissue engineering applications. As Giemsa staining and no L/D assays for MSC, MG-63 and MLO-Y4 cells were performed, no statement can be made concerning the vitality of these cells. The present study did not explore the effects of varied pore size on the biological behavior of scaffolds in vivo. This was published previously for comparable scaffolds [25,27].