Surface roughness prediction of FFF-fabricated workpieces by artificial neural network and Box–Behnken method

Fused Filament Fabrication (FFF) or Fused Deposition Modelling (FDM) or three-dimension (3D) printing are rapid prototyping processes for workpieces. There are many factors which have a significant effect on surface quality, including bed temperature, printing speed, and layer thickness. This empirical study was conducted to determine the relationship between the above-mentioned factors and average surface roughness (Ra). Workpieces of cylindrical shape were fabricated by an FFF system with a Polylactic acid (PLA) filament. The surface roughness wasmeasured at five different positions on the bottom and top surface. A response surface (Box-Behnken) method was utilised to design the experiment and statistically predict the response. The total number of treatments was sixteen, while five measurements (Ra1, Ra2, Ra3, Ra4 and Ra5) were carried out for each treatment. The settings of each factor were as follows: bed temperature (80, 85, and 90 °C), printing speed (40, 80 and 120mm/s), and layer thickness (0.10, 0.25 and 0.40mm). The prediction equation of surface roughness was then derived from the analysis. The same set of data was also used as the inputs for a machine learning method, an artificial neural network (ANN), to construct the prediction equation of surface roughness. Rectified linear unit (ReLU) was utilised as the activation function of ANN. Two training algorithms (resilient backpropagation with weight backtracking and globally convergent resilient backpropagation) were applied to train multi-layer perceptrons. Moreover, the different number of neurons in each hidden layer was also studied and compared. Another interesting aspect of this study is that the ANN was based on a limited number of training samples. Finally, the prediction errors of each method were compared, to benchmark the prediction performance of the two methods: Box-Behnken and ANN.


Introduction
Surface quality is an important characteristic of workpieces fabricated by a Fused Filament Fabrication (FFF) system; the lowest surface roughness is desirable. There are many factors affecting the quality of the surface, e.g., bed temperature, layer thickness, print speed, and raster width. As a result, it is critical to numerically quantify the effect of these factors on the surface roughness. According to the literature, different methods have been used to study this effect and experimental design was one of them. Another method was the artificial neural network (ANN), which was based on a machine learning algorithm. However, the challenging issue of the ANN application is the limited number of training and testing data, since the data is obtained from experiments which consume resources. Another dilemma is the question regarding the number of neurons in each hidden layer. Moreover, the right selection of training algorithms also plays an important role in prediction capability.
Therefore, the main aim of this study is the assessment of surface roughness prediction by a response surface method, Box-Behnken, and ANN method. Another objective is the performance comparison of two ANN training algorithms: resilient backpropagation with weight backtracking (RPROP+) and globally convergent resilient backpropagation (GRPROP). Finally, the utilisation of different number of neurons in each hidden layer is also studied and compared.

Literature review
Several studies have focused on the application of different training algorithms of the ANN method to model the mechanical properties of workpieces fabricated by FFF systems.
Mahapatra and Sood [1] applied ANN to determine the relationship between five input factors (raster angle, raster width, layer thickness, orientation, and air gap) while the surface roughness at the top and bottom position were measured. The algorithm used to train the neural network was the Levenberg-Marquardt algorithm.
Sood, Ohdar, and Mahapatra [2] applied the ANN with resilient backpropagation, to model the compressive strength of specimens built by the FDM process. The data set used to train the network was obtained from the central composite design of the experiment.
Nagarajan et al. [3] had improved the product quality of the FDM manufactured parts by introducing the knowledge-based ANN, which was a combination of the classical ANN and dimensional analysis conceptual modelling (DACM). The classical ANN was based on the Levenberg-Marquardt training algorithm. Jiang et al. [4] utilised the Bayesian Regularization as the training algorithm, to develop the neural network and predict the printable bridge length (PBL) of specimens fabricated by FDM.
A neural network with multi-layer perceptrons was used by Boschetto, Giordani, and Veniali [5] to predict the surface roughness of test specimens made from Acrylonitrile Butadiene Styrene (ABS). The experiment was conducted to study the effect of different activation functions and numbers of neurons. The accuracy of the models was based on these indexes: root mean square error (RMSE), mean absolute percentage error (MAPE), and coefficient of determination (R 2 ). Noriega et al. [6] combined the ANN and optimisation method to determine the optimal dimensional values of FDM-fabricated workpieces. The results point to the fact that the hybrid approach leads to significant reductions in the external and internal dimensions. Vijayaraghavan et al. [7] applied an improved approach of multi-gene genetic programming (Im-MGGP) to formulate the functional relationship between input factors, layer thickness, orientation, raster angle, raster width, air gap, and an output, wear strength.
Panda, Bahubalendruni and Biswal [8] compared two soft computing methods, multi-gene genetic programming (MGGP) and a general regression neural network (GRNN), to find the relationship between three input factors (layer thickness, orientation, and raster angle) and an output (the compressive strength of an FDM fabricated prototype). The study led to the conclusion that a model generated by GRNN had a better fit than the MGGP model. Bayraktar et al. [9] applied the ANN method to optimise the tensile strength of test specimens fabricated by the FDM method; a PLA filament was used as the material. The training algorithm used was the backpropagation (BP) with Levenberg-Marquardt algorithm (LMA). The experiments were conducted using three melt temperatures, four-layer thicknesses and three raster pattern orientations. Vahabli and Rahmati [10] applied the radial basis function neural network (RBFNN) to estimate the surface roughness of specimens. The imperialist competitive algorithm (ICA) was used to optimise the effective parameters of the effective variables. The results signify that the prediction error (MAPE) of the combined ICA-RBFNN was lower than that of RBFNN alone.
Other machine learning algorithms (e.g., random forests (RFs), support vector regression (SVR), ridge regression (RR), and least absolute shrinkage and selection operator (LASSO)), were utilised by Wu, Wei and Terpenny [11] to accurately predict the surface roughness of manufacturing parts. Li et al. [12] introduced the ensemble learning algorithm to determine the predictive model of surface roughness. There were two stages of calculation: offline and online predicting. For the offline, different learning algorithms, such as random forests (RF), AdaBoost, classification and regression trees (CART), SVR, RR, and random vector functional link (RVFL) networks were utilised for training the model. The model testing was conducted in this online stage. The input data was temperature and vibration data collected from the build-plate temperature and vibrations, the extruder temperature and vibration, and the temperature of the deposited material.
Jiang et al. [13] investigated the response of the FDM process and printable bridge length (PBL), in test specimens fabricated by the FDM method. The data set to train the ANN network was obtained from the experiment which followed the orthogonal analysis L 32 . BP with Bayesian Regularization was used as the training algorithm. According to the results, the application of BPNN leads to the accurate prediction of a test specimen's PBL. Barrios and Romero [14] deployed the decision-tree method, J48, random forest, and random tree, to predict the surface roughness of FDM parts. The dataset was categorised into two groups: training and testing. The design matrix followed the Taguchi design, L 27 , and there were five responses: layer height, print temperature, print speed, print acceleration, and flow rate. Yadav et al. [15] studied the relationship between three inputs (material density, infill density, and extrusion temperature) and an output (tensile strength). The material used was a polymer blend between ABS and Polyethylene Terephthalate Glycol (PETG). The ANN and genetic algorithm (GA) were utilised to predict the tensile strength.
On the other hand, the response surface and design of the experimental method were applied in order to screen and optimise the surface roughness. Plaza et al. [16] studied the effect of build orientation, feed rate, and layer thickness on the surface roughness of workpieces fabricated from PLA filament. The design of experiment method used was the factorial design (3 levels for build orientation, 4 levels of layer thickness, and 3 levels of feed rate). Three replica workpieces were fabricated for each treatment. The results showed that the layer thickness had a significant effect on the surface roughness while the feed rate had no effect at all. Moreover, the ANN, with feedforward BP (training algorithm was LM), was applied to predict the surface roughness of workpieces. The central composite design (CCD) method was utilised by Yang et al. [17], to optimise the tensile strength, surface roughness, and build time of fabricated workpieces. The results indicated that nozzle diameter, filling velocity, and layer thickness had a significant impact on the surface roughness. Kandananond [18] applied the Box-Behnken design to optimise the surface roughness of workpieces fabricated from ABS filament. The optimal settings of the input factors, nozzle temperature, bed temperature, and printing speed were recommended, to minimise the surface roughness.
The response surface method (RSM) and machine learning were also utilised to study the relationship between inputs and other outputs of the FDM system. The comparison study of RSM and machine learning method, ANN and fuzzy inference system (FIS), as well as GA, was conducted by Peng, Xiao and Yue [19] to benchmark the prediction accuracy of the dimensional error, warp deformation, and built time when the inputs were: line width compensation, extrusion velocity, filling velocity, and layer thickness. Mohamed, Masood, and Bhowmik [20] compared the performance of two methods (the definitive screening design and ANN), to predict two parameters, creep compliance and recoverable compliance, of the FDM fabricated test specimens. The model of ANN used was multi-layer perceptron. The performance of ANN method was measured by the prediction error.

FFF system and workpiece model
The FDM or FFF unit used is a consumer-grade system equipped with a single nozzle, which has a maximum extrusion temperature of 220°C and maximum bed temperature of 90°C. The highest printing speed is 120 mm/s. The shape of the workpiece in this study was cylindrical and the dimensions were as follows: diameter = 40 mm and height = 10 mm (adapted from the work by Pérez et al. [21]). The infill density was ten percent and the structure was honeycomb. The material used was PLA filament, which is naturally degradable. Regarding the workpiece, the raster angle was set at 45°and the pattern was rectilinear. The solid model and its dimensions are shown in Figure 1 and the sample workpiece is depicted in Figure 2.  The surface roughness was measured at the top and bottom surface. At each surface, the measurement of surface roughness was conducted at five different locations (as shown in Fig. 3), four positions with a 90°separation and the centre of the workpiece. The average surface roughness at these locations was denoted Ra 1 , Ra 2 , Ra 3 , Ra 4 , and Ra 5 . The surface roughness tester was used to measure the surface roughness.

Experimental design
In this study, the Box-Behnken design was applied to create the design matrix. As a result, the levels of treatments were divided into three levels. There are three input factors, bed temperature, print speed, and layer thickness, while the output is surface roughness. The bed temperature is the heating temperature of the print bed, so the filament material from the extruder is slowly cooled. The lowest bed temperature for this experiment was 80°C and it increased by five degrees, so the second and third levels were 85 and 90°C. The printing speed indicates how fast the workpiece is fabricated, the slowest speed was 40 mm/s. The second level of printing speed was 80 mm/s, which is a moderate speed, and the fastest speed was 120 mm/s. The layer thicknesses were the layer heights of each successive addition of material, which were stacked. The thinnest layer was 0.1 mm, while the thickest was 0.4 mm. For the Box-Behnken design, the designated levels are categorised as coded levels (À1 = low, 0 = middle and 1 = high). The experimental factors and their levels, in the form of actual and coded values, are shown in Table 1.
Because the measurements were conducted at five positions, the response at each position can be considered as a block. As a result, the blocking effect is also quantified and analysed in this study.

Results
Since the Box-Behnken design is based on the three levels of each factor, this design requires sixteen design points. At each design point, there are five measurements at different positions, expressed as Ra 1 -Ra 5 (mm). The results are divided into two cases: top and bottom surface, as shown in Tables 2 and 3.

Top surface
The top surface roughness of the workpieces was measured and is presented in Table 2.

Bottom surface
The bottom surface roughness of workpieces and their design matrix are shown in Table 3.

Analysis
Two prediction methods, Box-Behnken RSM and ANN, were applied to quantify the effect of each input factor and analyse the relationship between input factors and surface roughness. Afterwards, the prediction equation was derived to determine the surface roughness.

Box-Behnken RSM
Box-Behnken is an RSM focusing on the optimisation of the response. This method was introduced by Box and Behnken [22]. To achieve the desired target, the experiment was conducted to collect the data at different treatments including low, middle, and high levels. Moreover, the response data at the centre points was also collected, to test the lack-of-fit. Some design points of this response surface design are graphically represented in Figure 4.

Artificial neural network (ANN)
Normally, after the experiment was carried out, the obtained data was used to train and test the network. The structure of the network consists of the input and output layers, while the hidden layer connects both layers together. In this study, the training inputs are bed temperature, printing speed, and layer thickness while the training output is surface roughness. The structure of neural network is illustrated in Figure 5. The algorithms, RPROP+ and GRPROP, were utilised to achieve the optimal weight so that the neural network was trained to appropriately map inputs to an output. According to the ANN structure, the inputs were fed into each hidden neuron. At each hidden neuron, the weighted sum of the inputs, including bias, were added together. The summation of weighted inputs and bias to the hidden neuron h 11 is shown in equation (1).
where Input h11 is the summation of the product of the weight and inputs to the hidden neuron h 11 , including bias; w 11 is the weight of the first input to the hidden neuron h 11 ;  w 21 is the weight of the second input to the hidden neuron h 11 ; w 31 is the weight of the third input to the hidden neuron h 11 ; b 11 is the bias to the hidden neuron h 11 . Next, the weighted sum was then passed through the activation function to limit the amplitude of the output. In this study, the activation function was set at ReLU as elaborated in equation (2). The activation function is a function that is applied to smooth the product of the weights and inputs including bias. Afterwards, the available data is divided into two groups for training and testing purposes. The training dataset is a set of examples used to fit the parameters of the model (e.g., the weights of the connections between neurons in ANNs).
The graphical representation of the activation function (ReLU) is shown in Figure 6.
For the network training process, the weights are adjusted by using different learning algorithms and the objective is the minimisation of prediction errors. Two algorithms, RPROP+ and GRPROP, were utilised in this study.
According to [23] and [24], RPROP+ is the adaptation of the resilient backpropagation (RPROP) algorithm. Basically, the update weight of RPROP is dependent on the sign of the partial derivative of error function as shown in equation (3). If the sign of derivative is positive, the updated weight will decrease. However, the update weight will increase when the sign of derivative is negative. The process of weight update stops when ∂E ∂w jk ¼ 0. The calculation of update-value D jk follows the adaptation rule in equation (4) indicating how the update-value changes according to the sign of product of two partial derivatives. If the product is positive (signs of two consecutive partial derivatives do not change), the update value will increase with the factor h + . On the other hand, if the signs of two partial derivatives change (negative product implying that the weight skips the optimal value), the update-value will decrease with the factor hÀ.  where w jk is the weight connecting neuron j to k; D jk is the update-value, h + and h À are the increase and decrease factors and 0 < h À < 1 < h + , ∂E ðtÀ1Þ ∂w ðtÀ1Þ jk is the partial derivative of the error function with respect to w jk at the (t À 1) th iteration. However, for RPROP+, when the product of two partial derivatives is negative (indicating that the optimal weight is skipped), the step of weight backtracking is added to the RPROP algorithm to revert the weight-update. GRPROP is another modified version of RPROP and its objective is to improve the convergence of the optimal weight [25]. At each step, the descent direction is identified by utilizing RPROP algorithm to calculate the individual direction component of each weight dimension [26]

Analysis results
The surface analysis is divided into two cases: bottom and top surface roughness, For the Box-Behnken analysis, a DOE software (Design Expert version 12) was used to analyse the results, while the programming in R was applied to determine the ANN structures for both responses.

Box-Behnken analysis of bottom surface
The Box-Behnken analysis of bottom surface roughness is illustrated in the analysis of variance (ANOVA), as shown in Table 4. The ANOVA indicates that the main sources of variation combine a quadratic effect from printing speed (B 2 ) and the interaction between bed temperature and printing speed (AB). The main effect from printing speed (B) and layer thickness (C) is also included in the model. It is interesting to note that the bed temperature (A) is not excluded from the model because of the hierarchical effect, although the main effect of A is not significant (p-value = 0.3447). The measuring position or block does not seem to have a significant effect on the surface roughness, since the p-value is high. The Box-Behnken design also separates the residual into the lack-of-fit and pure error terms. The lack-of-fit is the variability that the prediction model cannot explain. On the other hand, the pure error represents the variability that is random and it is obtained from the replicated runs. Since the result derived from the ANOVA indicates that the lack-of-fit is not significant (p-value = 0.0844), the prediction model is capable of explaining the bottom surface roughness. The prediction model for bottom surface roughness is expressed in equation (6).
The normal probability plot (Fig. 7) illustrates that the points follow the pattern of a straight-line, so the residual is normally distributed. Therefore, the surface roughness can be explained by the proposed model.
The main effect of layer thickness (see Fig. 8) indicates that the surface roughness at a layer thickness of 0.1 mm is significantly lower than surface roughness at a layer thickness of 0.4 mm.
The interaction plot shown in Figure 9 illustrates that the surface roughness could be minimised in two scenarios: low bed temperature and high printing speed or high bed temperature and low printing speed. After the analysis was carried out, a workpiece was fabricated for validation purposes as indicated in Table 5. The bottom surface roughness was measured (Ra = 4.13) and compared to the Ra predicted by the proposed equation (3.9775). The results prove that measured Ra and predicted Ra are not significantly different.

ANN analysis of bottom surface
At the beginning, the whole data set was normalised and the min-max method was chosen to scale the data in the interval between 0 and 1. After this, the data was randomly split into two groups: training and testing. The next step was to fit the neural network to the training data. The number of hidden layers and neurons as well as the activation function (ReLU) were identified at this stage. In this study, the number of hidden layers was set at two while the number of neurons in the first layer was between 3 and 5. The number of neurons in the second layer was between 2 and 4. The algorithm of neural network training was coded and run in R, version 3.5.1. The system used has the following specifications: the processor speed of 1.66 GHz and 4 GB RAM. The training algorithms used were RPROP+ and GRPROP. The stopping criteria (based on the threshold of the partial derivatives of the error function) was set at 0.01. The number of neurons in each hidden layer was determined by varying the number of neurons for different ANN structures. After the appropriate neural network was fit to the training data, the prediction error, in the form of mean squared error (MSE) based on the testing data, whereŷ i is the predicted values for test data; y i is the actual values for test data; n is the total number of test data. The MSE was calculated and is shown in Tables 6 and 7, depending on the algorithms utilised. According to Tables 6 and 7     Its neural network is shown in Figure 10. This training process required 3728 steps until all absolute partial derivatives were smaller than 0.01. The error of training process is 0.1447.

Box-Behnken analysis of top surface
The ANOVA of the top surface roughness is shown in Table 8 Fig. 11. Normal probability plot of residuals (top surface roughness).  The prediction equation for top surface roughness is shown in equation (8).

See equation (8) above.
The normal probability plot in Figure 11 signifies that the residual from the prediction model is normally distributed and the above model can be efficiently used to predict the response.
According to Figure 12, there is an interaction between bed temperature (A) and layer thickness (C). The interaction plot shows that the layer thickness should be set at the lowest level (0.1 mm), so the top surface roughness would be minimised.
The interaction plot (Fig. 13) between printing speed (B) and layer thickness (C) indicates that the printing speed should be higher than 70 mm/s and the layer thickness should be as low as possible, so the top surface roughness could be minimised.
A workpiece was fabricated at a different condition, according to Table 9. The result shows that the top surface roughness measured for this piece was not significantly different from those determined from the prediction equation. As a result, the proposed models are capable of predicting the top surface roughness of the workpieces.

ANN analysis of top surface
Neural networks with different numbers of neurons in each layer were fitted to the data. Their prediction errors were calculated and are shown in Tables 10 and 11. According to Table 10, the application of the RPROP+ algorithm led to the lowest MSE (0.163) while S = {3, 5, 4, 1}. When the GRPROP algorithm was applied with S = {3, 5, 3, 1}, the lowest MSE was equal to 0.147. The ANN structure is illustrated in Figure 14. For the training process, 1600 steps were required before the partial derivatives of error functions was below the threshold value (0.01). The error of training process is 0.106.

Conclusions
The surface quality of workpieces fabricated from the FFF system is the primary focus of this study. The prediction models of surface roughness were determined from two different methods, a response surface method (Box-Behnken design) and machine learning (ANN). The experiment is based on the design matrix retrieved from the Box-Behnken design. Since there are three inputs and a response, there are sixteen treatments with four centre points. Each run consists of five blocks, so there is a total of 80 runs. The measurement positions are the bottom and top surface. One important research question is the performance capability of the Box-Behnken and ANN method. Another question to be investigated is the prediction performance of ANN, based on different attributes, training algorithms (RPROP+ or GRPROP) and the number of neurons in the hidden layers. The performance index used comprises the prediction errors, in terms of MSE.
According to Table 12, the prediction errors of each method (Box-Behnken and ANN) were recapped and compared, and the results show that the prediction capability of ANN outperforms the Box-Behnken method, although the number of data is limited to only 80 for both bottom and top surface roughness. When further analysis was conducted, the results showed that GRPROP delivers a promising result and has a better performance than RPROP+ (bottom Ra: 0.1598 VS 0.2094, top Ra: 0.147 VS 0.163). The number of neurons in the first hidden layer should be set at 5 and the number of neurons in the second layer are between 3 and 4. Based on these findings, both hobbyists and industrial users will be able to use the fabricating solutions so they can operate the FFF system to achieve the best surface finish. The application of ANN to accurately predict the surface roughness is also recommended in this study.