Optimal Design of Bubble Deck Concrete Slabs: Sensitivity Analysis and Numerical Homogenization

The use of layered or hollow floors in the construction of buildings obviously reduces the self-weight of the slab, and their design requires some expertise. In the present work, a sensitivity analysis and numerical homogenization were used to select the most important characteristics of bubble deck floors that have a direct or indirect impact on their load capacity. From the extensive case study, conclusions were drawn regarding the optimal selection of geometry, materials, and the arrangement and size of air voids in such a way as to ensure high stiffness of the cross-section and at the same time maximally reduce the self-weight of the slabs. The conducted analyses showed that the height of the slab and the geometry of the voids had the greatest impact on the load-bearing capacity. The concrete class and reinforcement used are of secondary importance in the context of changes in load-bearing capacity. Both the type of steel and the amount of reinforcement has a rather small or negligible influence on the bubble deck stab stiffness. Of course, the geometry of the voids and their arrangement and shape have the greatest influence on the drop in the self-weight of the floor slabs. Based on the presented results of the sensitivity analysis combined with numerical homogenization, a set of the most important design parameters was ordered and selected for use in the optimization procedure.


Introduction
The history of ceilings as a covering dates back to the Stone Age. It was during the Paleolithic era that the vaulting of rocks was used as a barrier against atmospheric agents such as precipitation or solar radiation. As early as 5450-4400 BC, Egyptian civilization used palm trunks stacked tightly next to each other as a supporting structure for the vault. Additionally, in ancient Greece, stone architraves were used as support structures for wooden beams and were covered with planks or terracotta slabs [1,2]. The further evolution of stone ceilings in Greece contributed to the development of coffered ceilings. During subsequent eras, other materials began to be used to make different coverings.
The beginning of documented knowledge of reinforced concrete is generally considered to have occurred in 1867, when the Frenchman Joseph Monier patented a concrete pot reinforced with wire mesh. A few years later, in 1892, Johann Franz Kleine, a master mason from Essen, patented a ceiling consisting of steel beams and an infill of ceramic slabs reinforced with cooper [3]. With the beginning of the 20th century, structural thought developed and ceilings with concrete slabs on steel beams and ceilings with monolithic reinforced concrete beams appeared. The constant search for new structural solutions to meet innovative architectural visions led to the development of monolithic prestressed floors [4][5][6]. Due to the greater strength of the floors made of these materials, the design Table 1. Differences between floors with HDPE inserts [25].

Bubble Deck Cobiax
One spatial structure, upper and lower reinforcement combined Single baskets with balls combined with only one (upper or lower) steel mesh Reinforcement to match the cartridges (balls) Bottom and top reinforcement independent of inserts (balls) Plate thicknesses 25-60 cm Thicknesses from 20-60 cm There is a danger of damaging the bullets No danger of damaging the bullets Such ceilings are usually designed in an analytical manner using the available standards and guidelines contained in manufacturers' catalogs. However, this approach is very time consuming and impractical because, in most cases, the engineering solution must be obtained quickly and must be repeatable for different load schemes and structures. Additionally, the use of full three-dimensional FE models of such floors in analyses involves a lot of work and is uneconomical. The solution to the above problem may be the use of simplified models obtained as a result of homogenization.
The homogenization technique allows scholars to replace the complex multilayer composite cross-section with a model with equivalent parameters. Therefore, the behavior of the replacement model is very similar to the initial (reference) model. Issues related to the homogenization of complex cross-sections of structures have been raised by many researchers for several decades. It is possible to distinguish, e.g., the method of periodic homogenization presented by Baunnic [26], the method based on an inverse analysis [27], or the homogenization method based on strain energy [28]. The differences between static and dynamic homogenization are presented in [29]. On the other hand, a list and comparison of different homogenization methods for sandwich and composite panels with a corrugated core can be found in [30]. Xin et al. proposed an approach based on the mechanics of the genome structure [31,32]. Biancolini proposed the homogenization method based on equivalent strain energy [33], in 2021, this method was improved by Garbowski and Gajewski [34]. Homogenization based on strain energy equivalence is a universal method that can be used not only for cardboard homogenization. It was used, e.g., for the homogenization of the cross-section of prefabricated Filigran slabs [35] and the analysis of thin-walled beams [36].
This paper presents the particular application of the numerical homogenization method based on the equivalence of strain energy between the full 3D RVE model and the simplified model of a single-layer coating. Initially, the method was used only for the analysis of shell structures [33,34]. This paper presents the application of the method to solid elements with internal and regular holes as well as truss elements. The application of the homogenization technique is shown on the example of a bubble deck ceiling with evenly spaced openings over the entire surface. The adopted homogenization method was used here together with a standard sensitivity analysis to automate the optimization of the numerical model of bubble deck reinforced concrete slabs. The presented methodology allows scholars to easily define both the objective function and fully automatically control all the design parameters of such ceilings. The novelty of the presented work is the algorithm, which was prepared to automatically optimize any reinforced concrete ceiling with regular openings, and the original homogenization method (adapted to the problem analyzed in this paper), which significantly improves the process of optimally designing ceilings constructed with bubble deck slabs.

Representative Volume Element
The numerical examples considered in this article were reinforced, concrete bubble deck slabs. The analyses were carried out for two selected versions of the above ceiling (BD340 and BD230) using the numerical homogenization method, the details of which are described in Section 2.2. The influence of particular floor slab parameters on the change in its stiffness and weight was examined. For this purpose, it was necessary to correctly define the stiffness matrix RVE (representative volume element) of the floor element. The RVE was created as a result of separating a repetitive fragment of the structure from the entire 3D plate model. In the first model, a typical BD340 bubble deck floor was analyzed. The basic RVE of the slab was 0.34 m high, 0.30 m wide, and 0.30 m long, see Figure 1a. The width and length were adopted in such a way as to maintain the minimum required distance between the bubbles for this type of ceiling, which was 0.03 m. The void was of a spherical shape with a diameter of 0.27 m. The slab was modeled as a solid 3D structure made of concrete. The bidirectional reinforced plate had a degree of bending reinforcement of 0.66% in both directions. This gave the reinforcement mesh the ability to move in two directions, up and down, and the mesh consisted of bars with a diameter of φ = 12 mm spaced every 0.10 m, see Figure 1. These elements were modeled as steel wire structures embedded in the concrete slab. In the model, the concrete cross-section was divided into C3D10 solid elements, i.e., 10-node quadric tetrahedron elements with three degrees of freedom with dimensions of 0.05 m in each RVE direction. The steel reinforcement was meshed with truss elements T3D2 (2-node truss elements with three degrees of freedom according to [37]). The discretization of the BD340 model gives a total number of 5289 elements and 9198 nodes, which translates to 27,594 degrees of freedom (since all the elements have three degrees of freedom in each node).
In the numerical model 2, a BD230 bubble deck slab was considered. An RVE of 0.25 × 0.25 × 0.23 m (width × length × height) was modeled, see Figure 2a. The plate had an elliptical bubble with dimensions of 0.23 × 0.18 m. The distance between the hole was equal to 0.02 m. As in model 1, the reinforcement was modeled as wire structures and was made of steel, while the plate was a concrete 3D solid structure. The degree of reinforcement of the slab equal to 0.79% in two directions was assumed. In the upper and lower reinforcement mesh, bars with a diameter = 12 mm were spaced every 0.125 cm in two directions (Figure 2b). The concrete slab was divided into a numerical model by using C3D10 elements with dimensions of 0.05 m. The steel bar was meshed with T3D2-type elements. The total number of elements in this model equaled 1866, which gave 3486 nodes and 10,458 degrees of freedom in total.
In the analyzed examples, two materials were used to describe the mechanical behavior of the individual elements of the models: steel and concrete class C30/37. The steel was described by the isotropic linear elasticity relations between stress and strain. For concrete, the linear stress-strain relations were used. In Table 2, the engineering parameters of the materials used in the models are shown, where is the Young's modulus, is the Poisson′s ratio, and is the density.  In the model, the concrete cross-section was divided into C3D10 solid elements, i.e., 10-node quadric tetrahedron elements with three degrees of freedom with dimensions of 0.05 m in each RVE direction. The steel reinforcement was meshed with truss elements T3D2 (2-node truss elements with three degrees of freedom according to [37]). The discretization of the BD340 model gives a total number of 5289 elements and 9198 nodes, which translates to 27,594 degrees of freedom (since all the elements have three degrees of freedom in each node).
In the numerical model 2, a BD230 bubble deck slab was considered. An RVE of 0.25 × 0.25 × 0.23 m (width × length × height) was modeled, see Figure 2a. The plate had an elliptical bubble with dimensions of 0.23 × 0.18 m. The distance between the hole was equal to 0.02 m. As in model 1, the reinforcement was modeled as wire structures and was made of steel, while the plate was a concrete 3D solid structure. The degree of reinforcement of the slab equal to 0.79% in two directions was assumed. In the upper and lower reinforcement mesh, bars with a diameter φ = 12 mm were spaced every 0.125 cm in two directions ( Figure 2b). The concrete slab was divided into a numerical model by using C3D10 elements with dimensions of 0.05 m. The steel bar was meshed with T3D2-type elements. The total number of elements in this model equaled 1866, which gave 3486 nodes and 10,458 degrees of freedom in total.

Numerical Homogenization Algorithm
All calculations in this article were made based on the principle of numerical homogenization, which uses the energy equivalence between the full three-dimensional model and the simplified shell model. In the full model, the complete RVE geometry (with all details) is defined, which is then discretized by applying a finite element mesh, and the constitutive relationships for each material are also defined. Since the equivalent stiffness In the analyzed examples, two materials were used to describe the mechanical behavior of the individual elements of the models: steel and concrete class C30/37. The steel was described by the isotropic linear elasticity relations between stress and strain. For concrete, the linear stress-strain relations were used. In Table 2, the engineering parameters of the materials used in the models are shown, where E is the Young's modulus, ν is the Poisson's ratio, and ρ is the density.

Numerical Homogenization Algorithm
All calculations in this article were made based on the principle of numerical homogenization, which uses the energy equivalence between the full three-dimensional model and the simplified shell model. In the full model, the complete RVE geometry (with all details) is defined, which is then discretized by applying a finite element mesh, and the constitutive relationships for each material are also defined. Since the equivalent stiffness is determined for a certain element whose surface and dimensions are defined by the projection of the RVE onto the x-y surface, in the homogenization process, only the nodes on the vertical RVE surfaces are used (see Figure 3). Therefore, the overall stiffness of the entire model is condensed to those external nodes by using the static condensation method.

Numerical Homogenization Algorithm
All calculations in this article were made based on the principle of numerical homogenization, which uses the energy equivalence between the full three-dimensional model and the simplified shell model. In the full model, the complete RVE geometry (with all details) is defined, which is then discretized by applying a finite element mesh, and the constitutive relationships for each material are also defined. Since the equivalent stiffness is determined for a certain element whose surface and dimensions are defined by the projection of the RVE onto the x-y surface, in the homogenization process, only the nodes on the vertical RVE surfaces are used (see Figure 3). Therefore, the overall stiffness of the entire model is condensed to those external nodes by using the static condensation method.
The resulting (equivalent) model is reduced to a shell that is described by the tensile, bending, and transversal stiffnesses in both the x and y directions. In this process, both the geometry and the materials (which here are the concrete plate, steel bars, and voids) are homogenized into equivalent structural stiffnesses, on the basis of which, of course, the equivalent material parameters and the effective thickness of such a shell can be calculated. However, in the example analyzed in this paper, this is not necessary; therefore, only the basic steps leading to obtaining equivalent structural stiffnesses will be briefly discussed here. The resulting (equivalent) model is reduced to a shell that is described by the tensile, bending, and transversal stiffnesses in both the x and y directions. In this process, both the geometry and the materials (which here are the concrete plate, steel bars, and voids) are homogenized into equivalent structural stiffnesses, on the basis of which, of course, the equivalent material parameters and the effective thickness of such a shell can be calculated. However, in the example analyzed in this paper, this is not necessary; therefore, only the basic steps leading to obtaining equivalent structural stiffnesses will be briefly discussed here.
The generalized displacements at each node on the RVE surface are related to the generalized strains. Therefore, the relationship between the generalized constant strains and the position of the external nodes on the RVE boundary is expressed by the following transformation: where u is the node displacement vector and is the strain vector. Here, for a single node (x i = x, y i = y, z i = z), the W i matrix adopted for the RVE model can be easily derived: Matrix W i determines the relationship between displacements and effective deformations that are applied to nodes in boundary conditions to which the stiffness of the entire model is condensed. The total energy of elastic deformation is: where K is the global stiffness matrix. Taking into account the fact that the finite element model is subjected to bending, tension, and transverse shear, for the shell (or plate), the internal energy is: The stiffness matrix for a homogenized composite is easy to extract from a discrete matrix because: The matrix W k is the ABDR matrix, which can be saved as: where A contains the tensile and shear stiffnesses, B contains the combination of tension and bending stiffness, D contains the bending and torsional stiffness, and R contains the transverse shear stiffness.

Sensitivity Analysis
The sensitivity analysis of the bubble deck slabs was examined in terms of changes in the weight of the RVE and its individual stiffnesses, i.e., bending, shearing, and tension/compression stiffness. The parameters of the model were (a) the class of concrete; (b) the height of the slab; (c,d) the dimensions of the bubble or the distance between the holes; and (e,f) the diameter and spacing of the reinforcement. In each subsequent model, one parameter was changed, and the remaining parameters were identical to the reference model (BD340 and BD230). Then, parametric models were used to build the representative volume element that is necessary for numerical homogenization analyses. An individual RVE was built in the same way as the base models described in Section 2.1.
The diagram of the conducted analyses is shown in Figure 4. Although the presented algorithm does not look very complicated, it requires the use of many advanced techniques. After defining the set of input parameters, a python script was launched, and its task was to automatically create any RVE geometry based on the defined parameters. As a reminder, these parameters were concrete class, diameter and spacing of bars, slab height, and diameter and spacing of holes. Then, the defined model was automatically opened in the ABAQUS program, where only the finite element mesh was generated and the global stiffness matrix was assembled. On the basis of the saved data, i.e., node coordinates and the global stiffness matrix, in the next step, static condensation and the W i matrix were computed. These calculations were performed in the MATLAB program, which controlled the whole procedure. On the basis of the calculated W i matrix and the condensed K matrix, the W k matrix was determined (see Equation (6)) for each analyzed geometry. All the scripts written in different environments created a coherent and fully automatic procedure that, based on the input parameters, generated the ABDR stiffnesses of any RVE model, which allowed for the efficient performance of the sensitivity analysis.
The diagram of the conducted analyses is shown in Figure 4. Although the presented algorithm does not look very complicated, it requires the use of many advanced techniques. After defining the set of input parameters, a python script was launched, and its task was to automatically create any RVE geometry based on the defined parameters. As a reminder, these parameters were concrete class, diameter and spacing of bars, slab height, and diameter and spacing of holes. Then, the defined model was automatically opened in the ABAQUS program, where only the finite element mesh was generated and the global stiffness matrix was assembled. On the basis of the saved data, i.e., node coordinates and the global stiffness matrix, in the next step, static condensation and the matrix were computed. These calculations were performed in the MATLAB program, which controlled the whole procedure. On the basis of the calculated matrix and the condensed matrix, the matrix was determined (see Equation (6)) for each analyzed geometry. All the scripts written in different environments created a coherent and fully automatic procedure that, based on the input parameters, generated the ABDR stiffnesses of any RVE model, which allowed for the efficient performance of the sensitivity analysis.

Results
This section presents the results obtained from the analysis of the bubble deck slab using both the numerical homogenization and the sensitivity analysis. The influence of the selected cross-section parameters on the change in the stiffness and weight of the floor was analyzed. The analyses were carried out for two variants of the BD 340 and BD 230 slab. First, the stiffnesses for the reference models described in Section 2.1 were calculated. The results obtained using the homogenization method (Section 2.2) are presented in Table 3. Then, cases with a change in subsequent cross-section parameters were analyzed. The calculated changes in stiffness are presented in the tables with comparison to the reference values. Table 3 compiles the stiffness submatrices, namely A, D, and the transverse stiffness matrix R (represented in Equation (6)) for both of the RVE reference models.

Results
This section presents the results obtained from the analysis of the bubble deck slab using both the numerical homogenization and the sensitivity analysis. The influence of the selected cross-section parameters on the change in the stiffness and weight of the floor was analyzed. The analyses were carried out for two variants of the BD 340 and BD 230 slab. First, the stiffnesses for the reference models described in Section 2.1 were calculated. The results obtained using the homogenization method (Section 2.2) are presented in Table 3. Then, cases with a change in subsequent cross-section parameters were analyzed. The calculated changes in stiffness are presented in the tables with comparison to the reference values.  Table 3 compiles the stiffness submatrices, namely A, D, and the transverse stiffness matrix R (represented in Equation (6)) for both of the RVE reference models.
The slab weight is 5.66 kN/m 2 for the BD 340 variant and 3.90 kN/m 2 for the BD 320 variant.

BD 340
This subsection presents the results obtained from the sensitivity analysis based on the numerical homogenization of the bubble ceiling, variant BD340, due to the change in the individual geometrical and material parameters of the cross-section. In each case, only one characteristic of the RVE was changed. Table 4 shows the stiffness results of the BD340 slab obtained as a result of analyses after increasing the selected cross-section or material parameters. The class of concrete was changed to C35/45 (E = 34 GPa), i.e., the Young's modulus of the material increased by 6.25% compared to the reference value. The perturbation of the other parameters was not the same; for example, the section height (H), hole diameter (D), and rebar spacing were increased by 5.0%, obtaining, respectively, a height equal to H = 0.357 m, a bubble diameter D = 0.284 m, and 0.105 m for the spacing of the rebar. In order to keep realistic values, the diameter of the bars was assumed as φ = 0.014 m (therefore, the dimension increased by 16.67% compared to the initial RVE model). When increasing the class of the concrete, the weight of the slab did not change and was equal to 5.66 kN/m 2 . For the other variants of the bubble deck BD340, the weight of the plate changed and was 6.07 kN/m 2 (for H = 0.357 m), 5.22 kN/m 2 (D = 0.284 m), 5.75 kN/m 2 (for φ = 0.014 m), and 5.65 kN/m 2 for bar spacings every 0.105 m. Table 5 presents the stiffness resulting from the reduction in the individual crosssection and material parameters for the BD340 bubble slab. The Young's modulus of concrete decreased by 6.25% compared to the reference model (concrete class C25/30-E = 30 GPa was assumed). The height of the cross-section, the spacing of the reinforcement bars, and the diameter of the bubble were, respectively, H = 0.323 m, 0.095 m, and D = 0.256 m. In each of these cases, the given parameter decreased again by 5.0%. In the case of the reinforcement diameter, the reduction was 16.67%, which resulted in using the bars φ = 0.010 m in the analyses.
In order to better illustrate the influence of a given parameter, the weight of the slab was also calculated. As before, changing the concrete class to C25/30 had no effect on the weight of the slab, and it was 5.66 kN/m 2 . Reducing the section height by 5% reduced the slab weight to 5.25 kN/m 2 . Assuming a smaller diameter of the bars (φ = 0.010 m), the mass of the bubble plate was equal to 5.59 kN/m 2 , while the reduction in the bubble diameter by 5% compared to the reference model increased the weight of the plate and amounted to 6.05 kN/m 2 . Changing the spacing of the bars (0.095 m) did not significantly affect the weight of the bubble plate, which was 5.67 kN/m 2 .

BD 230
In the next step, the sensitivity analysis based on numerical homogenization was used to compute the BD230 bubble slab type. As for the BD340 variant, only one parameter was increased/decreased. The height was increased by 5% to H = 0.242 m. The dimensions of the elliptical bubble changed to D 1 = 0.189 m and D 2 = 0.242 m (increase in dimensions by 5.0% compared to the reference model). The concrete class was changed to C35/45 (increase in Young's modulus by 6.25%). However, the diameter and spacing of the reinforcing bars increased by 5.0% and 16.67%, respectively, to 0.014 m and 0.131 m. The comparison of the stiffnesses obtained from the plates with a changed selected cross-sectional parameter is presented in Table 6. The weight of the slab did not change in the case of increasing the concrete class and was equal to 3.90 kN/m 2 . Increasing the dimensions of the holes caused a reduction in the mass of the bubble plate to 3.81 kN/m 2 (for D 1 = 0.189 m) and 3.71 kN/m 2 (for D2 = 0.242 m). For height H = 0.242 m, the weight of the plate was equal to 4.18 kN/m 2 ; for reinforcement with diameter φ = 0.014 m, the weight was equal to 4.01 kN/m 2 ; and the weight was equal to 3.88 kN/m 2 when the bars were spaced 0.131 m.
The change in stiffness due to a decrease in the value of the selected parameters of cross-sections of the BD230 bubble plate type is presented in Table 7. As previously mentioned, a change in the six parameters was applied, i.e., a change in the concrete class to C25/30 (a reduction in Young's modulus by 6.25%) and a decrease in the slab height, bubble dimensions, and bar spacing by 5.0%, which resulted in H = 0.219 m, D 1 = 0.171 m, D 2 = 0.219 m, respectively, and the distance between the reinforcements being equal to 0.119 m. The diameter of the rods was 0.010 m, which was a decrease in the value of the parameter by 16.67% compared to the initial model. As before, the change in the concrete class did not affect the value of the slab weight and was exactly the same as in the reference model

Discussion
The presented results are limited only to changes in the bending, tension/compression, and transversal shear stiffness of plates with holes due to the perturbations of the selected parameters. The results were obtained for the two types of bubble plates using a sensitivity analysis and numerical homogenization method based on strain energy equivalence. The analyzed slabs had reinforcement and spherical or elliptical bubbles evenly distributed over their entire surface. The conducted analyses allowed us to obtain information on the influence of individual parameters of cross-sections and materials on the plate stiffness for the considered plate. In addition, the influence of the selected parameters on the change in weight was examined. This was a basic step before the optimization of such panels, as it showed which parameters were worth changing and how their change was reflected in the calculated stiffness and weight of the panel. Table 3 shows the stiffnesses of the bubble plates for the reference models obtained from the sensitivity analysis and numerical homogenization. In turn, Tables 4-7 present the stiffness of the slabs depending on the change in the parameter, i.e., concrete class, cross-section height, bubble size, diameter of reinforcing bars, and their spacing. The stiffnesses caused by the increase in the values of the individual parameters are shown in Tables 4 and 6, while those determined after reducing the given feature are shown in Tables 5 and 7. For the considered cases, it can be seen that the change in the spacing of the reinforcing bars by 5.0% did not significantly affect the difference in slab stiffness, and we even obtained the same value as in the reference model. Additionally, the diameter of the reinforcement did not affect the shear stiffness, and its effect on the values of compression/tension and bending stiffness was in the range of 2.0-4.0% for the analyzed examples (assuming that the diameter of the reinforcement varied +/−16.67%). The height of the cross-section had the greatest impact on the differences in the obtained results between the reference model and the modified model. When increasing the height of the plate by 5.0%, the bending stiffness increased in the range of about 17.5-19.3% compared to the original model for BD340 and 18.0-20.8% for BD230. For a height reduction of 5.0%, the bending stiffness was reduced by 15.7-17.0% for the first case and 16.0-18.0% for the second case. In addition, it can be seen that the greatest change in tensile/compressive stiffness occurred for different section heights.
The use of a higher class of concrete (C35/45) caused an increase in all stiffnesses by about 5.8-6.6%, while the use of concrete of a lower class (C25/30) reduced the stiffness by about 6.0%. The obtained results showed that the shape of the bubble used in the slab affected the stiffness values in a variable way. Based on the data in Tables 3-7, it can be seen that the spherical bore had the largest change in shear stiffness and the smallest change in bending stiffness. In the case of an elliptical hole, changing the vertical dimension of the bubble affected the stiffness to a lesser extent than its horizontal dimension. The diameter of the reinforcing bars did not change the shear stiffness, and the maximum increase/decrease in bending stiffness was 5% (assuming that the bar diameter changed +/−16.67% compared to the reference model).
In addition, on the basis of the analysis, it can be noted that the change in the class of concrete did not affect the change in the weight of the slab. The greatest impact on the change in the weight of the bubble plate was the increase/decrease in the height of the slab for the analyzed cases and the spherical opening in the case of BD340. A smaller change in the weight of the plate occurred when elliptical bubbles were used than when spherical bubbles were used. The spacing of the bars caused a slight change in the weight of the slab compared to the initial model. In order to make a fair comparison of all the observations presented here, the following normalization procedure was performed. First, the central difference of each change in the measured stiffness value due to the perturbation (decrease and increase) of each analyzed parameter was calculated, and then the results were normalized using the reference value of each stiffness and the actual perturbation step. The normalized sensitivities were computed using the following formula: where R + i is the i-th stiffness computed for the i-th positive perturbation h i , R − i is the i-th stiffness computed for the i-th negative perturbation (−h i ), and R 0 i is i-th reference stiffness value. Tables 8 and 9 present the normalized stiffnesses ADR computed with Equation (7) for the BD340 and BD230 concrete bubble deck slab.
Many researchers recently analyzed the structural, thermal, and acoustic properties of reinforced concrete slabs with balls [38]. In this paper, a procedure is presented by means of which only the structural-strength aspects of hollow plates can be analyzed. Of course, the procedure proposed here can also be extended to an acoustic and thermal analysis of concrete slabs by adding appropriate modules (in Figure 4) in which the necessary analyses would be carried out. General guidelines for the optimal design of reinforced concrete bubble deck slabs can be presented in three main steps: (1) determining the necessary stiffness value for the analyzed floor; (2) specifying additional constraints in the form of, for example, the thermal and acoustic conductivity or the dead weight of the ceiling; and (3) running any optimization algorithm based on the numerical homogenization procedure presented in this paper.  Although the procedure presented in this work has not been experimentally verified, the homogenization method itself, used in the presented procedure, has already been verified on various engineering examples in our previous works [35,36,[39][40][41]. The proposed algorithm can of course be implemented in a completely different computing environment; however, the procedure should be consistent with the diagram presented in Figure 4 and the general design guidelines presented in the previous paragraph.

Conclusions
This sensitivity analysis combined with numerical homogenization can be used in the computer-aided design of optimal bubble deck slabs in construction engineering, etc. The use of the homogenization method adopted for concrete reinforce boards with periodic holes can reduce the calculation time, which is always important in many problems, e.g., strength analysis and slab optimization. Based on the analyses carried out in this work, very important answers were obtained regarding the sensitivity of the basic stiffnesses of the bubble deck slab model and the basic design parameters. The observed relationships clearly indicated that the most important design parameters that play the most important role in the optimization process are the ceiling height and the arrangement and shape of the voids. The class of concrete is twice as less important as the ceiling height, and the class of steel and the diameter of the bars have little effect on the stiffness of the slab cross-section.
The presented algorithm is only the initial step in the optimization procedure, where one can simply change any input parameter and automatically calculate the required stiffness, whether bending or tensile/compressive, along with the change in the basic weight. The algorithmic procedure presented here allows for a quick and comprehensive optimization, which we are going to present in our next article. The sensitivities calculated here will constitute the basic information to appropriately select the optimization method and algorithms, as well as to eliminate the parameters whose impact on individual stiffnesses is negligible.

Conflicts of Interest:
The authors declare no conflict of interest.