Simulating Large-Scale ENPS Models by Means of GPU

Simulating Large-Scale ENPS Models by Means of GPU Manuel Garc´ıa–Quismondo1, Ana Brˆandu¸sa Pavel2, and Mario J. P´erez–Jim´enez1 1Research Group on Natural Computing Dpt. of Computer Science and Artificial Intelligence, University of Sevilla Avda. Reina Mercedes s/n. 41012 Sevilla, Spain E-mail:{mgarciaquismondo,marper}@us.es 2Department of Automatic Control and Systems Engineering, Politehnica University of Bucharest Splaiul Independent¸ei, Nr. 313, sector 6, 090042, Bucharest, Romania E-mail: [email protected] Summary. Enzymatic Numerical P Systems (ENPS), an extension of Numerical P Systems, have been successfully applied to model robot controllers. GPGPU is an innovative technological paradigm which applies the parallel architecture of graphic cards to solve parallel, general–purpose problems. In previous work, a GPU simulator for ENPS was introduced. In this paper, a performance analysis on the simulator is performed in order to experimentally measure the speed-up factors resulting from the simulations. Keywords: Enzymatic Numerical P Systems, GPU, simulation 1 Enzymatic Numerical P Systems Membrane computing is an interdisciplinary field which studies computational models inspired by the compartmental structure of biological cells. There exist many types of membrane systems [12], also known as P systems after mathematician Gh. P˘aun, who introduced them. Numerical P systems (NPS) are a type of P systems in which numerical variables evolve inside the compartments by means of programs; a program (or rule) is composed of a production function and a repartition protocol [11]. The variables have a given initial value and the production function is a multivariate polynomial. The value of the production function for the current values of the variables is distributed among variables in certain compartments according to a repartition protocol. A formal definition of NPS can be found in [11], where this type of P system is introduced with possible applications in economics. Enzymatic numerical P systems (ENPS) represent an extension of NPS, proposed and used in the context of modeling robot controllers [13], [14]. ENPS is a more powerful modelling tool than NPS, as it is proven in several articles [13], [3], [17]. ENPS models allow the existence of more than one rule per membrane than NPS while keeping the deterministic behavior. By using a special type of variables referred as enzymes, ENPS models provide a selection mechanism of the active rules during the computational 138 M. Garc´ıa–Quismondo et al. process. Therefore, ENPS are a flexible and efficient modelling framework that can be successfully used for modeling robot behaviors like obstacle avoidance, localization, follower, etc. [13], [4]. An ENPS model of degree m, m ≥1 is formally defined as follows: Π= (H, µ, (V ar1, E1, P r1, V ar1(0)),...,(V arm, Em, P rm, V arm(0))) (1) where: •His an alphabet that contains msymbols (the labels of the membranes); •µis a membrane structure; •V ariis the set of variables from compartment i, and the initial values for these variables are V ari(0); •Eiis a set of enzyme variables from compartment i,Ei⊆V ari •P riis the set of programs (rules) from compartment i. Programs process variables and have two components; a production function and a repartition protocol. In ENPS models, programs can have one of the two following forms: 1. Non-enzymatic form, which is exactly like the one from standard NPS: P rj,i = (Fj,i(x1,i,...,xki,i), cj,1|v1+... +cj,ni|vni) (2) 2. Enzymatic form P rj,i = (Fj,i(x1,i,...,xki,i), ej,i, cj,1|v1+... +cj,ni|vni) (3) where: –Fj,i(x1,i, ..., xki,i) is the production function; –ej,i ∈Ei, is the enzyme–like variable associated to ej,i; –kirepresents the number of variables in membrane i; –cj,1|v1+. . . +cj,ni|vniis the repartition protocol; –nirepresents the number of variables contained in membrane i, plus the number of variables contained in the parent membrane of i, plus the number of variables contained in the children membranes of i. In ENPS models, all active rules are executed in parallel on each computational step. A rule is always active if it is in the non-enzymatic form. Otherwise, if it is in the enzymatic form, a rule is active only if the associated enzyme–like variable has a greater value than the minimum of the absolute values of the variables involved in the production function. The repartition protocol works like in classical NPS [11] Both NPS and ENPS models have been successfully used for modelling robot controllers [4], [13], [14]. The advantages of using ENPS for this kind of applications are pointed out in [3]. For testing the robot controllers, a Java implementation of a numerical P systems simulator was implemented and used. This Java simulator, SimP, simulates ENPS models. Since ENPS are an extension of NPS, SimP can be used to simulate classical NPS as well. SimP was proposed in [15] and it is available as a free executable version (for a free executable version of SimP, please contact [email protected]). SimP allows the computation of rational production functions and not only polynomials. This is useful for implementing more complex models required for robotics applications. An example of model of a rational production function is presented in subsection 4.3. In this paper, the performance of a parallel and distributed simulator which computes ENPS structures is analysed. The simulator was first proposed in [7]. In order to test the Simulating Large-Scale ENPS Models by Means of GPU 139 parallel simulator and analyse its performance, a replication mechanism of membranes was used. For instance, in swarm robotics applications, the robots in a group may have similar behaviors which need to run in parallel. For example, each robot needs to run an obstacle avoidance behavior. Therefore, the ENPS structure for obstacle avoidance [13] must be executed in parallel for all the robots in the swarm (figure 1). Those applications in which each robot needs to run one or even more membrane controllers in parallel and all robots must work in parallel as well take real advantage of the parallel and distibuted ENPS simulator, which will be further discussed. Fig. 1. Replication of ENPS models can be used for swarm robotics applications 2 An introduction to GPU Computing Graphic cards, also known as Graphic Processing Units (GPUs) are devices whose main task is to solve image rendering problems. These problems are usually massive parallel problems, as they can commonly be reduced to render pixels and vertices in a parallel fashion. As a result, the industry has turned these devices into powerful highly–parallel computers with a large number of processors. However, they have been of limited use to scientist for quite a long time. The reason is that GPUs were only suitable for image– related problems, with little application in the scientific world out of image processing itself. Therefore, the amount of effort required to translate general–purpose problems into their graphical interpretations made scientists use other parallel architectures such as computer clusters and FPGAs [1] [9] [16]. However, NVIDIA changed this landscape by providing a toolkit for general–purpose GPU computing named Compute Unified Device Architecture (CUDA) [19]. This API permitted developers to solve scientific, non–graphical problems on GPUs. From then on, the adoption of GPU computing by the scientific community has gone widespread. As a result numerous scientific papers have shown moderate to impressive performance improvements on a GPU over a CPU [2]. 140 M. Garc´ıa–Quismondo et al. 2.1 The CUDA programming model Nowadays, the number of processors in a GPU can reach up to 448 processor cores and 1.536 processing units per core, thus resulting in a total number of 448×1.536 = 688.128 processing units [19]. Thus, GPUs are structured as a grid of multiprocessor cores or streaming multiprocessors. Each one of these multiprocessors (figure 2) is composed of several simpler processing units known as streaming processors. This hierarchy allows programmers to structure their code in a distributed way, so they can couple processes with a high communication bandwidth with each other, thus clustering and allocating them within the same multiprocessor in the GPU. All processing units execute the same code at the same time, hence applying a computational paradigm known as Single Instruction Multiple Program (SIMP) [16]. CUDA provides an abstract model of GPU Fig. 2. Simplified hardware block diagram of an Nvidia Graphic Card architecture. This model is known as the CUDA programming model. The idea is to make use of an abstraction of the specific graphic card on which the code runs, so as not to force the code depend on particular devices. This model is composed of a grid of elements known as blocks, which are abstractions of the streaming multiprocessors mentioned above. Each one of these blocks is composed of computing elements known as threads [5]. Thus, the CUDA programming model is a multidimensional three–levelled architecture , as the dimension of grids and blocks can be 1, 2 or 3, depending on the configuration options set by the developer. Figure 3 describes graphically this programming model. 2.2 Programming on CUDA–C The firsts steps towards obtaining languages for general–purpose GPU computing involved wrapping currently existing languages for graphics processing with mainstream general–purpose languages, such as C and Fortran. However, these approaches were still hindered by the limitations of graphic cards on treating general–purpose data structures and data flows as pixels and vertices. However, the CUDA programming Simulating Large-Scale ENPS Models by Means of GPU 141 Fig. 3. The CUDA programming model model, combined with more advanced graphic hardware with support for general–purpose parallel programming allowed to extend mainstream languages with specific primitives for GPU computing [16]. In this direction, Nvidia issued the first version of the CUDA–C language in 2006 [19]. This language was an extension of the C language with primitives and operations for general–purpose parallel problems [16]. A CUDA–C program is composed of two different parts: the host code and the device code. The host code is the part of the code executed on the CPU. This code contains calls to pieces of the device code, which is executed on the GPU. The device code is composed of CUDA–C functions which are executed on the GPU, known as kernels. The SIMP paradigm defines that each kernel is executed on all threads at the same time. When it comes to execution, threads are bundled in packages known as warps. Threads in the same warp communicate via a fast on–chip memory. However, the communication with threads in different warps is performed by using a slow off–chip memory. That is the reason why threads in the same warp should have a high communication bandwith with each other, as well as a low communication bandwith with threads out of the warp [5]. In practice, the CUDA programming model claims that each block is assigned a sequence of warps. Therefore, each block takes a warp from its sequence and executes it, assigning each warp thread to a different block thread [16]. The order of this sequence, as well as the block to which each warp is associated, is not controlled by the developer. This entails that the problem should not depend on the order in which the warps are executed, as this order cannot be guaranteed [5]. A more thorough description of the CUDA programming model can be found in [10]. 142 M. Garc´ıa–Quismondo et al. 3 A GPU Simulator for Enzymatic Numerical P Systems 3.1 A brief description of the simulator In a previous work [7], a GPU simulator for Enzymatic Numerical P Systems was introduced. This simulator takes as input an XML file defining a description of an ENPS model and simulates it for a number of cycles defined. This number of cycles can be defined in the XML file or as an optional argument in the call sentence. The simulator does not need to check for errors in the input file. The reason is that the XML format accepted by the simulator is also accepted by SimP, which is a Java simulator for ENPS proposed in [15]. This way, errors in the specification of the ENPS model can be checked on SimP. Therefore, if a specification is regarded as error–free by SimP, then it can be given as input to the GPU simulator. This simulator will be published under GNU GPL version 3 license [8], and it is currently available by contacting the authors. In [7], some open problems were proposed. One of this problem had to do with the comprehensive evaluation of the performance of the simulator on large–scale models composed of a considerably high number of programs. In these cases, the overheads related to the workload distribution and the setup operations needed to run the simulations were supposed to be minimal, in comparison to the speed–up factor obtained from the parallel application of the programs in the model. This is not the case of models with a small numbers of programs. In these cases, the performance gain obtained as a result of the parallel application of programs is suffocated by the considerably time–consuming task of setting up the CUDA programming model elements [7]. 3.2 Functioning of the GPU simulator The functioning of the GPU simulator described in [7] consists of two stages. The first stage initializes the model to simulate. The second stage simulates a computational step of the model. This stage is repeated for each computational step simulated. These stages are described below: Initialization stage: First, all operations concerning the setup of the GPU device itself are carried out. These operations include allocations in the GPU memory and transactions between the CPU memory and the GPU memory. Then, this stage normalizes the coefficients found in the repartition protocols in the input model. In practical terms, this normalization substitutes each coefficient cl,s by cl,s Pkl,i j=0 cl,j . This normalization can be performed on each simulation step, but performing it only once at the beginning of the simulation spares computational time. Computation stage: This stage simulates checks the programs in the model and applicates the checked programs. It consists of four sub–steps, which are: 1. For each program, set its activation. That is, check whether this program is active, that is, is going to be applied on the current computational step. 2. Calculate the production function of the active programs. 3. Set to 0 the values consumed by the active programs, that is, the values of variables such that there is any active program which depends on its value. Simulating Large-Scale ENPS Models by Means of GPU 143 4. Multiply the results of the active production functions by the normalized coefficients calculated on the initialization stage. 5. Add these results to the variables contributed by the active programs, according to their repartition protocols. 4 Performance analysis of the simulator In order to carry out an analysis of the performance and runtimes obtained from the GPU simulator, we have developed a sequential simulator for Enzymatic Numerical P Systems in C language. This sequential simulator has been developed to compare simulation times obtained from a sequential, low–level simulator to those obtained from the GPU simulator. Nevertheless, it can also be used for the efficient simulation of Enzymatic Numerical P Systems in those environments in which no Nvidia card is available. Besides, this simulator takes as input a file which describes an ENPS in the same format that the GPU simulator and SimP. Therefore, the same files can be used for all three simulators, hence sparing time on translations between formats. For a fair comparison between execution times, no memory allocation is performed after the setup stage. This feature is compulsory on the GPU simulator, because all computation steps are implemented by means of kernel calls and all memory in the GPU can only be allocated from the host code [19]. The importance of this feature rises from the fact that memory allocation in C is a time–consuming instruction. Therefore, if there were a significant number of memory allocations on each computational step the performance of the C sequential simulator would be severely hindered. This would result on an even higher speed–up factor due to a bad design of the sequential simulator, instead of a good design of its GPU counterpart. 4.1 Performance comparison with SimP Apart from comparing the CUDA–C simulator with a C–based one, we have also compared the GPU simulation times with the ones obtained from SimP [15]. SimP is an ENPS simulator in Java language. Java programs are executed on a virtual machine which does as a middleware between the actual device and the software. This virtual machine is known as Java Virtual Machine (JVM) [18]. JVM ensures that Java programs can be executed on any device in which JVM is installed, thus guaranteeing complete compatiblity among different hardware architectures. However, this virtual machine approach comes at a cost. Firstly, the programmer loses control of the way in which the memory is managed. For instance, memory objects cannot be freed directly. Instead, an execution thread named garbage collector checks which objects are not referenced anymore in the program and frees the allocated memory. Secondly, the translatiosn from JVM instructions to assembly instructions are performed in runtime. Thus, an overhead in the execution time is produced as a result of these translations. All in all, the programmer cannot control directly the execution flow of Java programs. Therefore, in cases where efficiency is required, Java is not, in most cases, a true rival to low–level languages such as C. 144 M. Garc´ıa–Quismondo et al. 4.2 Generation of input models In [7], some problems about an extensive analysis of the simulator are discussed. One of them has to do with the fact that the existing ENPS models have too few programs for parallel simulations to pay off. The reason is that the setup operations computed at the beginning of parallel simulations take a long time, in comparison to the whole sequential computation runtime. In order to overcome this difficulty, we have taken some ENPS models as reference and replicated them several times. This way, we can have an on–demand number of programs per model. Therefore, the more times the seed models are replicated, the more programs will compose the resulting models. When the number of replications given as input is high enough, the GPU simulation does pay off in terms of execution time. The algorithm used for performing the model replication takes the following inputs: N: The number of copies to perform. Π: The seed model to replicate: Π= (H, µ, (V ar1, P r1, E1, V ar1(0)) . . . (V arm, P rm, Em, V arm(0))). For these inputs, the algorithm takes the following steps: 1. Replicate Πa number of times N. 2. Associate an index o, 1≤o≤Nto each of the copies of Π. As a result, a set Φof ENPS models of degree mis obtained. Formally speaking, Φ={Πo= (Ho, µo,(V ar1,o, E1,o, P r1,o, V ar1,o(0)) . . . (V arm,o, Em,o, P rm,o, V arm,o(0)))}, 1≤o≤N, where: • ∀o|1≤o≤N, Ho={{1, o}. . . {m, o}} • ∀i, o|1≤i≤m, 1≤o≤N, V ari,o ={x1,i,o . . . xki,i,o} • ∀i, o|1≤i≤m, 1≤o≤N, Ei,o ⊆V ari,o is the set of all enzyme–like variables associated to programs in P ri,o. • ∀i, o|1≤i≤m, 1≤o≤N, P ri,o =P r1,i,o . . . P rqi,i,o, where: –P rl,i,o = (Fl,i,o(x1,i,o,...,xki,i,o)→cl,1,o|vo,1+. . . +cl,ni,o|vo,ni) if P rl,i is in non–enzymatic form. –P rl,i,o = (Fl,i,o(x1,i,o,...,xki,i,o)(el,i,o →) cl,1,o|vo,1+. . . +cl,ni,o|vo,ni) if P rl,i is in enzymatic form. • ∀i, o|1≤i≤m, 1≤o≤N, V ari,o(0) = {λ1,i,o . . . λki,i,o} 3. ∀i, o, j, l, (1 ≤i≤m),(1 ≤o≤N),(1 ≤j≤ki),(1 ≤l≤qi), assign a different random value in {1,...,10}to λj,i,o, cl,ni,o and each one of the numerical constants in Fl,i,o. 4. Return a new ENPS model Πr= (Hr, µr,(V ar1,1, E1,1, P r1,1, V ar1,1(0)) . . . (V arm,N , Em,N P rm,N , V arm,N (0)), (V arskin, Eskin, P rskin, V ar(0)skin)), where: •Hr=∪N o=1Ho∪ {skin} •µr= [µ1. . . µN]skin •V arskin =∅ •Eskin =∅ •P rskin =∅ •V ar(0)skin =∅ The replication process has been performed by using Java [18]. For the purposes of parsing the seed model and writing the resulting models, an extension of P–Lingua [6] Simulating Large-Scale ENPS Models by Means of GPU 145 has been developed. Specifically, a new input and a new output format has been included into the P–Lingua framework. The input format delegates the parsing process to SimP [15]. The output format generates an XML description of the model. This description is encoded on the common format accepted by all ENPS simulators (Java, C and CUDA– C). This extension proves the versatility of P–Lingua as a useful, assisting tool for a wide variety of Membrane Computing–related tasks; in our case; for analysing the performance of a GPU simulator. 4.3 Case study In our simulator, two seed models have been replicated. After these replications, the resulting expanded models have been simulated. The first seed model is a dummy one with no particular purpose apart from this performance analysis. This model is an ENPS composed of 2 membranes: Π1= (H, µ, (V ar1, E1, P r1,1, V ar1(0)),(V ar2, E2, P r1,2, V ar2(0))), where: •H= 1,2 •µ= [[]2]1 •V ar1={x1,1, x2,1, x3,1} •E1={x3,1} •P r1,1={3·x1,1(x3,1→)2|x1,1+ 1|x2,1} •V ar1(0) = {1,2,3}, •V ar2={x1,2} •E2=∅ •P r1,2={2·x1,2→2|x1,2} •V ar2(0) = {1} Fig. 4. Dummy ENPS used as seed for replication On the other hand, the second seed model performs a function approximation. Mathematical functions like trigonometric functions, exponential functions, etc. are often used in control algorithms in robotics. Therefore, in the following example an ENPS model which computes exis presented. The proposed GPU simulator also allows computation of rational production functions as well as polynomial functions, which is an