extfieldsofstudy:
Computer Science
Medicine

text:
High Three-Dimensional Detection Accuracy in Piezoelectric-Based Touch Panel in Interactive Displays by Optimized Artificial Neural Networks

High detection accuracy in piezoelectric-based force sensing in interactive displays has gained global attention. To achieve this, artificial neural networks (ANN)—successful and widely used machine learning algorithms—have been demonstrated to be potentially powerful tools, providing acceptable location detection accuracy of 95.2% and force level recognition of 93.3% in a previous study. While these values might be acceptable for conventional operations, e.g., opening a folder, they must be boosted for applications where intensive operations are performed. Furthermore, the relatively high computational cost reported prevents the popularity of ANN-based techniques in conventional artificial intelligence (AI) chip-free end-terminals. In this article, an ANN is designed and optimized for piezoelectric-based touch panels in interactive displays for the first time. The presented technique experimentally allows a conventional smart device to work smoothly with a high detection accuracy of above 97% for both location and force level detection with a low computational cost, thereby advancing the user experience, and serviced by piezoelectric-based touch interfaces in displays.

Introduction
Touch-based interactivity has become a must-have function in smartphones and has created abundant applications in the last decade, improving users human-machine interactivity (HMI) in a convenient and highly efficient manner. The two-dimensional touch sensing, implemented by capacitive and resistive architectures, has traditionally dominated the market. However, the rapid development of information technology requires more data to be exchanged between the user and the end-terminal nowadays, boosting the popularity of three-dimensional force touch sensing.
Current force sensing in commercialized smartphones is supported by capacitive (e.g., iPhone 6S) and piezoresistive (e.g., iPhone X) techniques. The former integrates a layer of capacitive sensors into the backlight of the display to measure the distance shift due to the applied force between the cover glass and the backlight. In contrast, the latter utilizes the force-induced resistance change at the inserted electrodes to interpret the force level. To obtain the force sensing functionality without affecting the capacitive touch sensing, both the capacitive and the piezoresistive techniques require

Literature Review
A typical structured piezoelectric force touch panel is shown in Figure 1. It can be observed that three layers are used to construct the touch panel: A thin piezoelectric film layer in the middle for responding to force touch events and two electrode layers, one on top and one underneath the piezoelectric layer, for conveying the force-induced electric signals [10][11][12][13][14][15]. When a force touch is applied at the surface of the piezoelectric touch panel, the polarization of the piezoelectric layer becomes positive or negative according to the direction of the applied force, hence attracting the charges at the electrodes to form an electric potential. The relationship between the applied force and the induced polarization can be expressed as follows [16,17]: where P i is the induced polarization, and σ jk and d ijk denote the stress and piezoelectric coefficient, respectively. The coefficients remain the same for direct and inverse piezoelectric effects. The coefficients d ijk are symmetric in j and k [10]. Thus, Equation (1) can be simplified as follows: P i = d ij σ j with i = 1, 2, 3 and j = 1, 2, . . . , 6 (2) and systems. The work presented here reveals the relationship between ANN structure, detection accuracy, and computational cost of the piezoelectric-based force touch panel for the first time, providing design considerations for both the academic and the industrial community. This paper is organized as follows: Section 2 briefly reviews the working principle of piezoelectric-based force touch panels and their limitations when being used as interactive interfaces for smartphones; Section 3 explains the methodology; Section 4 demonstrates and discusses the experimental results. A typical structured piezoelectric force touch panel is shown in Figure 1. It can be observed that three layers are used to construct the touch panel: A thin piezoelectric film layer in the middle for responding to force touch events and two electrode layers, one on top and one underneath the From Equations (1) and (2), one can observe that the induced polarization is determined by both the piezoelectric coefficient and the stress vector. The piezoelectric coefficient does not change considerably after fabrication [18]. However, the same force touch can induce different stress vectors when applied at different locations of the touch panel, which is due to the mechanical properties and boundary conditions of the touch panel. As a result, the force-polarization responsivity is not consistent, which downgrades the detection accuracy [8,9,19].

Methodology
In this section, the methodology used in this work is described. Details of the experimental setup and data acquisition technique are provided, followed by a description of the pre-processing method and the datasets for training and validation. Finally, the ANN for classification is presented.

Experimental Setup and Data Acquisition
The experiment is carried out with a previously developed multi-layered prototype utilized in Reference [19], consisting of a top cover layer, nine evenly distributed electrodes, a piezoelectric thin film layer, a ground electrode layer, and a bottom cover layer. Figure 2 depicts the topology, together with its geometries. piezoelectric layer, for conveying the force-induced electric signals [10][11][12][13][14][15]. When a force touch is applied at the surface of the piezoelectric touch panel, the polarization of the piezoelectric layer becomes positive or negative according to the direction of the applied force, hence attracting the charges at the electrodes to form an electric potential. The relationship between the applied force and the induced polarization can be expressed as follows [16,17]: where Pi is the induced polarization, and σjk and dijk denote the stress and piezoelectric coefficient, respectively. The coefficients remain the same for direct and inverse piezoelectric effects. The coefficients dijk are symmetric in j and k [10]. Thus, Equation (1) can be simplified as follows: From Equations (1) and (2), one can observe that the induced polarization is determined by both the piezoelectric coefficient and the stress vector. The piezoelectric coefficient does not change considerably after fabrication [18]. However, the same force touch can induce different stress vectors when applied at different locations of the touch panel, which is due to the mechanical properties and boundary conditions of the touch panel. As a result, the force-polarization responsivity is not consistent, which downgrades the detection accuracy [8,9,19].

Methodology
In this section, the methodology used in this work is described. Details of the experimental setup and data acquisition technique are provided, followed by a description of the pre-processing method and the datasets for training and validation. Finally, the ANN for classification is presented.

Experimental Setup and Data Acquisition
The experiment is carried out with a previously developed multi-layered prototype utilized in Reference [19], consisting of a top cover layer, nine evenly distributed electrodes, a piezoelectric thin film layer, a ground electrode layer, and a bottom cover layer. Figure 2 depicts the topology, together with its geometries. Four volunteers to perform force touch events are employed for the experiment, and their physical body conditions are listed in Table 1. Each of them is required to tap 1350 times at nine specified locations (as shown in Figure 2) of the touch panel with three different force levels. Hence, a touch set is defined as tapping 50 times per location per force level. Because of the difference in physical body conditions and touch habits of the volunteers, the strength of the three force levels are highly dependent on the individual [19][20][21]. Four volunteers to perform force touch events are employed for the experiment, and their physical body conditions are listed in Table 1. Each of them is required to tap 1350 times at nine specified locations (as shown in Figure 2) of the touch panel with three different force levels. Hence, a touch set is defined as tapping 50 times per location per force level. Because of the difference in physical body conditions and touch habits of the volunteers, the strength of the three force levels are highly dependent on the individual [19][20][21]. The voltages across the nine piezoelectric-based capacitors are continuously monitored. The retrieved raw data for one touch subset are given in Figure 3 as an example, including the raw data from the nine channels; each channel illustrates the relationship between voltage and time. Each peak in the raw data indicates a force touch event. Since the force-induced stress at one location can propagate to other areas, all nine channels produce an output signal when a touch event occurs. However, the amplitudes of the signals are different and depend on the stress distribution on the touch panel. Furthermore, as a complete force touch event consists of "press" and "release" actions, the polarizations of the piezoelectric layer change twice, resulting in positive and negative peaks, as illustrated in Figure 3. The "press" action induces a positive electric signal and is considered to represent the force level of the force touch event.  The voltages across the nine piezoelectric-based capacitors are continuously monitored. The retrieved raw data for one touch subset are given in Figure 3 as an example, including the raw data from the nine channels; each channel illustrates the relationship between voltage and time. Each peak in the raw data indicates a force touch event. Since the force-induced stress at one location can propagate to other areas, all nine channels produce an output signal when a touch event occurs. However, the amplitudes of the signals are different and depend on the stress distribution on the touch panel. Furthermore, as a complete force touch event consists of "press" and "release" actions, the polarizations of the piezoelectric layer change twice, resulting in positive and negative peaks, as illustrated in Figure 3. The "press" action induces a positive electric signal and is considered to represent the force level of the force touch event.  To obtain the information carried in the press-induced electric signals, a pre-processing algorithm has been developed through four steps to extract the voltage peaks from each channel, as

Pre-Processing and Dataset Preparation
To obtain the information carried in the press-induced electric signals, a pre-processing algorithm has been developed through four steps to extract the voltage peaks from each channel, as shown in Figure 4. First, the direct current (DC) offset of each channel is removed by subtracting the mean voltage of each channel. Second, envelope detection is conducted to each channel by the Hilbert Transform method, as the voltage variations are both positive and negative. Third, we assume that the voltage variation peaks of different channels appear at the same time and calculate the mean value of the nine different channels to suppress random noise. Finally, the voltage peaks are detected by finding the local maximum voltage values. After the voltage peaks for each channel have been acquired, two ANNs are employed to classify locations and force levels. Their details will be given in the next sub-section.
After pre-processing, nine voltage peaks for each force touch event are obtained and used as features to classify the touch locations and force levels. As mentioned before, 50 touches are performed at each location at the same force amplitude, so 1350 touch signals are obtained in total from each volunteer. Since one touch provides one output data point (including nine peak values from all nine channels), the dataset for a single volunteer consists of 50·9·3 output data points. Then, for the training, validation, and testing of the ANN, the dataset is further separated into a training set, a validation set, and a test set. The training set has 40 force touch events for each location and each force level; hence, it has 40·9·3 data points, which makes up 80% of the data; the validation and test sets have 10 events each (5·9·3 data points, which makes up 10% each).  To obtain the information carried in the press-induced electric signals, a pre-processing algorithm has been developed through four steps to extract the voltage peaks from each channel, as

Pre-Processing Multi-Layers Neural Networks for Classification
Convolution neural networks (CNN) [22][23][24][25][26][27][28][29] and recurrent neural networks (RNN) [30][31][32][33][34][35] are conventionally used to extract features for image and audio applications, and a fully connected network [36][37][38] can be used as an output layer of a CNN or RNN. However, a fully connected network can also be used alone for classification. In our case, the signal features have already been extracted with the help of the pre-processing step explained above, and, thus, we only employed the fully connected network for location and force level classification. shown in Figure 4. First, the direct current (DC) offset of each channel is removed by subtracting the mean voltage of each channel. Second, envelope detection is conducted to each channel by the Hilbert Transform method, as the voltage variations are both positive and negative. Third, we assume that the voltage variation peaks of different channels appear at the same time and calculate the mean value of the nine different channels to suppress random noise. Finally, the voltage peaks are detected by finding the local maximum voltage values. After the voltage peaks for each channel have been acquired, two ANNs are employed to classify locations and force levels. Their details will be given in the next sub-section. After pre-processing, nine voltage peaks for each force touch event are obtained and used as features to classify the touch locations and force levels. As mentioned before, 50 touches are performed at each location at the same force amplitude, so 1350 touch signals are obtained in total from each volunteer. Since one touch provides one output data point (including nine peak values from all nine channels), the dataset for a single volunteer consists of 50•9•3 output data points. Then, for the training, validation, and testing of the ANN, the dataset is further separated into a training set, a validation set, and a test set. The training set has 40 force touch events for each location and each force level; hence, it has 40•9•3 data points, which makes up 80% of the data; the validation and test sets have 10 events each (5•9•3 data points, which makes up 10% each).

Pre-Processing Multi-Layers Neural Networks for Classification
Convolution neural networks (CNN) [22][23][24][25][26][27][28][29] and recurrent neural networks (RNN) [30][31][32][33][34][35] are conventionally used to extract features for image and audio applications, and a fully connected network [36][37][38] can be used as an output layer of a CNN or RNN. However, a fully connected network can also be used alone for classification. In our case, the signal features have already been extracted with the help of the pre-processing step explained above, and, thus, we only employed the fully connected network for location and force level classification. The fully connected ANN has three kinds of layers with multiple nodes embedded, including an input layer, an output layer, and multiple hidden layers ( Figure 5). The number of the input and output layers' nodes are identical to the number of features and object classes, which are nine and nine for location classification, and nine and three, respectively, for force level classification. Alternatively, the number of hidden layers and their embedded nodes are two independent hyperparameters, which have a significant effect on variance and computational complexity. In general, an increase in the value of these two numbers will have a negative effect on the computational cost; however, ANN becomes more powerful to deal with complex classification scenarios. Thus, the tradeoff between the performance and computational cost must be taken into careful account when designing a fully connected ANN. The fully connected ANN has three kinds of layers with multiple nodes embedded, including an input layer, an output layer, and multiple hidden layers ( Figure 5). The number of the input and output layers' nodes are identical to the number of features and object classes, which are nine and nine for location classification, and nine and three, respectively, for force level classification. Alternatively, the number of hidden layers and their embedded nodes are two independent hyperparameters, which have a significant effect on variance and computational complexity. In general, an increase in the value of these two numbers will have a negative effect on the computational cost; however, ANN becomes more powerful to deal with complex classification scenarios. Thus, the tradeoff between the performance and computational cost must be taken into careful account when designing a fully connected ANN. The computational time cost includes training time and classification time costs, both of which increase as the layer number and node number increase. The training time cost only occurs once, and that is during the initial training of the ANN. However, the classification time cost occurs every time the user presses the touch panel. Therefore, a high tolerance is usually set for the training time cost, and a low tolerance for the classification time cost. In this study, the training time cost is defined experimentally, while the classification time cost is analyzed theoretically. We implemented our ANNs with different layer numbers and node numbers using Keras. Regarding the training time cost, we trained the ANNs using a laptop (Lenovo ThinkPad X1) with a CPU (Intel(R) Core(TM) i5-3427U @ 1.8GHz), and the time cost of training was recorded. With regard to the classification time cost, it is positively correlated to the number (W) of total weight parameters, which can be calculated by: where F is the feature number (nine for both force level and location classification), M and N are the node number and layer number, respectively, and C is the number of object classes. The cost function used for training can also greatly affect the performance [39][40][41][42][43][44]. Three commonly used cost functions have been chosen: The mean-squared-error loss function, the categorical cross-entropy loss function, and the binary cross-entropy loss function, and their classification accuracy was compared.
For the minimization of the loss function during the training process, a stochastic gradient-based optimizer (Adam 2015 [45]) is used to optimize the values of weights and bias of the network, due to its computational efficiency and its low memory requirements [46]. The active function used for the input and hidden layers is 'ReLU', since it has high time efficiency [47,48]. In order to avoid overfitting, we make use of L2 regularization with a regularization factor of 0.001, which is selected experimentally [49]. The learning rate is experimentally selected as 0.0005.

Results and Discussion
The pre-processing results after each step are plotted in Figure 6. One can observe that the DC offset at the original output signal is about 0.8 V (Figure 6a). After successfully removing the DC offset, the mean of the signal is close to 0 V (Figure 6b). The envelope detection algorithm is then applied to the offset-free signal, since there are positive and negative components during a force touch event (Figure 6c). Based on the assumption that voltage variation peaks of different channels appear at the same time, the mean value of the nine channels is calculated, as Figure 6d depicts. Finally, the peaks are located and marked by circles in the detected peaks figure. As Figure 6e confirms, the proposed pre-processing method is capable of appropriately detecting the peaks of the signal. Figure 7 shows the response voltages from one volunteer at the same location with different force levels. As shown, our prototype touch panel has different ranges of response voltage for different force levels.
After the pre-processing, the training of the force level classification ANN is performed using three different cost functions, the mean square error cost function, the categorical cross-entropy cost function, and the binary cross-entropy cost function. The corresponding accuracy and loss results are provided in Figure 8. Based on this figure, one can not only observe that the accuracies of both the training set and the validation set are identical, but also that their losses are close to each other. Hence, it is demonstrated that there is no obvious overfitting in this training. Among these three cost functions, the binary cross-entropy cost function provides the best accuracy (98.5%), hence, it is chosen as the cost function for force level classification.
network, due to its computational efficiency and its low memory requirements [46]. The active function used for the input and hidden layers is 'ReLU', since it has high time efficiency [47,48]. In order to avoid overfitting, we make use of L2 regularization with a regularization factor of 0.001, which is selected experimentally [49]. The learning rate is experimentally selected as 0.0005. The pre-processing results after each step are plotted in Figure 6. One can observe that the DC offset at the original output signal is about 0.8 V (Figure 6a). After successfully removing the DC offset, the mean of the signal is close to 0 V (Figure 6b). The envelope detection algorithm is then applied to the offset-free signal, since there are positive and negative components during a force touch event (Figure 6c). Based on the assumption that voltage variation peaks of different channels appear at the same time, the mean value of the nine channels is calculated, as Figure 6d depicts. Finally, the peaks are located and marked by circles in the detected peaks figure. As Figure 6e confirms, the proposed pre-processing method is capable of appropriately detecting the peaks of the signal. Figure 7 shows the response voltages from one volunteer at the same location with different force levels. As shown, our prototype touch panel has different ranges of response voltage for different force levels.   The three cost functions used above are also applied for training of the location classification ANN, and the results are shown in Figure 9. Unlike force level classification, over-fitting appears in two of the three cost functions: The mean-squared error function and categorical cross-entropy cost function. The over-fittings can be seen by observing that the accuracies of the training sets are higher than those of the validation sets. In contrast, the accuracy of the binary cross-entropy cost function for the training set is close to that of the validation set (Figure 9e). In fact, it is 97.8% higher than that of their counterparts. Therefore, the binary cross-entropy cost function is selected for location classification.   After the pre-processing, the training of the force level classification ANN is performed using three different cost functions, the mean square error cost function, the categorical cross-entropy cost function, and the binary cross-entropy cost function. The corresponding accuracy and loss results are provided in Figure 8. Based on this figure, one can not only observe that the accuracies of both the training set and the validation set are identical, but also that their losses are close to each other. Hence, it is demonstrated that there is no obvious overfitting in this training. Among these three cost functions, the binary cross-entropy cost function provides the best accuracy (98.5%), hence, it is chosen as the cost function for force level classification. The three cost functions used above are also applied for training of the location classification ANN, and the results are shown in Figure 9. Unlike force level classification, over-fitting appears in two of the three cost functions: The mean-squared error function and categorical cross-entropy cost function. The over-fittings can be seen by observing that the accuracies of the training sets are higher than those of the validation sets. In contrast, the accuracy of the binary cross-entropy cost function for the training set is close to that of the validation set ( Figure 9e). In fact, it is 97.8% higher than that of their counterparts. Therefore, the binary cross-entropy cost function is selected for location classification.  Figure 10 shows the correlation between the achieved accuracy and the fully connected ANN structure for different layer numbers and node numbers. Figure 10a shows the accuracy of force level classification. For layer numbers 1 to 3 and node numbers 4 to 32, the accuracy increases as the layer number or node number increases. The reason for this is that the ANN structure becomes more complex, hence, it can fit complex classification functions. However, when the layer number becomes larger than 3, and the node number greater than 32, no significant change in the accuracy is observed. Therefore, a layer number of 3 and a node number of 32 are selected for the force level classification in this paper. Figure 10b shows the accuracy of the location classification, which is slightly lower than that of the force level classification. As the layer number increases from 1 to 4 and the node number from 4 to 64, the accuracy also increases. However, for a layer number larger than 3 and a node number larger than 32, the accuracy does not change significantly. Thus, a layer number of 3 and a node number of 64 are chosen for location classification in this paper. Then, the accuracy of the test set is calculated as 97.7% and 97.0% for force and location classification, respectively.

Results and Discussion
The three cost functions used above are also applied for training of the location classification ANN, and the results are shown in Figure 9. Unlike force level classification, over-fitting appears in two of the three cost functions: The mean-squared error function and categorical cross-entropy cost function. The over-fittings can be seen by observing that the accuracies of the training sets are higher than those of the validation sets. In contrast, the accuracy of the binary cross-entropy cost function for the training set is close to that of the validation set ( Figure 9e). In fact, it is 97.8% higher than that of their counterparts. Therefore, the binary cross-entropy cost function is selected for location classification.  Figure 10 shows the correlation between the achieved accuracy and the fully connected ANN structure for different layer numbers and node numbers. Figure 10a shows the accuracy of force level classification. For layer numbers 1 to 3 and node numbers 4 to 32, the accuracy increases as the layer number or node number increases. The reason for this is that the ANN structure becomes more complex, hence, it can fit complex classification functions. However, when the layer number becomes larger than 3, and the node number greater than 32, no significant change in the accuracy is The training time costs for location and force level classification ANNs are given in Figure 11. The time cost remarkably increases when the node number is greater than 64. Also, it rises as the layer number increases. When the ANN has 5 hidden layers and 512 nodes, the maximum training time cost for location and force level classification are 418. 2  observed. Therefore, a layer number of 3 and a node number of 32 are selected for the force level classification in this paper. Figure 10b shows the accuracy of the location classification, which is slightly lower than that of the force level classification. As the layer number increases from 1 to 4 and the node number from 4 to 64, the accuracy also increases. However, for a layer number larger than 3 and a node number larger than 32, the accuracy does not change significantly. Thus, a layer number of 3 and a node number of 64 are chosen for location classification in this paper. Then, the accuracy of the test set is calculated as 97.7% and 97.0% for force and location classification, respectively. The training time costs for location and force level classification ANNs are given in Figure 11. The time cost remarkably increases when the node number is greater than 64. Also, it rises as the layer number increases. When the ANN has 5 hidden layers and 512 nodes, the maximum training time cost for location and force level classification are 418.2 and 424.8 ms/epoch, respectively.  Figure 12 gives the parameter number for force classification ANN, which is shown on a logarithmic scale, since it grows exponentially. The parameter number for the location classification ANN is a little larger than that of the force level classification ANN, because the output layer node number increases from three-nine. The classification time cost increases as the parameter number increases. Thus, the classification cost also grows exponentially. The location classification time measured for one touch is 0.98 ms for an ANN with 5 hidden layers and 512 nodes for each layer. The training time costs for location and force level classification ANNs are given in Figure 11. The time cost remarkably increases when the node number is greater than 64. Also, it rises as the layer number increases. When the ANN has 5 hidden layers and 512 nodes, the maximum training time cost for location and force level classification are 418.2 and 424.8 ms/epoch, respectively.  Figure 12 gives the parameter number for force classification ANN, which is shown on a logarithmic scale, since it grows exponentially. The parameter number for the location classification ANN is a little larger than that of the force level classification ANN, because the output layer node number increases from three-nine. The classification time cost increases as the parameter number increases. Thus, the classification cost also grows exponentially. The location classification time measured for one touch is 0.98 ms for an ANN with 5 hidden layers and 512 nodes for each layer.

Conclusions
ANNs have been used in the literature to address the inconsistent force-voltage responsivity issue in piezoelectric-based force touch panels for interactive displays. However, the high computational complexity and the relatively low detection accuracy limit their successful use in conventional consumer products. The work presented in this paper bridges this gap by investigating the relationship between the ANN hyperparameters, detection accuracy, and computational cost. After carefully designing the ANN on the basis of the specific characteristics of the piezoelectric-based touch events, both a high detection accuracy of above 97% and a low computational cost of less than 0.98 ms have been achieved experimentally.


Performance of the Belle II Conditions Database

The Belle II experiment at KEK observed its first collisions in the summer of 2018. Processing the large amounts of data that will be produced requires conditions data to be readily available to systems worldwide in a fast and efficient manner that is straightforward for both the user and maintainer. This was accomplished by relying on industrystandard tools and methods: the conditions database is built as an HTTP REST service using tools such as Swagger for the API interface development, Payara for the Java EE application server, and Squid for the caching proxy. This article presents the design of the Belle II conditions database production environment as well as details about the capabilities and performance during both Monte Carlo campaigns and data reprocessing. 1 The Belle II Experiment The Belle II experiment is part of a broad-based search for new physics in the Intensity Frontier, focused on precisely measuring and comparing with theory branching fractions, angular distributions, CP asymmetries, forward-backward asymmetries, and a host of other observables. It is the successor to the highly-successful Belle experiment, which provided the experimental foundation for the 2008 Nobel Prize in Physics in CP violation [1]. The SuperKEKB accelerator upgrade will provide asymmetric electron-positron beams tuned to produce large quantities of B/anti-B meson pairs at 40x the instantaneous luminosity of KEKB and 50x the data taken with Belle. Initial “Phase 2” data-taking started in May 2018, concurrent with beam-tuning efforts, and “Phase 3” full running starts in 2019. Figure 1 presents the various interactions that the conditions database (DB) needs to support. The primary purpose of the conditions database is to support data processing across the international Belle II computing grid, shown as the WAN connection at lower left that connects to the KEK computing center and Tier 1 computing centers, as well as interactive users. Data enters the database via the calibration cycle, which uses data collected by the Belle II DAQ at KEK and processes it into the appropriate format for the * Corresponding author: lynn.wood@pnnl.gov © The Authors, published by EDP Sciences. This is an open access article distributed under the terms of the Creative Commons Attribution License 4.0 (http://creativecommons.org/licenses/by/4.0/). EPJ Web of Conferences 214, 04050 (2019) https://doi.org/10.1051/epjconf/201921404050 CHEP 2018

The Belle II Experiment
The Belle II experiment is part of a broad-based search for new physics in the Intensity Frontier, focused on precisely measuring and comparing with theory branching fractions, angular distributions, CP asymmetries, forward-backward asymmetries, and a host of other observables. It is the successor to the highly-successful Belle experiment, which provided the experimental foundation for the 2008 Nobel Prize in Physics in CP violation [1]. The SuperKEKB accelerator upgrade will provide asymmetric electron-positron beams tuned to produce large quantities of B/anti-B meson pairs at 40x the instantaneous luminosity of KEKB and 50x the data taken with Belle. Initial "Phase 2" data-taking started in May 2018, concurrent with beam-tuning efforts, and "Phase 3" full running starts in 2019. Figure 1 presents the various interactions that the conditions database (DB) needs to support. The primary purpose of the conditions database is to support data processing across the international Belle II computing grid, shown as the WAN connection at lower left that connects to the KEK computing center and Tier 1 computing centers, as well as interactive users. Data enters the database via the calibration cycle, which uses data collected by the Belle II DAQ at KEK and processes it into the appropriate format for the conditions database. Currently this is done offline, but is envisioned to become significantly more automated as Belle II matures.

Fig. 1.
Interactions of the Belle II conditions database with production nodes, interactive users, and the Belle II data acquisition system.

Concepts
A conditions database holds time-dependent status of detectors for data processing and reprocessing, consisting of either constants or time-varying parameters such as calibration settings, geometry, alignment, etc.
User access to the Belle II database is through a REST API, meaning that all accesses take the form of HTTP requests. For example, the following HTTP request would return a list of payloads for a given global tag, experiment and run number: This type of interface provides several key advantages to Belle II: it restricts DB access to allowed operations only, and makes it easy to scale with industry standard (HTTP) tools. Responses can be formatted in either XML or JSON format.
Payloads are kept as external files outside of the database; a request for a payload returns a partial URL of the payload file, which can be combined with a hostname to access a file server containing the payload. This method avoids a database bottleneck for large payload files, and allows for different scenarios for file storage and transfer.

Global Tags and Intervals of Validity
The data in the Belle II conditions database is organized using several structures. Intervals of Validity (IOVs) specify the starting and ending experiments and runs for which a payload is valid. A given payload may have multiple IOVs assigned to it if it is valid for different periods of time.
A Global Tag contains a list of IOV-payload relationships, and is used to select a complete set of conditions for a given reprocessing effort. The IOVs and payloads can be reused in different global tags, and the global tags can be classified by type (development, release, online, data, etc.) and have evolving status over time (new, published, invalid, etc.).

Service Layers
When operating at the expected full processing capacity, the Belle II computing grid will generate ~100 database requests per second, and can burst to significantly higher rates. To attain the required performance, various features of commercial tools designed for scalability are used.
As examples: Squid HTTP caches [2] are used to support repeated requests for the same query. These are configured as reverse proxy, i.e. optimized for many clients and few servers. These serve to cache the most common global tags in Belle II. Requests that are not yet found in the Squid caches reach the "b2s" Belle II service layer. This layer is implemented in Java using Payara Micro [3], a microservice JavaEE server, which translates REST requests into SQL queries. The REST API is built using Swagger [4], which can auto-generate client and server code, as well as provide interactive online documentation of the API. The storage back-end behind these layers is the PostgreSQL [5] relational database, which resides on the IBM General Parallel File System (GPFS) [6].

Belle II Database Server Architecture
As mentioned previously, the database and payloads are kept on separate servers for improved performance and support. Both servers are built with multiple instances of all components, using Docker [7] containers to provide easy support, greater reliability, and improved performance. As shown on the top in Figure 2, the database service consists of multiple nodes behind an HTTP load balancer. Each node has its own dedicated Squid cache and multiple instances of the b2s service which has a connection to query the single database service. The database service contains both a primary database as well as a hot standby backup database to minimize any downtime.
As shown at the bottom in Figure 2, the payload file is designed similarly, with multiple nodes behind an HTTP load balancer. Each node contains multiple instances of the NGINX [8] high-performance HTTP server, and the payload files themselves are stored in the IBM GPFS accessible by all nodes.
The primary interface is through the Belle II Analysis Software Framework [9], known as basf2. To simplify efforts, several assumptions are made about the conditions data. First, payload IOV granularity is assumed to be at the run level, as opposed to time or event granularity. This is expected to be sufficient for the majority of the conditions data in Belle II, but if sub-run granularity is required basf2 provides support for storing event numberbased granularity dependence inside of a single run payload. Another assumption made is that all payloads are stored as ROOT [10] objects. This is not a requirement of the database system itself since the payloads are simply files stored on a file server, but makes the analysis software more manageable.
The DBStore service inside of basf2 manages storage and automatic updates of payloads as data is processed. This is handled through C++ template classes that always provide the user with a pointer to the currently-correct payload. The DBStore service ensures that the pointer is updated when a new IOV becomes valid for a payload, making this rundependence transparent to the user. basf2 also supports accessing IOVs and payloads from either a remote database or a local directory, allowing development and testing of new conditions data payloads before loading them into the global database. This is handled by creating a "chain" of database locations to search for a given payload in the basf2 steering file. The relationship between the conditions database components in basf2 is shown in Figure 3.
Due to the easy accessibility of the HTTP API, several other interfaces are also available. The Swagger tool used to develop the API provides online documentation of all queries, and allows the user to interactively send queries and view responses. A command-line client written in Python [11] is also available, and is currently the primary method of adding new payloads to the database.

Service Layers
The Belle II conditions database is a sophisticated system with many interactions between the different service layers. Extensive monitoring of all layers is implemented using Prometheus [12] and Grafana [13]; example parameters to be monitored include CPU and RAM usage, disk latency, input HTTP request rate, HTTP response types, and so on.

Database Performance
The current Monte Carlo campaigns and the data reprocessing for data taken in 2018 do not put a significant load on the DB and payload file servers. The current biggest stressor is users submitting large numbers of non-grid-based jobs, which cause a significant spike in requests in a short period of time. To ensure that the database is performing to the required level, directed tests are necessary to confirm that performance is sufficient for the expected full processing load.
Multiple types of tests are performed to ensure proper database behavior via Monte Carlo job submissions to the grid, scripting, and the Gatling HTTP load tester [14]. These tests include verification of basf2 software releases as well as specific scripts to generate "problem" workflows, such as the non-grid user submissions mentioned above.
In one example of performance testing, the Gatling tool was configured to generate 400 requests per second for a period of 15 minutes. The requests were randomly selected from a list of 21,000 actual request types recorded from HTTP logs. As the test ran, the network bandwidth and CPU usage was monitored, as well as the actual rates of requests submitted and responses received. No performance issues were seen during this test, and all parameters showed steady values, indicating sufficient resources were available for this type of workload.

Upcoming Development
While already in use by the Belle II experiment, a number of improvements are planned and in progress for the conditions database.
Authentication and authorization is currently not present for database access. The expectation is that all reads will be open, but write access needs to be user authenticated and authorized. Several methods of implementing this functionality are being investigated, including leveraging the existing X.509 certification authentication used for Belle II grid processing. It has been proposed to support three roles for database users: "user", which would allow read-only access; "developer", which would support adding data to existing global tags; and "coordinator", who can create and manipulate global tags. These access rights would also vary with the type of global tag, to better control the update process.
While the implemented servers are fully capable of supporting Belle II during full production, the loss of networking access to the single location at Brookhaven National Laboratory would prevent proper operation of data reprocessing. To this end, replication of the system to other sites is being investigated. An intermediate solution being used currently is replicating the database contents and payload files on a distributed file system (CVMFS [15]), but as the database contents grow that may become a burden on the distributed grid site nodes. PostgreSQL supports standard replication such as streaming or log-based methods, although replication to remote sites is more complicated. Because of this, persistent remote Squid caching is being considered as a solution as well, where long cache invalidation times are implemented.


Data-driven monitoring of the gearbox using multifractal analysis and machine learning methods

Data-driven diagnostic methods allow to obtain a statistical model of time series and to identify deviations of recorded data from the pattern of the monitored system. Statistical analysis of time series of mechanical vibrations creates a new quality in the monitoring of rotating machines. Most real vibration signals exhibit nonlinear properties well described by scaling exponents. Multifractal analysis, which relies mainly on assessing local singularity exponents, has become a popular tool for statistical analysis of empirical data. There are many methods to study time series in terms of their fractality. Comparing computational complexity, a wavelet leaders algorithm was chosen. Using Wavelet Leaders Multifractal Formalism, multifractal parameters were estimated, taking them as diagnostic features in the pattern recognition procedure, using machine learning methods. The classification was performed using neural network, k-nearest neighbours’ algorithm and support vector machine. The article presents the results of vibration acceleration tests in a demonstration transmission system that allows simulations of assembly errors and teeth wear.

Introduction
Monitoring the operation of a machine, we assume that detection, isolation and identification of a damage requires a comparison of the current state to a state taken as a reference. Diagnosis can therefore be considered a process in line with the paradigm of pattern recognition. Properly defined purpose, data collection, data cleaning and normalisation, data transformation and selection of relevant diagnostic features, as well as classification based on the chosen exploration technique are required. In view of the intrinsic dynamic nature and external excitations, most of the real signals coming from the monitoring of mechanical vibrations of rotating machines has nonlinear, non-stationary and multiscale properties both in time and space [1]. Traditional methods of signal analysis in engineering applications are based on the assumption of stationarity. Such methods can not reveal local features in time and frequency domains.
One of more frequently used nonlinear methods is the study of point regularity of time series, represented by the Holder scaling exponent. The estimation of scaling exponents leads to determination of a multifractal spectrum, which is a statistical model of the monitored signal [2]. This fact was used to determine the diagnostic features of the test gear being tested and to identify assembly errors and wear. Small fluctuations of shaft rotational speed occurred during the test, but the proposed method does not require synchronisation and resampling of the data, as in the case of e.g. Fast Fourier Transform (FFT). In the studies of vibration signals in terms of their multifractality, methods of multifractal version of the detrending fluctuation analysis (MF-DFA) [3 -9] and various methods based on wavelet analysis [10 -12] are applied.
The article analyses the possibility of using wavelet leaders to select features and neural network, k-nearest neighbours and support vector machines in the process of fault classification of recorded data from the pattern of the monitored system and inferring the occurrence of damage.

Wavelet leaders multifractal formalism
This paragraph presents the idea of quantities called wavelet leaders, which are determined on the basis of wavelet coefficients and are representative of Holder exponents [13 -15]. Multifractal analysis based on leaders is not dependent on wavelet selection, provided that the number of zero wavelet moments is greater than the largest Holder exponent of the signal.
The wavelet leaders were calculated based on the discrete wavelet transform (DWT). The point Holder exponent of the function/signal x(t) in point t denotes the number h defined as the supremum of all exponents satisfying, for a certain C > 0, the condition: where ( − 0 ) is a polynomial of < ℎ degree.
Information on variability of regularity x(t) along t can be described by the multifractal spectrum, estimated using wavelet leaders multifractal formalism.
Let us define discrete wavelet transform of the function x(t) and mother wavelet 0 ( ) with the compact support using its coefficients: Moreover, let us introduce the dyadic intervals: and the union of such neighbour intervals: The quantities: illustrated in Figure 1 are referred to as wavelet leaders. If ( ) has a Holder exponent h at a point 0 = 2 : then the corresponding wavelet leaders decay as power laws of scales: The structure functions defined as sample mean estimators for the ensemble averages of q-power of wavelet leaders exhibit power law behaviours with respect to the analysis scale: The scaling exponents ( ) can be expanded using cumulants : ∀ ≥ 1 satisfy the relation: where ( , ) are the ≥ 1 order cumulants of the ln ( , . ). Moreover 1 describes the location of the multifractal spectrum maxima and 2 , 3 describe the level of multifractality: spectrum width and spectrum asymmetry, respectively.
This leads to a Wavelet Leader based Multifractal Formalism Under some mild regularity conditions (ℎ) acts as an upper bound for multifractal spectrum of analysed signal.

Neural Network NN
Artificial neural networks are relatively crude electronic networks of neurons based on the neural structure of the brain. They process records one at a time and learn by comparing their classification of the record (i.e., largely arbitrary) with the known actual classification of the record. The errors from the initial classification of the first record is fed back into the network and used to modify the networks algorithm for further iterations.
In the training phase, the correct class for each record is known (termed supervised training), and output nodes can be assigned correct values. The typical backpropagation network has an input layer, an output layer, and at least one hidden layer. There is no theoretical limit on the number of hidden layers but typically there are just one or two [16,17].

k-Nearest Neighbours k-NN
The k-NN algorithm is among the simplest of all machine learning algorithms. In k-NN classification, the output is a class membership. An object is classified by a majority vote of its neighbours, with the object being assigned to the class most common among its k nearest neighbours (k is a positive integer, typically small). If k = 1, then the object is simply assigned to the class of that single nearest neighbour. The neighbours are taken from a set of objects for which the class is known. This can be thought of as the training set for the algorithm, though no explicit training step is required [18,19].
The training examples are vectors in a multidimensional feature space, each with a class label. The training phase of the algorithm consists only of storing the feature vectors and class labels of the training samples.
In the classification phase, k is a user-defined constant, and an unlabelled vector (a query or test point) is classified by assigning the label which is most frequent among the k training samples nearest to that query point.
A commonly used distance metric for continuous variables is Euclidean distance. In the Euclidean plane, if p = (p1, p2) and q = (q1, q2) then the distance is given by Various metrics can be used to determine the distance, for example Mahalanobis, city block, Minkowski, Chebychev, cosine, correlation, Hamming, Jaccard, Spearman distance.

Support Vector Machines SVM
In machine learning, support vector machines are supervised learning models with associated learning algorithms that analyse data used for classification. Given a set of training examples, each marked as belonging to one or the other two categories, an SVM training algorithm builds a model that assigns new examples to one category or the other, making it a non-probabilistic binary linear classifier [20,21].
A support vector machine constructs a hyperplane or set of hyperplanes in a high-or infinite-dimensional space, which can be used for classification. Intuitively, a good separation is achieved by the hyperplane that has the largest distance to the nearest training-data point of any class (so-called functional margin), since in general the larger the margin the lower the generalisation error of the classifier.
Given a training vector { | ∈ } =1 and an indication vector { | ∈ {−1,1}} =1 , here the training vector contains two different categories and the element yi in the indication vector serves to indicate which category the sample xi in the training vector belongs to. Then, an SVM is trained to find a maximum-margin hyperplane for distinguishing the samples marked with yi = −1 from those marked with yi = 1. If there is no hyperplane which can separate these two categories, a soft margin method can be employed to determine a hyperplane which can separate these two categories, still keeping the most maximum distance to the closest clearly separated points.
The original SVM can be used only for solving a binary classification problem. Nonetheless, the multiclass classification is a frequently faced problem in the real world. To solve this problem, "one-against-one" and "one-against-all" approaches can be proposed [22].

Gearbox fault diagnosis
An experimental system, which photo is in Figure 2, makes up of the drive motor, one-stage planetary gearbox, and hydraulic brake. The main measurement path included B&K DeltaTron type 4519-002 acceleration sensor with a frequency range of 1 Hz-20 kHz, screwed to the bearing housing and DBK04 dynamic signal acquisition card equipped with a set of antialiasing filters cooperating with the main IQtech DaqBook 2005 acquisition module, which transmits data using an Ethernet connection to the computer. Data were registered with sampling frequency 10kHz during about 100 seconds (about 1 million samples).

Fig. 2. Experimental setup.
Three types of faults were simulated on the test bench: misalignment 1/3 o , enlarged clearance 0,4mm and worn gears. Figure 3 displays the vibration waveforms for faultfree and faulty states of the gearbox.  The motor rpm and load were slightly changing during the test, so the signal is not stationary. So, the Fast Fourier Transform (FFT) for signal analysis without resampling is not a good solution. In turn, resampling does not provide sufficient accuracy in on-line analysis. FFT spectra of raw vibration signal (without resampling) for different states of gearbox are shown in Figure 4. Slightly changing rotational speed causes the frequency spectrum to be modulated by variable frequency and difficult to compare.
Then for each gearbox condition, 100 pieces of data were gathered, each of the length of 8192 points. The Multifractal Wavelet Leaders method was adopted to analyse these gearbox vibration data. The multifractal spectra of sample data are demonstrated in Figure 5. Obtained spectra can be treated as multifractal dimensions (ℎ) associated with singularities ℎ representing local scaling of the measure at various places in the time series. The description of the dynamic properties of the system based on the multifractal spectrum of time series is permitted by: the singularity with the largest dimension, i.e. the most often met the time series singularity and dimensions span of the singularities subsets. These parameters are represented by first and second cumulants, respectively.
Next, for each piece of four gearbox condition, the cumulants were calculated. For classification purpose first two cumulants were chosen. Figure 6 shows a scatter plot second cumulant versus first cumulant. Since the scratch teeth and enlarged clearance states can be easily separated, separation of fault-free and misalignment states may be more complicated. Three methods of classification were compared for separation of the states: neural networks (NN), nearest neighbours (kNN) and Support Vector Machine (SVM) using Matlab environment.
Each of four sets registered for different fault states was divided into a training set (70% of samples) and a test set (30% of samples) and then subjected to machine learning. Using NN method a two-layer, feed-forward network with 10 sigmoid hidden neurons and linear output neurons was performed. The network was trained with Levenberg-Marquardt back propagation algorithm. In the k-NN method, k = 3 was assumed. In SVM method the "one-against-one" approach was employed to solve the problem of classification. The results of the classification are summarised in Table 1. The states: worn teeth and enlarged clearance were separated perfectly. Some of the fault-free sets were recognised as the misalignment sets and vice versa. This can be explained by the fact that during dynamic systems work, the assembly errors change and may temporarily give incorrect classification results. When analysing the time course of the vibration signal (see Figure 3), it can be noticed that some cycles of fault-free and misalignment states have a similar character. However, these are occasional errors occurring once in a dozen or so cases. If there is a permanent change in the system state, the mean value of the instant cumulants will move (see Figure 6).

Conclusions
The paper proposed Multifractal Wavelet Leader method for the automated recognition of the condition of the gearbox based on the monitored vibration signal. Using the wavelet leaders two first cumulants as the signal features were estimated. Then the methods of machine learning were employed for classification of the gearbox state: neural network, nearest neighbours and support vector machines. Regardless of the classification method, satisfactory results were obtained.
The best classification results were obtained for the NN method. The weakest results of the classification of sets: fault-free and misalignment are due to their overlap. This can be explained by the similarity of some fragments of the vibration signal for two states. Despite the complicated mathematical notation, the wavelet leaders method is characterised by low computational complexity and can be performed on-line.


extfieldsofstudy:
Computer Science
Physics

text:
Simulation of Wireless Network Throughput in Providing the Monitoring of Building Structures

The article shows the feasibility of using wireless networks in the field of construction monitoring. It is established that there is a need to increase the data transmission rate over the network in order to collect readings from distributed sensors in order to timely detect defects in building structures. An experiment with predefined parameters for the size of the transmitted packets was conducted using a software implementation of a model of network throughput assessment.

Introduction
More and more tasks in the modern world are being automated every day, and last but not least, the construction of buildings. Many aspects of construction have been already automatedyou can take a look at the semi-automatic marking of the terrain, electronic personnel accounting systems, digital construction models, and much more. Technology development plays a key role in the construction of buildings, but at the same time, new technologies are changing consumer requirements for building operation. Nowadays, it is difficult to imagine a world without mobile communications, a workplace or office without Internet access, including wireless onewhether it is Wi-Fi, LTE, or fifth-generation mobile networks (5G).
The concept of complete safety of building structures involves monitoring the state of building structures. The use of wireless networks is the developing technology in the field of construction monitoring [1].
Wireless networks are a promising area to collect readings from distributed sensors in order to timely detect defects in building structures.
With the development of wireless technologies in the field of construction monitoring and an increase in the number of network users, there is a growing need to increase the speed of data transmission over the network. In this regard, it is necessary to search for the necessary parameters at which the throughput of the data communication network will increase. Thus, the urgent task is to simulate the throughput of a Wi-Fi network in the distributed control mode.

Wireless Networks with Distributed Control
Distributed control in 802.11 wireless networks involves transmission of data packets in two ways: using the Basic Access mechanism and using the RTS/CTS mechanism, or as it is also called an access method with a preliminary connection setup.
We will use the model [7] to simulate the information transmission Model constants defined by IEEE 802.11 are:  SIFS, DIFS, and EIFS inter-frame gaps  the length of a RTS frame and a CTS frame. The model input parameters are:  the number of stations in the network -N;  the payload length of the transmitted data frames is the same and is equal to L bytes;  successful delivery confirmation frame length is L ACK bytes;  r is RTS/CTS threshold. Network throughput S is the output parameter, whose study and improvement is the aim of the simulation.
Suppose that the probability of transmission start in an arbitrarily chosen virtual slot does not depend on the history list or on the current states of other stations and is the same for all stations in the network. The duration of the empty slot = is constant and strictly defined by the protocol.
Then p = (1 − τ) is the probability that an arbitrarily selected virtual slot will be empty, p = Nτ(1 − τ) −1 is the probability that an arbitrarily selected virtual slot will be successful, p = 1 − p − p is the probability that an arbitrarily selected virtual slot will be collision, The durations of a successful slot and a collision slot are equal, consequently: and are the transmission times of the data frame and the delivery confirmation frame respectively, and SIFS, DIFS, and EIFS are the inter-frame gaps defined by the IEEE 802.11 protocol.
The length of one frame determines the average amount of information transmitted per successful slot, therefore: = .
The average network throughput S for the basic access method is calculated by the formula: To calculate the probability of the transmission start to an arbitrarily selected virtual slot, we use the simple consideration that this probability is equal to the ratio of the average number of attempts to transmit a frame to the total number of virtual slots that have passed from the moment of frame transmission beginning until it was successfully completed or all transmission attempts were used. Thus, there is a formula: = ̅ ̅ + ̅ , where ̅ is the average number of slots, in which the station is forced to delay the transmission of the current frame.
It is necessary to take into account that is the probability that at least one of the − 1 remaining stations transmits in a time slot in order to calculate the probability that a transmitted packet gets into a collision.
A numerical solution of this expression, knowing the value of , and, therefore, knowing the value of , is found.
In the case of using the RTS/CTS mechanism to transmit packets, whose length is not equal, the formula (1) to calculate the network throughput will be slightly changed: IOP Conf. Series: Journal of Physics: Conf. Series 1425 (2020) 012118 IOP Publishing doi:10.1088/1742-6596/1425/1/012118 where = + + is still the average duration of a virtual slot. And the probability that the packet has a length is equal to , therefore: The distribution function will be calculated by the formula: Since the RTS/CTS mechanism is used, the calculation of the duration of a successful and collision slot will depend on the RTS/CTS threshold. If the packet length exceeds it, then the RTS/CTS mechanism is used, otherwise the basic method is used.
Thus, the general formula to calculate the duration of a successful slot will look like: is packet header; is the number of bytes; is the channel bit rate.
Collision slot duration will be calculated by the formula: where is the probability that the stations selected this slot for transmission, if it is known that a collision has occurred, ( ) is the average duration of a collision, provided that n stations got into a collision.
After calculating the durations of successful and collision slots and substituting the obtained values in formula (2), we can find the throughput of the entire system as a whole.

Simulation of Wi-Fi Network Throughput in Distributed Control Mode
An R program was written in the RStudio integrated development environment, which allows you to read and display in the form of a matrix of values and in the form of a graph the network throughput with the specified parameters to simulate the Wi-Fi network throughput in distributed control mode.
Network throughput depends on the duration of the transmitted packet, since the size of the packet depends on how the transmission will occur: according to the Basic Method or using the RTS/CTS mechanism. A method to transmit a frame is selected by comparing the packet size with the threshold value of RTS/CTS. However, network throughput depends not only on the selected RTS/CTS threshold, but also on the size of the concurrent window. The channel usage efficiency depends on the duration of the CW concurrent window. The right choice of parameters for the delay mechanism greatly affects the throughput. Improper selection of CW parameters causes a decrease in throughput and the average packet delay can increase by several times.
Let's apply the developed mathematical model to asses the network efficiency with various parameters for 50 stations with various parameters of the transmitted data. The number of stations was selected as the most optimal one due to the frequency of its use. In this case, we will vary the values of the CW delay window, theRTS/CTS threshold values, and the maximum value of the delay stage m. IOP Conf. Series: Journal of Physics: Conf. Series 1425 (2020) 012118 IOP Publishing doi:10.1088/1742-6596/1425/1/012118 4 We will alsoconsider the values obtained for 15 and 100 stations for a detailed consideration of the obtained values.
The initial data transmission parameters were a 15-bit interval delay window size and an RTS/CTS threshold size of 200.
The initial data transmission parameters were a 15-bit interval delay window size and an RTS/CTS threshold size of 200.     The best value of throughput for operation of 100 stations is observed at the highest threshold value (R = 800). However, it is equal to only 2.617537 (Table 2). Let's consider how the value of throughput will change with an increase in the delay window. To do this, let's take different CW values with their maximum delay stages with the same RTS/CTS threshold value of 200 (Table 3).  Figures 5 and 6 show that the maximum value of the throughput begins to shift towards an increase in the number of stations capable to transmit packets with a given throughput.  It was previously noted that throughput of 100 stations network is better with a threshold RTS/CTS of 800. Let's consider how the network performance will change in this case (   Figures 7 and 8 show that although the maximum throughput is observed with the smallest delay window, in this case the network throughput decreases uniformly with an increase in the number of stations. Therefore, the largest value of the delay window will be the best value for the operation of a large number of station. Let's significantly increase the value of the delay window (CW = 255, 511, 1023) and consider the resulting network parameters (Table 5).  Figure 9 shows that, fora small number of stations (N = 15), throughput is the highest when the smallest delay window (CW = 255) is chosen. Figure 10 shows that the best value for the delay window for stations is CW = 1023. However, there are three intersections in the graphs, which show that each value of the delay window will be the best for a certain number of stations. Let's consider what parameters will be best to transmit data packets to 50 stations (Table 7). Despite the fact that the maximum value of throughput is observed whenCW delay window of 255 is chosen, when transmitting packets between 50 stations, the network throughput will be greater with CW of 511.  Figure 12 shows a comparison of the initial parameters of the system with the final selected values, which shows how much the latter improve the network throughput indicator for 50 stations Thus, the best parameters to transmit using 50 stations are the CW delay window of 511, the maximum delay stage m of 2, and the RTS/CTS threshold of 800. Thus, the network throughputwith various values of the concurrent window, its maximum delay stage and RTS/CTS threshold value was found.
In this case, the value of throughput for 15, 50, and 100 stations was considered. The optimal number of stations was chosen, which is equal to 50.
Software implementation of the developed mathematical model was used to find the necessary parameters with which throughputof the 50 stations network is the greatest. The best parameters to transmit using 50 stations are the CW concurrent window of 511, with maximum delay stage m of 2, and the RTS/CTS threshold of 800.

Conclusion
This article studies a model to asses the throughput of a wireless network in providing the high-speed data transmission, which takes into account influence of many factors. The distributed management mechanism in a wireless network running on the IEEE 802.11 protocol can be used in planning a wireless network to monitor the state of building structures in order to ensure the safety of buildings and structures. Data transmission in such a network occurs in two ways: using the Basic Access mechanism and the access method with a preliminary connection setup (RTS/CTS mechanism). An experiment with predefined parameters for the size of the transmitted packets and their distribution was conducted using a software implementation of a mathematical model of network throughputassessment. The network throughput directly depends on the duration of the transmitted packet, since the size of the packet depends on how the packet transmission will occur: according to the Basic Method or using the RTS/CTS mechanism. A method to transmit a frame is selected by comparing the packet length with the RTS/CTS threshold size. Network throughput depends not only on the selected RTS/CTS threshold, but also on the selected value of the concurrent window. The duration of the CW concurrent window affects the efficiency of the data transmission channelutilization. The network throughput depends on the correct choice of parameters of the delay mechanism. Improper selection of CW parameters causes a decrease in network throughput, and the average packet delay can increase by several times.


Using Covariance Matrix Adaptation Evolutionary Strategy to boost the search accuracy in hierarchic memetic computations

Many global optimization problems arising naturally in science and engineering exhibit some form of intrinsic ill-posedness, such as multimodality and insensitivity. Severe ill-posedness precludes the use of standard regularization techniques and necessitates more specialized approaches, usually comprised of two separate stages global phase, that determines the problem’s modality and provides rough approximations of the solutions, and a local phase, which refines these approximations. In this work, we attempt to improve one of the most efficient currently known approaches Hierarchic Memetic Strategy (HMS) by incorporating the Covariance Matrix Adaptation Evolutionary Strategy (CMA-ES) into its local phase. CMA-ES is a stochastic optimization algorithm that in some sense mimics the behavior of population-based evolutionary algorithms without explicitly evolving the population. This way, it avoids, to an extent, the associated cost of multiple evaluations of the objective function. We compare the performance of the HMS on relatively simple multi∗Corresponding author Email addresses: jsawicki@agh.edu.pl (Jakub Sawicki), marcin.los.91@gmail.com (Marcin Loś), smolka@agh.edu.pl (Maciej Smo lka), julen.alvarez.aramberri@gmail.com (Julen Alvarez-Aramberri) The work presented in this paper has been partially supported by Polish National Science Center grants No. DEC-2015/17/B/ST6/01867 and by the AGH statutory research grant. This research is supported by the Basque Government through the BERC2018-2021 program and by the Spanish Ministry of Science, Innovation and Universities: BCAM Severo Ochoa accreditation SEV-2017-0718. Preprint submitted to Journal of Computational Science July 9, 2019 modal benchmark problems and on an engineering problem. To do so, we consider two configurations: the CMA-ES and the standard SEA (Simple Evolutionary Algorithm). The results demonstrate that the HMS with CMA-ES in the local phase requires less objective function evaluations to provide the same accuracy, making this approach more efficient than the standard SEA.

Introduction
Many global optimization problems provided by environmental science, geophysics, medical diagnosis, economy, and technology are not correctly posed. For instance, the ambiguity appearing by lens design is described in [1], while [2] mentions the ill-conditioning of the heat flux estimation in irradiative dryer/furnace. Some examples of uncertain geophysical problems are [3] (calibration of conceptual Rainfall-Runoff Models) and [4] (investigation of oil and gas resources). The problems of cancer diagnosis (see, e.g. [5]) as well as restoring the high-resolution image from the low-level Electrical Impedance Tomography (EIT) images [6] might also be classified as ill-posed.
Generally, the ill-posedness mentioned above is not caused by neglecting the necessary mathematical preciseness but generally, by the insufficient knowledge about the problem's nature (e.g., lack of the total observability, incomplete and inaccurate measurements, computational errors or even poor formal model of the phenomenon under consideration only available).
Typical symptoms of such ill-conditioning are the objective multimodality or/and the regions of objective insensitivity, which might be called "lowlands" in the optimization landscape (see [7] for details). The traditional formulations and solving methods are difficult to apply for such problems. In particular, strategies leading only to find the global, stable minimizer and approximate it by a gradient based method starting from an arbitrary admissible point, are unable to deliver satisfactory information of solutions to such problems.
More reliable knowledge about the problem should contain the information about all local minimizers and the objective behavior in some neighborhood of each of them.
If the ill-conditioning of the problem is mild and has a form of small perturbation (e.g., the problem is conditionally well-posed in the sense of Tikhonov [8]), then the solution can be obtained by applying a local convex minimization method and some regularized objective supplement (see e.g. [9]).
Unfortunately, global optimization problems with a multimodal and an insensitive objective can rarely be solved by a single convex optimization algorithm. Typical strategies are composed of two phases: the global phase in which the objective modality is determined and the local phase, in which the solution's approximation is refined in the basin of attraction of each local minimizer (see, e.g., [10,11]).
The global phase can be performed by a stochastic search, sometimes of a complex structure, involving post processing of the random sample with an artificial intelligence method (see e.g. [12,13]). It is also possible to apply an evolutionary or memetic approach in the global phase (see [14], [15]). The local phase is traditionally performed by gradient-based methods (see, e.g., [16]), if the total strategy leads to find only one minimizer in each distinguished basin of attraction.
In this paper, we use the Hierarchic Memetic Strategy (HMS), which is a concurrent multi-phase global optimizer [17]. This is a complex strategy aimed at solving ill-posed inverse problems formulated as multi-modal global optimization tasks with expensive objective function evaluations. The goal of HMS is to provide possibly all solutions with satisfactory accuracy and in a reasonable time. The strategy was already tested on several real-world engineering inverse problems and proved its superiority over commonly-used global optimizers like the Simple Evolutionary Algorithm (SEA) and the multi-start method, cf. [18,4,17].
Roughly, HMS consists of running in parallel the tree of demes searching with growing accuracy and locality when passing from the root to leaves. The information about the promising search regions propagates down the tree to introduce new child-demes in the best positions found to date. The root and upper demes are responsible for the global phase, while the leaves are involved in the local HMS phase. It is possible to execute both phases concurrently.
While we have applied the Evolutionary Algorithm (EA) with a normal and α-stable mutations [19] in the global HMS phase, we have already checked the following algorithms in the local phase of HMS: • EA and steepest descent methods based on the accurate and approximated objective gradient (see e.g. [4]), • the steepest descent methods associated with the local Tikhonov regularization [20], • EA equipped with the multi-winner selection [21] and the local objective approximation [22], if the insensitivity regions have to be determined.
In this paper, we concentrate on improving the local phase of the hierarchical memetic search by introducing a Covariance Matrix Adaptation Evolutionary Strategy (CMA-ES) as the search engine in the leaves, similarly as in [14]. With this strategy, we aim to: • Increase the efficiency with respect to the EA strategy by reducing the total number of the most expensive fitness calls in the leaves.
• Reduce the number of parameters for tuning (since CMA-ES includes a self-tuning mechanism. • Obtain additional information about the local fitness behavior drawn from the parameters (e.g., the mean and the correlation matrix) of the CMA-ES sampling measure, especially in the case, when the basins of attraction of local minimizers significantly vary in shape and size.
The structure of the paper is as follows. Section 2.1 shortly describes the HMS strategy. In Section 2.3 we show how to apply CMA-ES as a deme evolutionary engine in HMS. Section 3 shows the results of the evaluation of HMS-CMA-ES composition on multi-modal benchmark functions and on an engineering problem. We provide a comparison of the new strategy with both the single-population CMA-ES and the HMS with SEA as the deme evolutionary engine. We finish the paper with some conclusions.

Hierarchic Memetic Strategy
In this section, we describe the algorithmic framework of our global search strategy. In particular, we focus on multi-population evolutionary algorithms.
The HMS is a mixture of a multi-deme evolutionary algorithm and some deterministic techniques, such as gradient-based optimization methods. The HMS algorithm is employed for the accuracy boosting, dynamic accuracy adjustment for the computation time reduction, sample clustering, and additional evolutionary components. Additionally, it is equipped with a multiwinner selection operator for the objective plateau discovery, cf. [17,21,22]. In this paper, we focus on the evolutionary core of HMS that was its straightforward predecessor and has been developed under its own name, i.e., Hierarchic Genetic Strategy (HGS), cf. [23,24,18,4].
As a multi-population evolutionary strategy, the HMS makes use of a particular data structure that provides an organization among the subpopulations (demes). It is a tree with a fixed maximal height and variable internal node degree (see Figure 1). The tree has the parent-child semantics. An outline of the whole strategy is presented in Algorithms 1 and 2. d.population ← pop 10: d.level ← level 11: d.searchAccuracy ← accuracy(level) 12: return d 13: end function on a particular level of the tree. In a typical setup, the HMS begins with a single parentless deme, i.e., the root. A single HMS global step (called a metaepoch) consists in executing a prescribed number of genetic epochs of each active deme. During a metaepoch, the tree structure remains unchanged. However, after that, each active deme located at a non-maximal tree level tries to create a new (next-level) deme by sampling some points around its current best individual. This operation is called sprouting and it can change the tree structure by adding new child nodes. Nevertheless, it is conditional: a node cannot sprout a new child if its current best is too close to other demes located at the target tree level. A node is active as long as its local stopping condition is not satisfied. The whole tree evolves until a global stopping condition starts to hold. The range and accuracy of the performed evolutionary search are different among the HMS demes. The search is more chaotic and inaccurate in the root population, and it becomes more and more focused and accurate while going towards the leaves. The reason is because the task of higher-level populations is to discover promising search areas that are subsequently explored more thoroughly by lower-level populations. The final solutions to the whole strategy are then discovered by the leaf demes. For the full description of the HGS and the HMS we refer the reader to papers [23,24,17,21,25,22].

Covariance Matrix Adaptation Evolutionary Strategy
CMA-ES is a stochastic, population-based method of optimization, well suited to the semi-local search in continuous domains D ⊂ R n , n ≥ 1. The strategy was introduced and initially tested in [26].
The main idea of CMA-ES is to represent the sampling measure as a multivariate normal distribution N (m k , σ 2 k C k ) appearing in each iteration k. It is parametrized by a tuple (m k , σ k , C k ); where m k is the mean, σ k the standard deviation scaling factor (step size), and C k the normalized covariance matrix.
CMA-ES starts its operation by creating an initial population P 0 by λtimes sampling with return using some initial parameters: the step size σ 0 , a distribution's covariance matrix C 0 of ones at the diagonal and zeros outside, and the mean m 0 chosen in the admissible domain. At each iteration, the tuple of parameters is modified based on the current evaluated population P k to better fit the central part of the basin of attraction of the nearest minimizer (see (1)).
We employ two stopping criteria: the algorithm stops when (a) the parameters (m k+1 , σ k+1 , C k+1 ) do not differ significantly from (m k , σ k , C k ), or when (b) the fitness of at least one individual is below the satisfactory level. The high efficiency of CMA-ES in the local search in irregular landscapes in comparison with the conventional EA comes mainly from the more synthetic and concise representation of the knowledge about the problem gathered during all iterations (contained in the parameters (m k , σ k , C k )). Moreover, its consecutive improvements need less objective evaluations.
The extraordinary capabilities shown by CMA-ES in the local gradientfree search are the result of an economical and flexible representation of a sampling measure density and the deterministic character of the adaptation technique used throughout the strategy. The topic is broadly covered in the literature, in particular in [26,27,28,14].

Improving HMS by the CMA-ES component
We use CMA-ES only in leaves of the HMS hierarchy. The higher levels of the computation tree perform a more global search for which the sampling distribution adaptation offered by SEA is better suited and can be controlled by the mutation and crossover parameters. We argue that CMA-ES with its local-like capabilities is going to be more effective than ordinary SEA at converging to local and global minimizers if the starting population is located in their basins of attraction.

Testing new component
The benchmarks shown here are meant to test the CMA-ES algorithm coupled with HMS applied to solve multimodal problems.
Benchmark functions. In tests 1-4 we use two benchmark functions defined on 2-and 4-dimensional domains. Both are composed of inverted Gaussian functions. The function is defined as: where x 0 ∈ D is a point in the domain, r ∈ R + is a constant, and S is a diagonal matrix of size n with diagonal elements equal to r −1 and zeros outside, and n is the number of dimensions of the domain D. The function g r x 0 has a minimum of 0 at point x 0 and raises to 1 exponentially away from x 0 . Function g r x 0 can be multiplied by another function of the same type to obtain multimodal expressions.
The four-dimension benchmark used in Test 4 is defined by: It has minima in (2, 2, 2, 2), (6, 1, 2, 3) and (7,8,2,3). The radii were scaled up by 2 compared to f 1 , in order to maintain a similar susceptible region. We describe the engineering problem that we employ in Test 6 in Section 3.6. Metrics. We assess the algorithms with the following four metrics: covered minima count (CMC), M covered minima evaluations (M -CME), mean minimum distance (MMD) and mean minimum distance fitness (MMDF).
The number of minima covered (CMC) is calculated by taking all the points generated by the demes (referenced from now on as deme points) and finding the closest one for each minimum. If the point is within a prescribed distance around the minimum, then we consider the minimum covered.
The M -CME metric provides the number of evaluations that an algorithm needs to cover M minima (in our case between one and three) using the idea of the first hitting time.
M -CME is equal to (M − 1)-CME plus the number of evaluations required to cover an additional minimum. The value of M -CME is undefined for a particular run of an algorithm. Besides, it is only included in the statistics when the run is unable to cover all M minima. Nevertheless, the statistics can be evaluated if needed.
MMD is computed as follows for each benchmark run: where P leaf is an aggregated set of points from the leaf-demes, m i is the center of the i th minimum and by · we mean Euclidean norm. Unlike the CMC, MMD is parameterless. MMDF is computed similarly to MMD: In this case, we also consider the points from the leaves that are closest to the minima, and we average their fitness.
Configurations. Table Table 1 includes the parameters used for each of the tests. We select these values based on previous works solving similar problems (see, e.g., [17]).

Test 1 -comparison of HMS-SEA, CMA-ES and HMS-CMA-ES
Here we test the performance of the CMA-ES, the HMS-CMA-ES and the HMS-SEA using the CMC metric and the total number of evaluations. We check the performance of the CMA-ES and we compare it to HMS versions in terms of the number of evaluations. It is the only test which takes the single-deme CMA-ES into account.
We run each configuration 50 times and we show the results in Table 2. CMA-ES alone is unable to find more than one minimum, as expectedit was included to compare the number of evaluations needed. HMS-SEA is inefficient at covering more than one minimum and it uses 9 times larger budget than HMS-CMA-ES. HMS-CMA-ES covers more minima than the HMS-SEA, however, it still is not able to do it with great confidence -it fails to localize any minima in some cases. The budget of HMS-CMA-ES is larger than a single-deme CMA-ES, but the overhead is justified in better performance as far as CMC metric goes.

Test 2 -M -CME comparison
We wanted to check how would using CMA-ES in the lowest level of the HMS hierarchy modify the performance of such strategy. We aim at studying the impact of using the CMA-ES strategy in the lowest level of the HMS.
We run each variant 100 times to collect statistics of the 1-, 2-and 3-CME metrics and the number of runs which succeeded in covering that many minima.
The results are presented in Table 3. The performance of both variants is comparable though. HMS-CMA-ES does better when min and max statistics are taken into account, but the median in the 3-CME is larger than for HMS-SEA variant.

Test 3 -MMD and MMDF in function of budget
To analyze the performance of the strategy depending on the budget, we set up two benchmarks and we run them with different budgets. They have global stopping conditions of their own, so they may actually have used    smaller number of evaluations. We used CMC, MMD and MMDF metrics to compare the variants. We consider a range of budgets from 200 to 5000 and for each budget we run 100 repetitions of each algorithm. Figure 3 shows a significant improvement for budgets larger than 1000 evaluations when the MMD and MMDF metrics are employed. The number of minima covered is slightly better for the CMA-ES.
The sub-optimal peak performance is not a problem with CMA-ES specifically. While CMA-ES is expected to improve the efficiency on the leaf-level, the global search capabilities of the strategy are mostly on the side of HMS. Further configuration of HMS may allow for better CMC results, but it is beyond the scope of this test.
The CMC metric for both variants is around 2.4. However, HMS-SEA required 4 × 10 3 evaluations to reach it while HMS-CMA-ES took only 1.5 × 10 3 -over 2.5 times less.

Test 4 -MMD and MMDF in function of budget for 4D benchmark
To study the strategy on a more demanding problem, we run this test on a 4D benchmark (for the function f 2 defined in (4)). Most of the parameters of the algorithms remain unaltered. Since the 4D domain is more difficult,   we set the size of the population to 100 to enhance the exploration. The budget ranges from 10 3 to 10 6 evaluations. We run 100 repetitions of each algorithm variant per budget.
The results, which include the same metrics as in the previous test, are shown in Figure 4. The first panel of this figure shows that HMS-CMA-ES does better than its counterpart, confidently beginning to cover all the minima around 10 5 evaluations. When the number of budgets is sufficiently large (10 4 or more), the gains of HMS-CMA-ES increase.
The CMC metric value for HMS-SEA stabilizes around 1.8 and takes at least 10 5 evaluations. HMS-CMA-ES achieved this accuracy in 10 4 , which is 6 times less than the SEA variant.

Test 5 -M -CME comparison for 4D benchmark
We run each variant 100 times to collect statistics of the 1-, 2-and 3-CME metrics and the number of runs which succeeded in covering that many minima. The setup comparison is the same as in Section 3.2.
We present the results in Table 4. In this case, HMS-CMA-ES variant is able to identify the minima for a significantly smaller budget than HMS-SEA.

Test 6 -Magnetotelluric (MT) oil exploration method
The MT method [29] is a geophysical prospecting method governed by Maxwell's equations and intended to determine the conductivity (or resistivity) distribution of the Earth's subsurface. This method employs natural sources, while the receivers are placed at the surface of the Earth or the seafloor. Being clean and inexpensive, it has been applied for different purposes, such as the mapping of active faults [30], the study of volcanoes [31], or the exploration of offshore hydrocarbons [32]. of this subsurface formation, the air, a natural source placed at the ionosphere, and the receivers located at the surface of the Earth. We truncate the computational domain with the Perfectly Matched Layer (PML) described in [33] and by setting Dirichlet homogeneous boundary conditions.
To solve the equations involved in this problem, in the present work we employ an adaptive Finite Element Method described in [34].
A thorough introduction of the problem has been performed by Smo lka et al. [17,Sec. 3]. The task is to identify the conductivity σ = {σ i } i=1,...,4 of the regions presented in Figure 5. The domain of the problem D = [0, 6] 4 is mapped to the actual conductivity by a function D The fitness values are transformed by a function R + y → log 10 (y/y 0 ) ∈ R + ∪ {0} , where y 0 is the original fitness value returned by the MT solver at the point (1, 2, 10, 3) in the conductivity space. We compare HMS-SEA and HMS-CMA-ES algorithm variants. Their performance is evaluated by measuring the best fitness that their leaf demes are able to identify. We run them 10 times with a fixed budget of 1000 evaluations.
The results are presented in Figure 6. Using the same number of evaluations, HMS-CMA-ES is able to approach the global optimum fitness value closer.

Conclusions
The paper presents a global search strategy (HMS), proven to be efficient when dealing with ill-posed inverse parametric problems, together with the CMA-ES algorithm.
The HMS operates at several levels of search accuracy. It is known that when the search in the leaves is efficient, the strategy can explore larger regions of the domain (in the search for solutions) at the same computational cost level. Therefore, the objective of this work is to improve the efficiency of the optimization at the finest-search level, that is, in the leaves of the HMS tree.
The HMS strategy (equipped with local gradient-based optimizers) was already tested on real-world engineering problems, cf. [17], and was shown capable of finding almost twice as many good-quality solutions than the multi-start method. At the same time, it was better at avoiding poor-quality solution basins of attraction. However, the SEA with the same budget was unable to recover good-quality solutions.
This paper shows with different numerical examples that the performance of the proposed strategy (the core of the HMS strategy with the CMA-ES in the leaves) is better than when using the SEA in the leaves. In particular, the search is more efficient, in some examples reducing the number of evaluations required to achieve similar performance 9 times.
During its run, the proposed strategy acquires information about the objective function, which is expressed in the CMA-ES sampling distribution. Therefore, an additional benefit that we will explore in the future is to use this information to analyze how the fitness function behaves around its optima. This way, we will be able to determine regions of insensitivity around the solutions of the inverse problem, obtaining, accordingly, more information about the ill-posedness.