Typical errors (as multiple contents, not necessary list of figures, too long jumping addresses, extra bold line before formulas, formulas are not on their place) do not need comment)
Table of Contents
1. Parallel Computations
15.1 Parallel architectures
15.1.1 SIMD architectures
15.1.2 Symmetric multiprocessors
15.1.3 Cache-coherent NUMA architectures:
15.1.4 Non-cache-coherent NUMA architectures:
15.1.5 No remote memory access architectures
15.1.6 Clusters
15.1.7 Grids
15.2 Performance in practice
15.3 Parallel programming 
15.3.1 MPI programming
15.3.2 OpenMP programming
15.3.3 Other programming models
15.4 Computational models
15.4.1 PRAM
15.4.2 BSP, LogP and QSM 
15.4.3 Mesh, hypercube and butterfly
15.5 Performance in theory
15.6 PRAM algorithms
15.6.1 Prefix
15.6.2 Ranking
15.6.3 Merge
15.6.4 Selection
15.6.5 Sorting
15.7 Mesh algorithms
15.7.1 Prefix on chain
15.7.2 Prefix on square
List of Figures
15.1. SIMD architecture.
15.2. Bus-based SMP architecture.
15.3. ccNUMA architecture.
15.4. Ideal, typical, and super-linear speedup curves.
15.5. Locality optimisation by data transformation.
15.6. A simple MPI program.
15.7. Structure of an OpenMP program.
15.8. Matrix-vector multiply in OpenMP using a parallel loop.
15.9. Parallel random access machine.
15.10. Types of parallel random access machines.
15.11. A chain consisting of six processors.
15.12. A square of size .
15.13. A 3-dimensional cube of size .
15.14. A 4-dimensional hypercube .
15.15. A butterfly model.
15.16. A ring consisting of 6 processors.
15.17. Computation of prefixes of 16 elements using Optimal-Prefix.
15.18. Input data of array ranking and the the result of the ranking.
15.19. Work of algorithm Det-Ranking on the data of Example 15.4.
15.20. Sorting of 16 numbers by algorithm Odd-Even-Merge.
15.21. A work-optimal merge algorithm Optimal-Merge.
15.22. Selection of maximal integer number.
15.23. Prefix computation on square.
Chapter 1. Parallel Computations
Table of Contents
15.1 Parallel architectures
15.1.1 SIMD architectures
15.1.2 Symmetric multiprocessors
15.1.3 Cache-coherent NUMA architectures:
15.1.4 Non-cache-coherent NUMA architectures:
15.1.5 No remote memory access architectures
15.1.6 Clusters
15.1.7 Grids
15.2 Performance in practice
15.3 Parallel programming 
15.3.1 MPI programming
15.3.2 OpenMP programming
15.3.3 Other programming models
15.4 Computational models
15.4.1 PRAM
15.4.2 BSP, LogP and QSM 
15.4.3 Mesh, hypercube and butterfly
15.5 Performance in theory
15.6 PRAM algorithms
15.6.1 Prefix
15.6.2 Ranking
15.6.3 Merge
15.6.4 Selection
15.6.5 Sorting
15.7 Mesh algorithms
15.7.1 Prefix on chain
15.7.2 Prefix on square
Parallel computations is concerned with solving a problem faster by using multiple processors in parallel. These processors may belong to a single machine, or to different machines that communicate through a network. In either case, the use of parallelism requires to split the problem into tasks that can be solved simultaneously.
In the following, we will take a brief look at the history of parallel computing, and then discuss reasons why parallel computing is harder than sequential computing. We explain differences from the related subjects of distributed and concurrent computing, and mention typical application areas. Finally, we outline the rest of this chapter.
Although the history of parallel computing can be followed back even longer, the first parallel computer is commonly said to be Illiac IV, an experimental 64-processor machine that became operational in 1972. The parallel computing area boomed in the late 80s and early 90s when several new companies were founded to build parallel machines of various types. Unfortunately, software was difficult to develop and non-portable at that time. Therefore, the machines were only adopted in the most compute-intensive areas of science and engineering, a market too small to commence for the high development costs. Thus many of the companies had to give up.
On the positive side, people soon discovered that cheap parallel computers can be built by interconnecting standard PCs and workstations. As networks became faster, these so-called clusters soon achieved speeds of the same order as the special-purpose machines. At present, the Top 500 list, a regularly updated survey of the most powerful computers worldwide, contains 42% clusters. Parallel computing also profits from the increasing use of multiprocessor machines which, while designed as servers for web etc., can as well be deployed in parallel computing. Finally, software portability problems have been solved by establishing widely used standards for parallel programming. The most important standards, MPI and OpenMP, will be explained in Subsections the section called “15.3.1 MPI programming” and the section called “15.3.2 OpenMP programming” of this book.
In summary, there is now an affordable hardware basis for parallel computing. Nevertheless, the area has not yet entered the mainstream, which is largely due to difficulties in developing parallel software. Whereas writing a sequential program requires to find an algorithm, that is, a sequence of elementary operations that solves the problem, and to formulate the algorithm in a programming language, parallel computing poses additional challenges: 
•	Elementary operations must be grouped into tasks that can be solved concurrently.
•	The tasks must be scheduled onto processors.
•	Depending on the architecture, data must be distributed to memory modules.
•	Processes and threads must be managed, i.e., started, stopped and so on.
•	Communication and synchronisation must be organised.
Of course, it is not sufficient to find any grouping, schedule etc. that work, but it is necessary to find solutions that lead to fast programs. Performance measures and general approaches to performance optimisation will be discussed in Section the section called “15.2 Performance in practice”, where we will also elaborate on the items above. Unlike in sequential computing, different parallel architectures and programming models favour different algorithms.
In consequence, the design of parallel algorithms is more complex than the design of sequential algorithms. To cope with this complexity, algorithm designers often use simplified models. For instance, the Parallel Random Access Machine (see Subsection the section called “15.4.1 PRAM”) provides a model in which opportunities and limitations of parallelisation can be studied, but it ignores communication and synchronisation costs.
We will now contrast parallel computing with the related fields of distributed and concurrent computing. Like parallel computing, distributed computing uses interconnected processors and divides a problem into tasks, but the purpose of division is different. Whereas in parallel computing, tasks are executed at the same time, in distributed computing tasks are executed at different locations, using different resources. These goals overlap, and many applications can be classified as both parallel and distributed, but the focus is different. Parallel computing emphasises homogeneous architectures, and aims at speeding up applications, whereas distributed computing deals with heterogeneity and openness, so that applications profit from the inclusion of different kinds of resources. Parallel applications are typically stand-alone and predictable, whereas distributed applications consist of components that are brought together at runtime.
Concurrent computing is not bound to the existence of multiple processors, but emphasises the fact that several sub-computations are in progress at the same time. The most important issue is guaranteeing correctness for any execution order, which can be parallel or interleaved. Thus, the relation between concurrency and parallelism is comparable to the situation of reading several books at a time. Reading the books concurrently corresponds to having a bookmark in each of them and to keep track of all stories while switching between books. Reading the books in parallel, in contrast, requires to look into all books at the same time (which is probably impossible in practice). Thus, a concurrent computation may or may not be parallel, but a parallel computation is almost always concurrent. An exception is data parallelism, in which the instructions of a single program are applied to different data in parallel. This approach is followed by SIMD architectures, as described below.
For the emphasis on speed, typical application areas of parallel computing are science and engineering, especially numerical solvers and simulations. These applications tend to have high and increasing computational demands, since more computing power allows one to work with more detailed models that yield more accurate results. A second reason for using parallel machines is their higher memory capacity, due to which more data fit into a fast memory level such as cache.
The rest of this chapter is organised as follows: In Section the section called “15.1 Parallel architectures”, we give a brief overview and classification of current parallel architectures. Then, we introduce basic concepts such as task and process, and discuss performance measures and general approaches to the improvement of efficiency in Section the section called “15.2 Performance in practice”. Next, Section the section called “15.3 Parallel programming ” describes parallel programming models, with focus on the popular MPI and OpenMP standards. After having given this general background, the rest of the chapter delves into the subject of parallel algorithms from a more theoretical perspective. Based on example algorithms, techniques for parallel algorithm design are introduced. Unlike in sequential computing, there is no universally accepted model for parallel algorithm design and analysis, but various models are used depending on purpose. Each of the models represents a different compromise between the conflicting goals of accurately reflecting the structure of real architectures on one hand, and keeping algorithm design and analysis simple on the other. Section the section called “15.4 Computational models” gives an overview of the models, Section the section called “15.5 Performance in theory” introduces the basic concepts of parallel algorithmics, Sections the section called “15.6 PRAM algorithms” and the section called “15.7 Mesh algorithms” explain deterministic example algorithms for PRAM and mesh computational model.
15.1 Parallel architectures
A simple, but well-known classification of parallel architectures has been given in 1972 by Michael Flynn. He distinguishes computers into four classes: SISD, SIMD, MISD, and MIMD architectures, as follows:
•	SI stands for “single instruction”, that is, the machine carries out a single instruction at a time.
•	MI stands for “multiple instruction”, that is, different processors may carry out different instructions at a time.
•	SD stands for “single data”, that is, only one data item is processed at a time.
•	MD stands for “multiple data”, that is, multiple data items may be processed at a time.
SISD computers are von-Neumann machines. MISD computers have probably never been built. Early parallel computers were SIMD, but today most parallel computers are MIMD. Although the scheme is of limited classification power, the abbreviations are widely used.
The following more detailed classification distinguishes parallel machines into SIMD, SMP, ccNUMA, nccNUMA, NORMA, clusters, and grids.
15.1.1 SIMD architectures
As depicted in Figure Figure 15.1, “ SIMD architecture.”, a SIMD computer is composed of a powerful control processor and several less powerful processing elements (PEs). The PEs are typically arranged as a mesh so that each PE can communicate with its immediate neighbours. A program is a single thread of instructions. The control processor, like the processor of a sequential machine, repeatedly reads a next instruction and decodes it. If the instruction is sequential, the control processor carries out the instruction on data in its own memory. If the instruction is parallel, the control processor broadcasts the instruction to the various PEs, and these simultaneously apply the instruction to different data in their respective memories. As an example, let the instruction be LD reg, 100. Then, all processors load the contents of memory address 100 to reg, but memory address 100 is physically different for each of them. Thus, all processors carry out the same instruction, but read different values (therefore “SIMD”). For a statement of the form if test then if_branch else else_branch, first all processors carry out the test simultaneously, then some carry out if_branch while the rest sits idle, and finally the rest carries out else_branch while the formers sit idle. In consequence, SIMD computers are only suited for applications with a regular structure. The architectures have been important historically, but nowadays have almost disappeared.
Figure 15.1.  SIMD architecture.

15.1.2 Symmetric multiprocessors
Symmetric multiprocessors (SMP) contain multiple processors that are connected to a single memory. Each processor may access each memory location through standard load/store operations of the hardware. Therefore, programs, including the operating system, must only be stored once. The memory can be physically divided into modules, but the access time is the same for each pair of a processor and a memory module (therefore “symmetric”). The processors are connected to the memory by a bus (see Figure Figure 15.2, “ Bus-based SMP architecture.”), by a crossbar, or by a network of switches. In either case, there is a delay for memory accesses which, partially due to competition for network resources, grows with the number of processors.
Figure 15.2.  Bus-based SMP architecture.

In addition to main memory, each processor has one or several levels of cache with faster access. Between memory and cache, data are moved in units of cache lines. Storing a data item in multiple caches (and writing to it) gives rise to coherency problems. In particular, we speak of false sharing if several processors access the same cache line, but use different portions of it. Since coherency mechanisms work at the granularity of cache lines, each processor assumes that the other would have updated its data, and therefore the cache line is sent back and forth.
15.1.3 Cache-coherent NUMA architectures:  architectures
NUMA stands for Non-Uniform Memory Access, and contrasts with the symmetry property of the previous class. The general structure of ccNUMA architectures is depicted in Figure Figure 15.3, “ ccNUMA architecture.”.
Figure 15.3.  ccNUMA architecture.

As shown in the figure, each processor owns a local memory, which can be accessed faster than the rest called remote memory. All memory is accessed through standard load/store operations, and hence programs, including the operating system, must only be stored once. As in SMPs, each processor owns one or several levels of cache; cache coherency is taken care of by the hardware.
15.1.4 Non-cache-coherent NUMA architectures:  architectures
nccNUMA (non cache coherent Non-Uniform Memory Access) architectures differ from ccNUMA architectures in that the hardware puts into a processor's cache only data from local memory. Access to remote memory can still be accomplished through standard load/store operations, but it is now up to the operating system to first move the corresponding page to local memory. This difference simplifies hardware design, and thus nccNUMA machines scale to higher processor numbers. On the backside, the operating system gets more complicated, and the access time to remote memory grows. The overall structure of Figure Figure 15.3, “ ccNUMA architecture.” applies to nccNUMA architectures as well.
15.1.5 No remote memory access architectures
NORMA (NO Remote Memory Acess) architectures differ from the previous class in that the remote memory must be accessed through slower I/O operations as opposed to load/store operations. Each node, consisting of processor, cache and local memory, as depicted in Figure Figure 15.3, “ ccNUMA architecture.”, holds an own copy of the operating system, or at least of central parts thereof. Whereas SMP, ccNUMA, and nccNUMA architectures are commonly classified as shared memory machines, SIMD architectures, NORMA architectures, clusters, and grids (see below) fall under the heading of distributed memory.
15.1.6 Clusters
According to Pfister, a cluster is a type of parallel or distributed system that consists of a collection of interconnected whole computers that are used as a single, unified computing resource. Here, the term “whole computer” denotes a PC, workstation or, increasingly important, SMP, that is, a node that consists of processor(s), memory, possibly peripheries, and operating system. The use as a single, unified computing resource is also denoted as single system image SSI. For instance, we speak of SSI if it is possible to login into the system instead of into individual nodes, or if there is a single file system. Obviously, the SSI property is gradual, and hence the borderline to distributed systems is fuzzy. The borderline to NORMA architectures is fuzzy as well, where the classification depends on the degree to which the system is designed as a whole instead of built from individual components.
Clusters can be classified according to their use for parallel computing, high throughput computing, or high availability. Parallel computing clusters can be further divided into dedicated clusters, which are solely built for the use as parallel machines, and campus-wide clusters, which are distributed systems with part-time use as a cluster. Dedicated clusters typically do not contain peripheries in their nodes, and are interconnected through a high-speed network. Nodes of campus-wide clusters, in contrast, are often desktop PCs, and the standard network is used for intra-cluster communication.
15.1.7 Grids
A grid is a hardware/software infrastructure for shared usage of resources and problem solution. Grids enable coordinated access to resources such as processors, memories, data, devices, and so on. Parallel computing is one out of several emerging application areas. Grids differ from other parallel architectures in that they are large, heterogeneous, and dynamic. Management is complicated by the fact that grids cross organisational boundaries.
15.2 Performance in practice
As explained in the introduction, parallel computing splits a problem into tasks that are solved independently. The tasks are implemented as either processes or threads. A detailed discussion of these concepts can be found in operating system textbooks such as Tanenbaum. Briefly stated, processes are programs in execution. For each process, information about resources such as memory segments, files, and signals is stored, whereas threads exist within processes such that multiple threads share resources. In particular, threads of a process have access to shared memory, while processes (usually) communicate through explicit message exchange. Each thread owns a separate PC and other register values, as well as a stack for local variables. Processes can be considered as units for resource usage, whereas threads are units for execution on the CPU. As less information needs to be stored, it is faster to create, destroy and switch between threads than it is for processes.
Whether threads or processes are used, depends on the architecture. On shared-memory machines, threads are usually faster, although processes may be used for program portability. On distributed memory machines, only processes are a priori available. Threads can be used if there is a software layer (distributed shared memory) that implements a shared memory abstraction, but these threads have higher communication costs.
Whereas the notion of tasks is problem-related, the notions of processes and threads refer to implementation. When designing an algorithm, one typically identifies a large number of tasks that can potentially be run in parallel, and then maps several of them onto the same process or thread.
Parallel programs can be written in two styles that can also be mixed: With data parallelism, the same operation is applied to different data at a time. The operation may be a machine instruction, as in SIMD architectures, or a complex operation such as a function application. In the latter case, different processors carry out different instructions at a time. With task parallelism, in contrast, the processes/threads carry out different tasks. Since a function may have an if or case statement as the outermost construct, the borderline between data parallelism and task parallelism is fuzzy.
Parallel programs that are implemented with processes can be further classified as using Single Program Multiple Data (SPMD) or Multiple Program Multiple Data (MPMD) coding styles. With SPMD, all processes run the same program, whereas with MPMD they run different programs. MPMD programs are task-parallel, whereas SPMD programs may be either task-parallel or data-parallel. In SPMD mode, task parallelism is expressed through conditional statements.
As the central goal of parallel computing is to run programs faster, performance measures play an important role in the field. An obvious measure is execution time, yet more frequently the derived measure of speedup is used. For a given problem, speedup is defined by 

where denotes the running time of the fastest sequential algorithm, and denotes the running time of the parallel algorithm on processors. Depending on context, speedup may alternatively refer to using processes or threads instead of processors. A related, but less frequently used measure is efficiency, defined by 

Unrelated to this definition, the term efficiency is also used informally as a synonym for good performance.
Figure Figure 15.4, “ Ideal, typical, and super-linear speedup curves.” shows ideal, typical, and super-linear speedup curves. The ideal curve reflects the assumption that an execution that uses twice as many processors requires half of the time. Hence, ideal speedup corresponds to an efficiency of one. Super-linear speedup may arise due to cache effects, that is, the use of multiple processors increases the total cache size, and thus more data accesses can be served from cache instead of from slower main memory.
Figure 15.4.  Ideal, typical, and super-linear speedup curves.

Typical speedup stays below ideal speedup, and grows up to some number of processors. Beyond that, use of more processors slows down the program. The difference between typical and ideal speedups has several reasons: 
•	Amdahl's law states that each program contains a serial portion that is not amenable to parallelisation. Hence, , and thus , that is, the speedup is bounded from above by a constant. Fortunately, another observation, called Gustafson-Barsis law reduces the practical impact of Amdahl's law. It states that in typical applications, the parallel variant does not speed up a fixed problem, but runs larger instances thereof. In this case, may grow slower than , so that is no longer constant.
•	Task management, that is, the starting, stopping, interrupting and scheduling of processes and threads, induces a certain overhead. Moreover, it is usually impossible, to evenly balance the load among the processes/threads.
•	Communication and synchronisation slow down the program. Communication denotes the exchange of data, and synchronisation denotes other types of coordination such as the guarantee of mutual exclusion. Even with high-speed networks, communication and synchronisation costs are orders of magnitude higher than computation costs. Apart from physical transmission costs, this is due to protocol overhead and delays from competition for network resources.
Performance can be improved by minimising the impact of the factors listed above. Amdahl's law is hard to circumvent, except that a different algorithm with smaller may be devised, possibly at the price of larger . Algorithmic techniques will be covered in later sections; for the moment, we concentrate on the other performance factors.
As explained in the previous section, tasks are implemented as processes or threads such that a process/thread typically carries out multiple tasks. For high performance, the granularity of processes/threads should be chosen in relation to the architecture. Too many processes/threads unnecessarily increase the costs of task management, whereas too few processes/threads lead to poor machine usage. It is useful to map several processes/threads onto the same processor, since the processor can switch when it has to wait for I/O or other delays. Large-granularity processes/threads have the additional advantage of a better communication-to-computation ratio, whereas fine-granularity processes/threads are more amenable to load balancing.
Load balancing can be accomplished with static or dynamic schemes. If the running time of the tasks can be estimated in advance, static schemes are preferable. In these schemes, the programmer assigns to each process/thread some number of tasks with about the same total costs. An example of a dynamic scheme is master/slave. In this scheme, first a master process assigns one task to each slave process. Then, repeatedly, whenever a slave finishes a task, it reports to the master and is assigned a next task, until all tasks have been processed. This scheme achieves good load balancing at the price of overhead for task management.
The highest impact on performance usually comes from reducing communication/synchronisation costs. Obvious improvements result from changes in the architecture or system software, in particular from reducing latency, that is, the delay for accessing a remote data item, and bandwidth, that is, the amount of data that can be transferred per unit of time.
The algorithm designer or application programmer can reduce communication/synchronisation costs by minimising the number of interactions. An important approach to achieve this minimisation is locality optimisation. Locality, a property of (sequential or parallel) programs, reflects the degree of temporal and spatial concentration of accesses to the same data. In distributed-memory architectures, for instance, data should be stored at the processor that uses the data. Locality can be improved by code transformations, data transformations, or a combination thereof. As an example, consider the following program fragment to be executed on three processors: 
       for (i=0; i<N; i++) in parallel
          for (j=0; j<N; j++)
             f(A[i][j]);
Here, the keyword “in parallel” means that the iterations are evenly distributed among the processors so that runs iterations , runs iterations , and runs iterations (rounded if necessary). The function is supposed to be free of side effects. The place of the formulas is not good
Figure 15.5.  Locality optimisation by data transformation.

With the data distribution of Figure Figure 15.5, “ Locality optimisation by data transformation.”a), locality is poor, since many accesses refer to remote memory. Locality can be improved by changing the data distribution to that of Figure Figure 15.5, “ Locality optimisation by data transformation.”b) or, alternatively, by changing the program into 
       for (j=0; j<N; j++) in parallel
          for (i=0; i<N; i++)
             f(A[i][j]);
The second alternative, code transformations, has the advantage of being applicable selectively to a portion of code, whereas data transformations influence the whole program so that an improvement in one part may slow down another. Data distributions are always correct, whereas code transformations must respect data dependencies, which are ordering constraints between statements. For instance, in
       a = 3; (1)
       b = a; (2)
a data dependence occurs between statements (1) and (2). Exchanging the statements would lead to an incorrect program.
On shared-memory architectures, a programmer does not the specify data distribution, but locality has a high impact on performance, as well. Programs run faster if data that are used together are stored in the same cache line. On shared-memory architectures, the data layout is chosen by the compiler, e.g. row-wise in C. The programmer has only indirect influence through the manner in which he or she declares data structures.
Another opportunity to reduce communication costs is replication. For instance, it pays off to store frequently used data at multiple processors, or to repeat short computations instead of communicating the result.
Synchronisations are necessary for correctness, but they slow down program execution, first because of their own execution costs, and second because they cause processes to wait for each other. Therefore, excessive use of synchronisation should be avoided. In particular, critical sections (in which processes/threads require exclusive access to some resource) should be kept at a minimum. We speak of sequentialisation if only one process is active at a time while the others are waiting.
Finally, performance can be improved by latency hiding, that is, parallelism between computation and communication. For instance, a process can start a remote read some time before it needs the result (prefetching), or write data to remote memory in parallel to the following computations.
Exercises
15.2-1 For standard matrix multiplication, identify tasks that can be solved in parallel. Try to identify as many tasks as possible. Then, suggest different opportunities for mapping the tasks onto (a smaller number of) threads, and compare these mappings with respect to their efficiency on a shared-memory architecture. 
15.2-2 Consider a parallel program that takes as input a number and computes as output the number of primes in range . Task of the program should determine whether is a prime, by systematically trying out all potential factors, that is, dividing by . The program is to be implemented with a fixed number of processes or threads. Suggest different opportunities for this implementation and discuss their pros and cons. Take into account both static and dynamic load balancing schemes. 
15.2-3 Determine the data dependencies of the following stencil code: 
       for (t=0; t<tmax; t++)
          for (i=0; i<n; i++)
             for (j=0; j<n; j++)
                a[i][j] += a[i-1][j] + a[i][j-1]
Restructure the code so that it can be parallelised. 
15.2-4 Formulate and prove the bounds of the speedup known as Amdahl law and Gustafson-Barsis law. Explain the virtual contradiction between these laws. What can you say on the practical speedup? 
15.3 Parallel programming 
Partly due to the use of different architectures and the novelty of the field, a large number of parallel programming models has been proposed. The most popular models today are message passing as specified in the Message Passing Interface standard (MPI), and structured shared-memory programming as specified in the OpenMP standard. These programming models are discussed in Subsections the section called “15.3.1 MPI programming” and the section called “15.3.2 OpenMP programming”, respectively. Other important models such as threads programming, data parallelism, and automatic parallelisation are outlined in Subsection the section called “15.3.3 Other programming models”.
15.3.1 MPI programming
As the name says, MPI is based on the programming model of message passing. In this model, several processes run in parallel and communicate with each other by sending and receiving messages. The processes do not have access to a shared memory, but accomplish all communication through explicit message exchange. A communication involves exactly two processes: one that executes a send operation, and another that executes a receive operation. Beyond message passing, MPI includes collective operations and other communication mechanisms.
Message passing is asymmetric in that the sender must state the identity of the receiver, whereas the receiver may either state the identity of the sender, or declare its willingness to receive data from any source. As both sender and receiver must actively take part in a communication, the programmer must plan in advance when a particular pair of processes will communicate. Messages can be exchanged for several purposes: 
•	exchange of data with details such as the size and types of data having been planned in advance by the programmer
•	exchange of control information that concerns a subsequent message exchange, and
•	synchronisation that is achieved since an incoming message informs the receiver about the sender's progress. Additionally, the sender may be informed about the receiver's progress, as will be seen later. Note that synchronisation is a special case of communication.
The MPI standard has been introduced in 1994 by the MPI forum, a group of hardware and software vendors, research laboratories, and universities. A significantly extended version, MPI-2, appeared in 1997. MPI-2 has about the same core functionality as MPI-1, but introduces additional classes of functions.
MPI describes a set of library functions with language binding to C, C++, and Fortran. With notable exceptions in MPI-2, most MPI functions deal with interprocess communication, leaving issues of process management such as facilities to start and stop processes, open. Such facilities must be added outside the standard, and are consequently not portable. For this and other reasons, MPI programs typically use a fixed set of processes that are started together at the beginning of a program run. Programs can be coded in SPMD or MPMD styles. It is possible to write parallel programs using only six base functions: 
•	MPI_Init must be called before any other MPI function.
•	MPI_Finalize must be called after the last MPI function.
•	MPI_Comm_size yields the total number of processes in the program.
•	MPI_Comm_rank yields the number of the calling process, with processes being numbered starting from 0.
•	MPI_Send sends a message. The function has the following parameters:
o	address, size, and data type of the message,
o	number of the receiver,
o	message tag, which is a number that characterises the message in a similar way like the subject characterises an email,
o	communicator, which is a group of processes as explained below.
•	MPI_Recv receives a message. The function has the same parameters as MPI_Send, except that only an upper bound is required for the message size, a wildcard may be used for the sender, and an additional parameter called status returns information about the received message, e.g. sender, size, and tag.
Figure Figure 15.6, “ A simple MPI program.” depicts an example MPI program.
Figure 15.6.  A simple MPI program.

Although the above functions are sufficient to write simple programs, many more functions help to improve the efficiency and/or structure MPI programs. In particular, MPI-1 supports the following classes of functions:
•	Alternative functions for pairwise communication: The base MPI_Send function, also called standard mode send, returns if either the message has been delivered to the receiver, or the message has been buffered by the system. This decision is left to MPI. Variants of MPI_Send enforce one of the alternatives: In synchronous mode, the send function only returns when the receiver has started receiving the message, thus synchronising in both directions. In buffered mode, the system is required to store the message if the receiver has not yet issued MPI_Recv.On both the sender and receiver sides, the functions for standard, synchronous, and buffered modes each come in blocking and nonblocking variants. Blocking variants have been described above. Nonblocking variants return immediately after having been called, to let the sender/receiver continue with program execution while the system accomplishes communication in the background. Nonblocking communications must be completed by a call to MPI_Wait or MPI_Test to make sure the communication has finished and the buffer may be reused. Variants of the completion functions allow to wait for multiple outstanding requests.MPI programs can deadlock, for instance if a process first issues a send to process and then a receive from ; and does the same with respect to . As a possible way-out, MPI supports a combined send/receive function.In many programs, a pair of processes repeatedly exchanges data with the same buffers. To reduce communication overhead in these cases, a kind of address labels can be used, called persistent communication. Finally, MPI functions MPI_Probe and MPI_Iprobe allow to first inspect the size and other characteristics of a message before receiving it.
•	Functions for Datatype Handling: In simple forms of message passing, an array of equally-typed data (e.g. float) is exchanged. Beyond that, MPI allows to combine data of different types in a single message, and to send data from non-contiguous buffers such as every second element of an array. For these purposes, MPI defines two alternative classes of functions: user-defined data types describe a pattern of data positions/types, whereas packaging functions help to put several data into a single buffer. MPI supports heterogeneity by automatically converting data if necessary.
•	Collective communication functions: These functions support frequent patterns of communication such as broadcast (one process sends a data item to all other processes). Although any pattern can be implemented by a sequence of sends/receives, collective functions should be preferred since they improve program compactness/understandability, and often have an optimised implementation. Moreover, implementations can exploit specifics of an architecture, and so a program that is ported to another machine may run efficiently on the new machine as well, by using the optimised implementation of that machine.
•	Group and communicator management functions: As mentioned above, the send and receive functions contain a communicator argument that describes a group of processes. Technically, a communicator is a distributed data structure that tells each process how to reach the other processes of its group, and contains additional information called attributes. The same group may be described by different communicators. A message exchange only takes place if the communicator arguments of MPI_Send and MPI_Recv match. Hence, the use of communicators partitions the messages of a program into disjoint sets that do not influence each other. This way, communicators help structuring programs, and contribute to correctness. For libraries that are implemented with MPI, communicators allow to separate library traffic from traffic of the application program. Groups/communicators are necessary to express collective communications. The attributes in the data structure may contain application-specific information such as an error handler. In addition to the (intra)communicators described so far, MPI supports intercommunicators for communication between different process groups.MPI-2 adds four major groups of functions:
•	Dynamic process management functions: With these functions, new MPI processes can be started during a program run. Additionally, independently started MPI programs (each consisting of multiple processes) can get into contact with each other through a client/server mechanism.
•	One-sided communication functions: One-sided communication is a type of shared-memory communication in which a group of processes agrees to use part of their private address spaces as a common resource. Communication is accomplished by writing into and reading from that shared memory. One-sided communication differs from other shared-memory programming models such as OpenMP in that explicit function calls are required for the memory access.
•	Parallel I/O functions: A large set of functions allows multiple processes to collectively read from or write to the same file.
•	Collective communication functions for intercommunicators: These functions generalise the concept of collective communication to intercommunicators. For instance, a process of one group may broadcast a message to all processes of another group.
15.3.2 OpenMP programming
OpenMP derives its name from being an open standard for multiprocessing, that is for architectures with a shared memory. Because of the shared memory, we speak of threads (as opposed to processes) in this section.
Shared-memory communication is fundamentally different from message passing: Whereas message passing immediately involves two processes, shared-memory communication uncouples the processes by inserting a medium in-between. We speak of read/write instead of send/receive, that is, a thread writes into memory, and another thread later reads from it. The threads need not know each other, and a written value may be read by several threads. Reading and writing may be separated by an arbitrary amount of time. Unlike in message passing, synchronisation must be organised explicitly, to let a reader know when the writing has finished, and to avoid concurrent manipulation of the same data by different threads.
OpenMP is one type of shared-memory programming, while others include one-sided communication as outlined in Subsection the section called “15.3.1 MPI programming”, and threads programming as outlined in Subsection the section called “15.3.3 Other programming models”. OpenMP differs from other models in that it enforces a fork-join structure, which is depicted in Figure Figure 15.7, “ Structure of an OpenMP program.”. A program starts execution as a single thread, called master thread, and later creates a team of threads in a so-called parallel region. The master thread is part of the team. Parallel regions may be nested, but the threads of a team must finish together. As shown in the figure, a program may contain several parallel regions in sequence, with possibly different numbers of threads.
Figure 15.7.  Structure of an OpenMP program.

As another characteristic, OpenMP uses compiler directives as opposed to library functions. Compiler directives are hints that a compiler may or may not take into account. In particular, a sequential compiler ignores the directives. OpenMP supports incremental parallelisation, in which one starts from a sequential program, inserts directives at the most performance-critical sections of code, later inserts more directives if necessary, and so on.
OpenMP has been introduced in 1998, version 2.0 appeared in 2002. In addition to compiler directives, OpenMP uses a few library functions and environment variables. The standard is available for C, C++, and Fortran.
Programming OpenMP is easier than programming MPI since the compiler does part of the work. An OpenMP programmer chooses the number of threads, and then specifies work sharing in one of the following ways:
•	Explicitly: A thread can request its own number by calling the library function omp_get_thread_num. Then, a conditional statement evaluating this number explicitly assigns tasks to the threads, similar as in SPMD-style MPI programs.
•	Parallel loop: The compiler directive #pragma omp parallel for indicates that the following for loop may be executed in parallel so that each thread carries out several iterations (tasks). An example is given in Figure Figure 15.8, “ Matrix-vector multiply in OpenMP using a parallel loop.”. The programmer can influence the work sharing by specifying parameters such as schedule(static) or schedule(dynamic). Static scheduling means that each thread gets an about equal-sized block of consecutive iterations. Dynamic scheduling means that first each thread is assigned one iteration, and then, repeatedly, a thread that has finished an iteration gets the next one, as in the master/slave paradigma described before for MPI. Different from master/slave, the compiler decides which thread carries out which tasks, and inserts the necessary communications.
•	Task-parallel sections: The directive #pragma omp parallel sections allows to specify a list of tasks that are assigned to the available threads.
Threads communicate through shared memory, that is, they write to or read from shared variables. Only part of the variables are shared, while others are private to a particular thread. Whether a variable is private or shared is determined by rules that the programmer can overwrite.
Figure 15.8.  Matrix-vector multiply in OpenMP using a parallel loop.

Many OpenMP directives deal with synchronisation that is necessary for mutual exclusion, and to provide a consistent view of shared memory. Some synchronisations are inserted implicitly by the compiler. For instance, at the end of a parallel loop all threads wait for each other before proceeding with a next loop.
15.3.3 Other programming models
While MPI and OpenMP are the most popular models, other approaches have practical importance as well. Here, we outline threads programming, High Performance Fortran, and automatic parallelisation.
Like OpenMP, threads programming or by Java threads uses shared memory. Threads operate on a lower abstraction level than OpenMP in that the programmer is responsible for all details of thread management and work sharing. In particular, threads are created explicitly, one at a time, and each thread is assigned a function to be carried out. Threads programming focuses on task parallelism, whereas OpenMP programming focuses on data parallelism. Thread programs may be unstructured, that is, any thread may create and stop any other. OpenMP programs are often compiled into thread programs.
Data parallelism provides for a different programming style that is explicitly supported by languages such as High Performance Fortran (HPF). While data parallelism can be expressed in MPI, OpenMP etc., data-parallel languages center on the approach. As one of its major constructs, HPF has a parallel loop whose iterations are carried out independently, that is, without communication. The data-parallel style makes programs easier to understand since there is no need to take care of concurrent activities. On the backside, it may be difficult to force applications into this structure. HPF is targeted at single address space distributed memory architectures, and much of the language deals with expressing data distributions. Whereas MPI programmers distribute data by explicitly sending them to the right place, HPF programmers specify the data distribution on a similar level of abstraction as OpenMP programmers specify the scheduling of parallel loops. Details are left to the compiler. An important concept of OpenMP is the owner-computes rule, according to which the owner of the left-hand side variable of an assignment carries out an operation. Thus, data distribution implies the distribution of computations.
Especially for programs from scientific computing, a significant performance potential comes from parallelising loops. This parallelisation can often be accomplished automatically, by parallelising compilers. In particular, these compilers check for data dependencies. that prevent parallelisation. Many programs can be restructured to circumvent the dependence, for instance by exchanging outer and inner loops. Parallelising compilers find these restructuring for important classes of programs.
Exercises
15.3-1 Sketch an MPI program for the prime number problem of Exercise 15.2-3. The program should deploy the master/slave paradigma. Does your program use SPMD style or MPMD style? 
15.3-2 Modify your program from Exercise 15.3-1 so that it uses collective communication. 
15.3-3 Compare MPI and OpenMP with respect to programmability, that is, give arguments why or to which extent it is easier to program in either MPI or OpenMP. 
15.3-4 Sketch an OpenMP program that implements the stencil code example of Exercise 15.2-3. 
15.4 Computational models
15.4.1 PRAM
The most popular computational model is the Parallel Random Access Machine (PRAM) which is a natural generalisation of the Random Access Machine (RAM).
The PRAM model consists of synchronised processors , a shared memory with memory cells and memories of the processors. Figure Figure 15.9, “ Parallel random access machine.”. shows processors and the shared random access memory 
There are variants of this model. They differ in whether multiple processors are allowed to access the same memory cell in a step, and in how the resulting conflicts are resolved. In particular the following variants are distinguished:
Figure 15.9.  Parallel random access machine.

Types on the base of the properties of read/write operations are
•	EREW (Exclusive-Read Exclusive-Write) PRAM,
•	ERCW (Exclusive-Read Concurrent-Write) PRAM,
•	CREW (Concurrent-Read Exclusive-Write) PRAM,
•	CRCW (Concurrent-Read Concurrent-Write) PRAM.
Figure Figure 15.10, “ Types of parallel random access machines.”(a) shows the case when at most one processor has access a memory cell (ER), and Figure Figure 15.10, “ Types of parallel random access machines.”(d) shows, when multiple processors have access the same cell (CW).
Figure 15.10.  Types of parallel random access machines.

Types of concurrent writing are common, priority, arbitrary, combined.
15.4.2 BSP, LogP and QSM 
Here we consider the models BSP, LogP and QSM.
Bulk-synchronous Parallel Model (BSP) describes a computer as a collection of nodes, each consisting of a processor and memory. BSP supposes the existence of a router and a barrier synchronisation facility. The router transfers messages between the nodes, the barrier synchronises all or a subset of nodes. According to BSP computation is partitioned into supersteps. In a superstep each processor independently performs computations on data in its own memory, and initiates communications with other processors. The communication is guaranteed to complete until the beginning of the next superstep.
is defined such that is the time that is takes to route an -relation under continuous traffic conditions. An -relation is a communication pattern in which each processor sends and receives up to messages. 
The cost of a superstep is determined as , where is the maximum number of communications initiated by any processor. The cost of a program is the sum of the costs of the individual supersteps.
BSP contains a cost model that involves three parameters: the number of processors , the cost of a barrier synchronisation and a characteristics of the available bandwidth .
LogP model was motivated by inaccuracies of BSP and the restrictive requirement to follow the superstep structure. 
While LogP improves on BSP with respect to reflectivity, QSM improves on it with respect to simplicity. In contrast to BSP, QSM is a shared-memory model. As in BSP, the computation is structured into supersteps, and each processor has its own local memory. In a superstep, a processor performs computations on values in the local memory, and initiates read/write operations to the shared memory. All shared-memory accesses complete until the beginning of the next superstep. QSM allows for concurrent reads and writes. Let the maximum number of accesses to any cell in a superstep be . Then QSM charges costs , with , and being defined in BSP. 
15.4.3 Mesh, hypercube and butterfly
Mesh also is a popular computational model. A -dimensional mesh is an sized grid having a processor in each grid point. The edges are the communication lines, working in two directions. Processors are labelled by -tuples, as .
Each processor is a RAM, having a local memory. The local memory of the processor is . Each processor can execute in one step such basic operations as adding, subtraction, multiplication, division, comparison, read and write from/into the local memory, etc. Processors work in synchronised way, according to a global clock.
The simplest mesh is the chain, belonging to the value . Figure Figure 15.11, “ A chain consisting of six processors.” shows a chain consisting of 6 processors.
Figure 15.11.  A chain consisting of six processors.

The processors of a chain are . is connected with , is connected with , the remaining processors are connected with and .
If , then we get a rectangle. If now , then we get a square. Figure Figure 15.12, “ A square of size .” shows a square of size .
Figure 15.12.  A square of size .

A square contains several chains consisting of processors. The processors having identical first index, form a row of processors, and the processors having the same second index form a column of processors. Algorithms running on a square often consists of such operations, executed only by processors of some rows or columns.
If , then the corresponding mesh is a brick. In the special case the mesh is called cube. Figure Figure 15.13, “ A 3-dimensional cube of size .” shows a cube of size .
Figure 15.13.  A 3-dimensional cube of size .

The next model of computation is the d-dimensional hypercube . This model can be considered as the generalisation of the square and cube: the square represented on Figure Figure 15.12, “ A square of size .” is a 2-dimensional, and the cube, represented on Figure Figure 15.13, “ A 3-dimensional cube of size .” is a 3-dimensional hypercube. The processors of can be labelled by a binary number consisting of bits. Two processors of are connected iff the Hamming-distance of their labels equals to 1. Therefore each processors of has neighbours, and the of is . Figure Figure 15.14, “ A 4-dimensional hypercube .” represents .
Figure 15.14.  A 4-dimensional hypercube .

The butterfly model consists of processors and edges. The processors can be labelled by a pair , where is the columnindex and is the level of the given processor. Figure Figure 15.15, “ A butterfly model.” shows a butterfly model containing 32 processors in 8 columns and in 4 levels.
Figure 15.15.  A butterfly model.

Finally Figure Figure 15.16, “ A ring consisting of 6 processors.” shows a ring containing 6 processors.
Figure 15.16.  A ring consisting of 6 processors.

15.5 Performance in theory
In the previous section we considered the performance measures used in the practice.
In the theoretical investigations the algorithms are tested using abstract computers called computation models.
The required quantity of resources can be characterised using absolute and relative measures. 
Let , resp. denote the time necessary in worst case to solve the problem of size by the sequential algorithm A, resp. parallel algorithm P (using processors).
In a similar way let , resp. the time necessary for algorithm A, resp. P in best case to solve the problem of size (algorithm P can use processors).
Let , resp. the time needed by any sequential, resp. parallel algorithm to solve problem of size (algorithm P can use processors). These times represent a lower bound of the corresponding running time.
Let suppose the distribution function of the problem of size is given. Then let , resp. the expected value of the time necessary for algorithm A, resp. P to solve problem of size (algorithm P uses processors).
In the analysis it is often supposed that the input data of equal size have equal probability. For such cases we use the notation , resp. and termin average running time.
The value of the performance measures and depend on the used computation model too. For the simplicity of notations we suppose that the algorithms determine the computation model.
Usually the context shows in a unique way the investigated problem. If so, then the parameter is omitted.
Among these performance measures hold the following inequalities: 
Equation 15.1. 

Equation 15.2. 

Equation 15.3. 

In a similar way for the characteristic data of the parallel algorithms the following inequalities are true: 
Equation 15.4. 

Equation 15.5. 

Equation 15.6. 

For the expected running time we have 
Equation 15.7. 

Equation 15.8. 

and 
Equation 15.9. 

Equation 15.10. 

These notations can be used not only for the running time, but also for any other resource, as memory requirement, number of messages, etc.
Now we define some relative performance measures.
Speedup shows, how many times is smaller the running time of a parallel algorithm, than the running time of the parallel algorithm solving the same problem. 
The speedup (or relative number of steps or relative speed) of a given parallel algorithm P, comparing it with a given sequential algorithm A, is defined as 
Equation 15.11. 

If for a sequential algorithm A and a parallel algorithm P holds 
Equation 15.12. 

then the speedup of P comparing with A is linear, if 
Equation 15.13. 

then the speedup of P comparing with A is sublinear, and if 
Equation 15.14. 

then the speedup of P comparing with A is superlinear.
In the case of parallel algorithms it is a very important performance measure the work , defined by the product of the running time and the number of the used processors: 
Equation 15.15. 

This definition is used even then if some processors work only in a small fraction of the running time. Therefore the real work can be much smaller, then given by the formula (Equation 15.15).
The efficiency is a measure of the fraction of time for which the processors are usefully employed; it is defined as the ratio of the work of the sequential algorithm to the work of the parallel algorithm P: 
Equation 15.16. 

One can observe, that the ratio of the speedup and the number of the used parallel processors results the same value. If the parallel work is not less than the sequential one, then efficiency is between zero and one, and the relatively large values are beneficial.
In connection with the analysis of the parallel algorithms the work-efficiency is a central concept. If for a parallel algorithm P and sequential algorithm A holds 
Equation 15.17. 

then algorithm P work-optimal comparing with A.
This definition is equivalent with the equality 
Equation 15.18. 

According to this definition a parallel algorithm is work-optimal only if the order of its total work is not greater, than the order of the total work of the considered sequential algorithm.
A weaker requirement is the following. If there exists a finite positive integer such that 
Equation 15.19. 

then algorithm P is work-efficient comparing with A.
If a sequential algorithm A, resp. a parallel algorithm P uses only , resp. units of a given resource, then A, resp. P is called—for the given resource and the considered model of computation—asymptotically optimal.
If an A sequential or a P parallel algorithm uses only the necessary amount of some resource for all possible size of the input, that is , resp. units, and so we have 
Equation 15.20. 

for A and 
Equation 15.21. 

for P, then we say, that the given algorithm is absolute optimal for the given resource and the given computation model. In this case we say, that is the accurate complexity of the given problem.
Comparing two algorithms and having 
Equation 15.22. 

we say, that the speeds of the growths of algorithms and asymptotically have the same order.
Comparing the running times of two algorithms A and B (e.g. in worst case) sometime the estimation depends on : for some values of algorithm A, while for other values of algorithm B is the better. A possible formal definition is as follows. If the functions and are defined for all positive integer , and for some positive integer hold 
1.	;
2.	,
then the number is called crossover point of the functions and .
For example multiplying two matrices according to the definition and algorithm of Strassen we get one crossover point, whose value is about 20.
Exercises
15.5-1 Suppose that the parallel algorithms P and Q solve the selection problem. Algorithm P uses processors and its running time is . Algorithm Q uses processors and its running time is . Determine the work, speedup and efficiency for both algorithms. Are these algorithms work-optimal or at least work-efficient? 
15.5-2 Analyse the following two assertions.
a) Running time of algorithm P is at least .
b) Since the running time of algorithm P is , and the running time of algorithm B is , therefore algorithm B is more efficient. 
15.5-3 Extend the definition of the crossover point to noninteger values and parallel algorithms. 
15.6 PRAM algorithms
In this section we consider parallel algorithms solving simple problems as prefix calculation, ranking of the elements of an array, merging, selection and sorting.
In the analysis of the algorithms we try to give the accurate order of the running time in the worst case and try to decide whether the presented algorithm is work-optimal or at least work-efficient or not. When parallel algorithms are compared with sequential algorithms, always the best known sequential algorithm is chosen.
To describe these algorithms we use the following pseudocode conventions.
         IN PARALLEL FOR  TO 
          DO 
              
             .
             .
             .
              
For PRAM ordered into a square grid of size the instruction begin with
         IN PARALLEL FOR  TO ,  TO 
             DO
For a -dimensional mesh of size the similar instruction begins with
         IN PARALLEL FOR  TO  TO 
             DO
It is allowed that in this commands represents a group of processors.
15.6.1 Prefix
Let be a binary associative operator defined over a set . We suppose that the operator needs only one set and the set is closed for this operation. 
A binary operation is associative on a set, if for all holds 
Equation 15.23. 

Let the elements of the sequence be elements of the set . Then the input data are the elements of the sequence , and the prefix problem is the computation of the elements . These elements are called prefixes.
It is worth to remark that in other topics of parallel computations the starting sequences of the sequence are called prefixes.
Example 15.1 Associative operations. If is the set of integer numbers, means addition and the sequence of the input data is , then the sequence of the prefixes is . If the alphabet and the input data are the same, but the operation is the multiplication, then . If the operation is the minimum (it is also an associative operation), then . In this case the last prefix is the minimum of the input data.
The prefix problem can be solved by sequential algorithms in time. Any sequential algorithm A requires time to solve the prefix problem. There are parallel algorithms for different models of computation resulting a work-optimal solution of the prefix problem.
In this subsection at first the algorithm CREW-Prefix is introduced, which solves the prefix problem in time, using CREW PRAM processors. 
Next is algorithm EREW-Prefix, having similar quantitative characteristics, but requiring only EREW PRAM processors.
These algorithms solve the prefix problem quicker, then the sequential algorithms, but the order of the necessary work is larger.
Therefore interesting is algorithm Optimal-Prefix, which uses only CREW PRAM processors, and makes only steps. The work of this algorithm is only , therefore its efficiency is , and so it is work-optimal. The speedup of this algorithm equals to .
For the sake of simplicity in the further we write usually instead of . 
A CREW PRAM algorithm.
As first parallel algorithm a recursive algorithm is presented, which runs on CREW PRAM model of computation, uses processors and time. Designing parallel algorithm it is often used the principle divide-and-conquer, as we we will see in the case of the next algorithm too
Input is the number of processors and the array , output data are the array . We suppose is a power of 2. Since we use the algorithms always with the same number of processors, therefore we omit the number of processors from the list of input parameters. In the mathematical descriptions we prefer to consider and as sequences, while in the pseudocodes sometimes as arrays. 
CREW-Prefix()
  1  IF  
  2    THEN  
  3       RETURN  
  4  IF  
  5    THEN  IN PARALLEL FOR  TO  
             DO    compute recursive ,
                the prefixes, belonging to 
            IN PARALLEL FOR  TO 
             DO compute recursive 
                the prefixes, belonging to 
  6      IN PARALLEL FOR  
             DO read  from the global memory and compute 
  7  RETURN  
Example 15.2 Calculation of prefixes of 8 elements on 8 processors. Let and . The input data of the prefix calculation are 12, 3, 6, 8, 11, 4, 5 and 7, the associative operation is the addition. 
The run of the recursive algorithm consists of rounds. In the first round (step 4) the first four processors get the input data 12, 3, 6, 8, and compute recursively the prefixes 12, 15, 21, 29 as output. At the same time the other four processors get the input data 11, 4, 5, 7, and compute the prefixes 11, 15, 20, 27.
According to the recursive structure and work as follows. and get and , resp. and get and as input. Recursivity mean for and , that gets and gets , computing at first and , then updates . After this computes and . 
While and , according to step 4, compute the final values and , and compute the local provisional values of and . 
In the second round (step 5) the first four processors stay, the second four processors compute the final values of and , adding to the provisional values 11, 15, 20 and 27 and receiving 40, 44, 49 and 56.
In the remaining part of the section we use the notation instead of and give the number of used processors in verbal form. If , then we usually prefer to use . 
Theorem 15.1 Algorithm CREW-Prefix uses time on p CREW PRAM processors to compute the prefixes of p elements.
Proof. The lines 4–6 require steps, the line 7 does steps. So we get the following recurrence: 
Equation 15.24. 

Solution of this recursive equation is .
CREW-prefix is not work-optimal, since its work is and we know sequential algorithm requiring only time, but it is work-effective, since all sequential prefix algorithms require time.
An EREW PRAM algorithm.
In the following algorithm we use exclusive write instead of the parallel one, therefore it can be implemented on the EREW PRAM model. Its input is the number of processors and the sequence , and its output is the sequence containing the prefixes. 
EREW-Prefix()
  1   
  2   IN PARALLEL FOR  TO  
  3    DO  
  4   
  5  WHILE  
  6    DO  IN PARALLEL FOR  TO  
  7       DO  
  8             
  9  RETURN  
Theorem 15.2 Algorithm EREW-Prefix computes the prefixes of elements on EREW PRAM processors in time.
Proof. The commands in lines 1–3 and 9 are executed in time. Lines 4–7 are executed so many times as the assignment in line 8, that is times.
A work-optimal algorithm.
Next we consider a recursive work-optimal algorithm, which uses CREW PRAM processors. Input is the length of the input sequence and the sequence , output is the sequence , containing the computed prefixes.
Optimal-Prefix()
  1   IN PARALLEL FOR  TO  
  2    DO compute recursive , 
             the prefixes of the following  input data
              
  3   IN PARALLEL FOR  TO  
  4    DO using CREW-Prefix compute , 
             the prefixes of the following  elements:
              
  5   IN PARALLEL FOR  TO  
  6    DO FOR  TO  
  7          DO  
  8   FOR  TO  
  9    DO  
 10  RETURN  
This algorithm runs in logarithmic time. The following two formulas help to show it: 
Equation 15.25. 

and 
Equation 15.26. 

where summing goes using the corresponding associative operation.
Theorem 15.3 (parallel prefix computation in time) Algorithm Optimal-Prefix computes the prefixes of elements on CREW PRAM processors in time. 
Proof. Line 1 runs in time, line 2 runs time, line 3 runs time.
This theorem imply that the work of Optimal-Prefix is , therefore Optimal-Prefix is a work-optimal algorithm.
Figure 15.17.  Computation of prefixes of 16 elements using Optimal-Prefix.

Let the elements of the sequence be the elements of the alphabet . Then the input data of the prefix computation are the elements of the sequence , and the prefix problem is the computation of the elements . These computable elements are called prefixes.
We remark, that in some books on parallel programming often the elements of the sequence are called prefixes.
Example 15.3 Associative operations. If is the set of integers, denotes the addition and the sequence of the input data is 3, -5, 8, 2, 5, 4, then the prefixes are 3, -2, 6, 8, 13, 17. If the alphabet and the input data are the same, the operation is the multiplication, then the output data (prefixes) are 3, -15, -120, -240, -1200, -4800. If the operation is the minimum (it is also associative), then the prefixes are 3, -5, -5, -5, -5, -5. The last prefix equals to the smallest input data.
Sequential prefix calculation can be solved in time. Any A sequential algorithm needs time. There exist work-effective parallel algorithms solving the prefix problem.
Our first parallel algorithm is CREW-Prefix, which uses CREW PRAM processors and requires time. Then we continue with algorithm EREW-Prefix, having similar qualitative characteristics, but running on EREW PRAM model too.
These algorithms solve the prefix problem quicker, than the sequential algorithms, but the order of their work is larger.
Algorithm Optimal-Prefix requires only CREW PRAM processors and in spite of the reduced numbers of processors requires only time. So its work is , therefore its efficiency is and is work-effective. The speedup of the algorithm is .
15.6.2 Ranking
The input of the list ranking problem is a list represented by an array : each element contains the index of its right neighbour (and maybe further data). The task is to determine the rank of the elements. The rank is defined as the number of the right neighbours of the given element.
Since the further data are not necessary to find the solution, for the simplicity we suppose that the elements of the array contain only the index of the right neighbour. This index is called pointer. The pointer of the rightmost element equals to zero.
Example 15.4 Input of list ranking. Let be the array represented in the first row of Figure Figure 15.18, “ Input data of array ranking and the the result of the ranking.”. Then the right neighbour of the element is , the right neighbour of is . is the last element, therefore its rank is 0. The rank of is 1, since only one element, is to right from it. The rank of is 4, since the elements and are right from it. The second row of Figure Figure 15.18, “ Input data of array ranking and the the result of the ranking.” shows the elements of in decreasing order of their ranks.
Figure 15.18.  Input data of array ranking and the the result of the ranking.

The list ranking problem can be solved in linear time using a sequential algorithm. At first we determine the head of the list which is the unique having the property that does not exist an index with . In our case the head of is . The head of the list has the rank , its right neighbour has a rank and finally the rank of the last element is zero.
In this subsection we present a deterministic list ranking algorithm, which uses EREW PRAM processors and in worst case time. The pseudocode of algorithm Det-Ranking is as follows. 
The input of the algorithm is the number of the elements to be ranked , the array containing the index of the right neighbour of the elements of , output is the array containing the computed ranks.
Det-Ranking()
  1  IN PARALLEL FOR  TO  
  2    DO IF  
  3       THEN  
  4       ELSE  
  5  FOR  TO  
  6    DO  IN PARALLEL FOR  TO  
  7       DO IF  
  8          THEN  
  9               
 10  RETURN  
The basic idea behind the algorithm Det-Ranking is the pointer jumping. According to this algorithm at the beginning each element contains the index of its right neighbour, and accordingly its provisional rank equal to 1 (with exception of the last element of the list, whose rank equals to zero). This initial state is represented in the first row of Figure Figure 15.19, “ Work of algorithm Det-Ranking on the data of Example 15.4.”.
Figure 15.19.  Work of algorithm Det-Ranking on the data of Example 15.4.

Then the algorithm modifies the element so, that each element points to the right neighbour of its right neighbour (if it exist, otherwise to the end of the list). This state is represented in the second row of Figure Figure 15.19, “ Work of algorithm Det-Ranking on the data of Example 15.4.”. 
If we have processors, then it can be done in time. After this each element (with exception of the last one) shows to the element whose distance was originally two. In the next step of the pointer jumping the elements will show to such other element whose distance was originally 4 (if there is no such element, then to the last one), as it is shown in the third row of Figure Figure 15.19, “ Work of algorithm Det-Ranking on the data of Example 15.4.”.
In the next step the pointer part of the elements points to the neighbour of distance 8 (or to the last element, if there is no element of distance 8), according to the last row of Figure Figure 15.19, “ Work of algorithm Det-Ranking on the data of Example 15.4.”.
In each step of the algorithm each element updates the information on the number of elements between itself and the element pointed by the pointer. Let , resp. the rank, resp. neighbour field of the element . The initial value of is 1 for the majority of the elements, but is 0 for the rightmost element ( in the first line of Figure Figure 15.19, “ Work of algorithm Det-Ranking on the data of Example 15.4.”). During the pointer jumping gets the new value (if ) gets the new value , if . E.g. in the second row of Figure Figure 15.19, “ Work of algorithm Det-Ranking on the data of Example 15.4.”) , since its previous rank is 1, and the rank of its right neighbour is also 1. After this will be modified to point to . E.g. in the second row of Figure Figure 15.19, “ Work of algorithm Det-Ranking on the data of Example 15.4.” , since the right neighbour of the right neighbour of is .
Theorem 15.4 Algorithm Det-Ranking computes the ranks of an array consisting of elements on EREW PRAM processors in time.
Since the work of Det-Ranking is , this algorithm is not work-optimal, but it is work-efficient.
The list ranking problem corresponds to a list prefix problem, where each element is 1, but the last element of the list is 0. One can easily modify Det-Ranking to get a prefix algorithm.
15.6.3 Merge
The input of the merging problem is two sorted sequences and and the output is one sorted sequence containing the elements of the input.
If the length of the input sequences is , then the merging problem can be solved in time using a sequential processor. Since we have to investigate all elements and write them into the corresponding element of , the running time of any algorithm is . We get this lower bound even in the case when we count only the number of necessary comparisons.
Merge in logarithmic time.
Let and be the input sequences. For the shake of simplicity let be the power of two and let the elements be different.
To merge two sequences of length it is enough to know the ranks of the keys, since then we can write the keys—using processors—into the corresponding memory locations with one parallel write operation. The running time of the following algorithm is a logarithmic, therefore it is called Logarithmic-Merge.
Theorem 15.5 Algorithm Logarithmic-Merge merges two sequences of length on CREW PRAM processors in time.
Proof. Let the rank of element be in (in ). If , then let . If we assign a single processor to the element , then it can determine, using binary search, the number of elements in , which are smaller than . If is known, then computes the rank in the union of and , as . If belongs to , the method is the same.
Summarising the time requirements we get, that using one CREW PRAM processor per element, that is totally processors the running time is .
This algorithm is not work-optimal, only work-efficient.
Odd-even merging algorithm.
This following recursive algorithm Odd-Even-Merge follows the classical divide-and-conquer principle.
Let and be the two input sequences. We suppose that is a power of 2 and the elements of the arrays are different. The output of the algorithm is the sequence , containing the merged elements. This algorithm requires EREW PRAM processors.
Odd-Even-Merge()
  1  IF  
  2    THEN get  by merging  and  with one comparison 
  3    RETURN  
  4  IF  
  5    THEN  IN PARALLEL FOR  TO  
  6       DO merge recursively  and 
  7            to get  
  8      IN PARALLEL FOR  TO  
  9       DO merge recursively  and 
 10            to get  
 11      IN PARALLEL FOR  TO  
 12       DO  
 13            
 14          IF  
 15             THEN  
 16                  
 17                  
 18  RETURN  
Example 15.5 Merge of twice eight numbers. Let = 1, 5, 8, 11, 13, 16, 21, 26 and = 3, 9, 12, 18, 23, 27, 31, 65. Figure Figure 15.20, “ Sorting of 16 numbers by algorithm Odd-Even-Merge.” shows the sort of 16 numbers. 
At first elements of with odd indices form the sequence and elements with even indices form the sequence , and in the same way we get the sequences and . Then comes the recursive merge of the two odd sequences resulting and the recursive merge of the even sequences resulting .
After this Odd-Even-Merge shuffles and , resulting the sequence : the elements of with odd indices come from and the elements with even indices come from . 
Finally we compare the elements of with even index and the next element (that is with , with etc.) and if necessary (that is they are not in the good order) they are changed.
Figure 15.20.  Sorting of 16 numbers by algorithm Odd-Even-Merge.

Theorem 15.6 (merging in time) Algorithm Odd-Even-Merge merges two sequences of length elements in time using EREW PRAM processors.
Proof. Let denote the running time of the algorithm by . Step 1 requires time, Step 2 time. Therefore we get the recursive equation 
Equation 15.27. 

having the solution .
We prove the correctness of this algorithm using the zero-one principle.
A comparison-based sorting algorithm is oblivious, if the sequence of comparisons is fixed (elements of the comparison do not depend on the results of the earlier comparisons). This definition means, that the sequence of the pairs of elements to be compared is given.
Theorem 15.7 (zero-one principle) If a simple comparison-based sorting algorithm correctly sorts an arbitrary 0-1 sequence of length n, then it sorts also correctly any sequence of length n consisting of arbitrary keys. 
Proof. Let A be a comparison-based oblivious sorting algorithm and let be such a sequence of elements, sorted incorrectly by A. Let suppose A sorts in increasing order the elements of . Then the incorrectly sorted sequence contains an element on the -th position in spite of the fact that contains at least keys smaller than . 
Let be the first (having the smallest index) such element of . Substitute in the input sequence the elements smaller than by 0's and the remaining elements by 1's. This modified sequence is a 0-1 sequence therefore A sorts it correctly. This observation implies that in the sorted 0-1 sequence at least 0's precede the 1, written on the place of . 
Now denote the elements of the input sequence smaller than by red colour, and the remaining elements by blue colour (in the original and the transformed sequence too). We can show by induction, that the coloured sequences are identical at the start and remain identical after each comparison. According to colours we have three types of comparisons: blue-blue, red-red and blue-red. If the compared elements have the same colour, in both cases (after a change or not-change) the colours remain unchanged. If we compare elements of different colours, then in both sequences the red element occupy the position with smaller index. So finally we get a contradiction, proving the assertion of the theorem.
Example 15.6 A non comparison-based sorting algorithm. Let be a bit sequence. We can sort this sequence simply counting the zeros, and if we count zeros, then write zeros, then ones. Of course, the general correctness of this algorithm is not guaranteed. Since this algorithm is not comparison-based, therefore this fact does not contradict to the zero-one principle.
But merge is sorting, and Odd-Even-Merge is an oblivious sorting algorithm.
Theorem 15.8 Algorithm Odd-Even-Merge sorts correctly sequences consisting of arbitrary numbers.
Proof. Let and sorted 0-1 sequences of length . Let the number of zeros at the beginning of . Then the number of zeros in equals to , while the number of zeros in is . Therefore the number of zeros in equals to and the number of zeros in equals to .
The difference of and is at most 2. This difference is exactly then 2, if and are both odd numbers. Otherwise the difference is at most 1. Let suppose, that (the proof in the other cases is similar). In this cases contains two additional zeros. When the algorithm shuffles and , begins with an even number of zeros, end an even number of ones, and between the zeros and ones is a short “dirty” part, 0, 1. After the comparison and change in the last step of the algorithm the whole sequence become sorted.
A work-optimal merge algorithm.
Algorithm Work-Optimal-Merge uses only processors, but solves the merging in logarithmic time. This algorithm divides the original problem into parts so, that each part contains approximately elements. 
Let and be the input sequences. Divide into parts so, that each part contain at most elements. Let the parts be denoted by . Let the largest element in be . 
Assign a processor to each element. These processors determine (by binary search) the correct place (according to the sorting) of in . These places divide to parts (some of these parts can be empty). Let denote these parts by . We call the subset corresponding to in (see Figure Figure 15.21, “ A work-optimal merge algorithm Optimal-Merge.”).
Figure 15.21.  A work-optimal merge algorithm Optimal-Merge.

The algorithm gets the merged sequence merging at first with , with and so on, and then joining these merged sequences. 
Theorem 15.9 Algorithm Optimal-Merging merges two sorted sequences of length in time on CREW PRAM processors.
Proof. We use the previous algorithm.
The length of the parts is , but the length of the parts can be much larger. Therefore we repeat the partition. Let an arbitrary pair. If , then and can be merged using one processor in time. But if , then divide into parts—then each part contains at most keys. Assign a processor to each part. This assigned processor finds the subset corresponding to this subsequence in : time is sufficient to do this. So the merge of and can be reduced to subproblems, where each subproblem is the merge of two sequences of length. 
The number of the used processors is , and this is at most , what is not larger then .
This theorem imply, that Optimal-Merging is work-optimal.
Corollary 15.10 Optimal-Merging is work-optimal.
15.6.4 Selection
In the selection problem elements and a positive integer are given and the -th smallest element is to be selected. Since selection requires the investigation of all elements, and our operations can handle at most two elements, so . 
Since it is known sequential algorithm A requiring only time, so A is asymptotically optimal.
The search problem is similar: in that problem the algorithm has to decide, whether a given element appears in the given sequence, and if yes, then where. Here negative answer is also possible and the features of any element decide, whether it corresponds the requirements or not.
We investigate three special cases and work-efficient algorithms to solve them. 
Selection in constant time using processors.
Let , that is we wish to select the largest key. Algorithm Quadratic-Select solves this task in time using CRCW processors. 
The input ( different keys) is the sequence , and the selected largest element is returned as .
Quadratic-Select()
  1  IF  
  2    THEN  
  3       RETURN  
  4   IN PARALLEL FOR  TO ,  TO  
          DO IF 
  5       THEN FALSE 
  6       ELSE TRUE 
  7   IN PARALLEL FOR  TO  
  8    DO TRUE 
  9   IN PARALLEL FOR  TO ,  TO  
 10    IF FALSE 
 11       THEN FALSE 
 12   IN PARALLEL FOR  TO  
 13    DO IF TRUE 
 14       THEN  
 15  RETURN  
In the first round (lines 4–6) the keys are compared in parallel manner, using all the processors. so, that processor computes the logical value . We suppose that the keys are different. If the elements are not different, then we can use instead of the pair (this solution requires an additional number of length bits. Since there is a unique key for which all comparison result FALSE, this unique key can be found with a logical OR operation is lines 7–11. 
Theorem 15.11 (selection in time) Algorithm Quadratic-Select determines the largest key of different keys in time using CRCW common PRAM processors.
Proof. First and third rounds require unit time, the second round requires time, so the total running time is .
The speedup of this algorithm is . The work of the algorithm is . So the efficiency is . It follows that this algorithm is not work-optimal, even it is not work-effective.
Selection in logarithmic time on processors.
Now we show that the maximal element among keys can be found, using even only common CRCW PRAM processors and time. The used technique is the divide-and-conquer. For the simplicity let be a square number.
The input and the output are the same as at the previous algorithm.
Quick-Selection()
  1  IF  
  2    THEN  
  3       RETURN  
  4  IF  
  5    THEN divide the input into groups  and 
             divide the processors into groups 
  6   IN PARALLEL FOR  TO  
  7    DO recursively determines the maximal element  of the group  
  8  Quadratic-Select() 
  9  RETURN  
The algorithm divides the input into groups so, that each group contains elements , and divides the processors into groups so, that group contains processors . Then the group of processors computes recursively the maximum of group . Finally the previous algorithm Quadratic-Select gets as input the sequence and finds the maximum y of the input sequence . 
Theorem 15.12 (selection in time) Algorithm Quick-Select determines the largest of different elements in time using common CRCW PRAM processors.
Let the running time of the algorithm denoted by . Step 1 requires time, step 2 requires time. Therefore satisfies the recursive equation 
Equation 15.28. 

having the solution .
The total work of algorithm Quick-Select is , so its efficiency is , therefore Quick-Select is not work-optimal, it is only work-effective.
Selection from integer numbers.
If the problem is to find the maximum of keys when the keys consist of one bit, then the problem can be solved using a logical OR operation, and so requires only constant time using processors. Now we try to extend this observation. Let be a given positive integer constant, and we suppose the keys are integer numbers, belonging to the interval . Then the keys can be represented using at most bits. For the simplicity we suppose that all the keys are given as binary numbers of length bits. 
The following algorithm Integer-Selection requires only constant time and CRCW PRAM processors to find the maximum.
The basic idea is to partition the bits of the numbers into parts of length . The -th part contains the bits , the number of the parts is . Figure Figure 15.22, “ Selection of maximal integer number.” shows the partition.
Figure 15.22.  Selection of maximal integer number.

The input of Integer-Selection is the number of processors and the sequence containing different integer numbers, and output is the maximal number .
Integer-Selection()
  1  FOR  TO  
  2    DO compute the maximum  of the remaining numbers on the base of 
             their -th part
  3       delete the numbers whose -th part is smaller than  
  4  one of the remaining numbers 
  5  RETURN  
The algorithm starts with searching the maximum on the base of the first part of the numbers. Then it delete the numbers, whose first part is smaller, than the maximum. Then this is repeated for the second, ..., last part of the numbers. Any of the non deleted numbers is maximal.
Theorem 15.13 (selection from integer numbers) If the numbers are integers drawn from the interval , then algorithm Integer-Selection determines the largest number among numbers for any positive in time using CRCW PRAM processors.
Proof. Let suppose that we start with the selection of numbers, whose most significant bits are maximal. Let this maximum in the first part denoted by . It is sure that the numbers whose first part is smaller than are not maximal, therefore can be deleted. If we execute this basis operation for all parts (that is times), then exactly those numbers will be deleted, what are not maximal, and all maximal element remain. 
If a key contains at most bits, then its value is at most . So algorithm Integer-Select in its first step determines the maximum of integer numbers taken from the interval . The algorithm assigns a processor to each number and uses common memory locations , containing initially . In one step processor writes into . Later the maximum of all numbers can be determined from memory cells using processors by Theorem 15.11 in constant time.
General selection.
Let the sequence contain different numbers and the problem is to select the th smallest element of . Let we have CREW processors.
General-Selection()
  1  divide the  processors into  groups  so, that group  
       contains the processors  and divide
       the  elements into  groups  so, that group 
       contains the elements 
  2   IN PARALLEL FOR  TO  
  3    DO determine  (how many elements of  are smaller, than ) 
  4   IN PARALLEL FOR  TO  
  5    DO using Optimal-Prefix determine  
          (how many elements of  are smaller, than )
  6   IN PARALLEL FOR  TO  
  7    DO IF  
  8       THEN RETURN  
Theorem 15.14 (general selection) The algorithm General-Selection determines the -th smallest of different numbers in time using processors.
Proof. In lines 2–3 works as a sequential processor, therefore these lines require time. Lines 4–5 require time according to Theorem 15.3. Lines 6–8 can be executed in constant time, so the total running time is .
The work of General-Selection is , therefore this algorithm is not work-effective. 
15.6.5 Sorting
Given a sequence the sorting problem is to rearrange the elements of e.g. in increasing order.
It is well-known that any A sequential comparison-based sorting algorithm needs comparisons, and there are comparison-based sorting algorithms with running time.
There are also algorithms, using special operations or sorting numbers with special features, which solve the sorting problem in linear time. If we have to investigate all elements of and permitted operations can handle at most 2 elements, then we get . So it is true, that among the comparison-based and also among the non-comparison-based sorting algorithms are asymptotically optimal sequential algorithms. In this subsection we consider three different sorting algorithm.
Sorting in logarithmic time using processors.
Using the ideas of algorithms Quadratic-Selection and Optimal-Prefix we can sort elements using processors in time. 
Quadratic-Sort()
  1  IF  
  2    THEN  
  3       RETURN  
  4   IN PARALLEL FOR  TO ,  TO  
          DO IF 
  5       THEN  
  6       ELSE  
  7  divide the processors into  groups  so, that group  contains 
          processors 
  8   IN PARALLEL FOR  TO  
  9    DO compute  
 10   IN PARALLEL FOR  TO  
 11    DO  
 12  RETURN  
In lines 4–7 the algorithm compares all pairs of the elements (as Quadratic-Selection), then in lines 7–9 (in a similar way as Optimal-Prefix works) it counts, how many elements of is smaller, than the investigated , and finally in lines 10–12 one processor of each group writes the final result into the corresponding memory cell.
Theorem 15.15 (sorting in time) Algorithm Quadratic-Sort sorts elements using CRCW PRAM processors in time.
Proof. Lines 8–9 require time, and the remaining lines require only constant time.
Since the work of Quadratic-Sort is , this algorithm is not work-effective.
Odd-even algorithm with running time.
The next algorithm uses the Odd-Even-Merge algorithm and the classical divide-and-conquer principle. The input is the sequence , containing the numbers to be sorted, and the output is the sequence , containing the sorted numbers.
Odd-Even-Sort()
  1  IF  
  2    THEN  
  3  IF  
  4    THEN let  and . 
  5      IN PARALLEL FOR  TO  
  6       DO sort recursively  to get  
  7      IN PARALLEL FOR  TO  
  8       DO sort recursively  to get  
  9      IN PARALLEL FOR  TO  
 10       DO merge  and  using Odd-Even-Merge() 
 11  RETURN  
The running time of this EREW PRAM algorithm is .
Theorem 15.16 (sorting in time) Algorithm Odd-Even-Sort sorts elements in time using EREW PRAM processors.
Proof. Let be the running time of the algorithm. Lines 3–4 require time, Lines 5–8 require time, and lines 9–10 require time, line 11 require time. Therefore satisfies the recurrence 
Equation 15.29. 

having the solution .
Example 15.7 Sorting on 16 processors. Sort using 16 processors the following numbers: 62, 19, 8, 5, 1, 13, 11, 16, 23, 31, 9, 3, 18, 12, 27, 34. At first we get the odd and even parts, then the first 8 processors gets the sequence , while the other 8 processors get . The output of the first 8 processors is , while the output of the second 8 processors is . The merged final result is .
The work of the algorithm is , its efficiency is , and its speedup is . The algorithm is not work-optimal, but it is work-effective.
Algorithm of Preparata with running time.
If we have more processors, then the running time can be decreased. The following recursive algorithm due to Preparata uses CREW PRAM processors and time. Input is the sequence , and the output is the sequence containing the sorted elements.
Preparata()
  1  IF  
  2    THEN sort  using  processors and Odd-Even-Sort 
  3    RETURN    Let RETURN begin there, where sort
  4  divide the  elements into  parts  so, that each part 
       contains  elements, and divide the processors into  groups
         so, that each group contains  processors
  5   IN PARALLEL FOR  TO  
  6    DO sort the part  recursively to get a sorted sequence  
  7       divide the processors into  groups  
          containing  processors
  8   IN PARALLEL FOR  TO  TO  
  9    DO merge  and  
 10  divide the processors into  groups  so, that each group 
       contains  processors
 11   IN PARALLEL FOR  TO  
 12    DO determine the ranks of the  element in  using the local ranks 
             received in line 9 and using the algorithm Optimal-Prefix
 13        the elements of  having a rank  
 14  RETURN  
This algorithm uses the divide-and-conquer principle. It divides the input into parts, then merges each pair of parts. This merge results local ranks of the elements. The global rank of the elements can be computed summing up these local ranks.
Theorem 15.17 (sorting in time) Algorithm Preparata sorts elements in time using CREW PRAM processors.
Proof. Let the running time be . Lines 4–6 require time, lines 7–12 together . Therefore satisfies the equation 
Equation 15.30. 

having the solution .
The work of Preparata is the same, as the work of Odd-Even-Sort, but the speedup is better: . The efficiency of both algorithms is .
Exercises
15.6-1 The memory cell of the global memory contains some data. Design an algorithm, which copies this data to the memory cells in time, using EREW PRAM processors. 
15.6-2 Design an algorithm which solves the previous Exercise 15.6-1 using only EREW PRAM processors saving the running time. 
15.6-3 Design an algorithm having running time and determining the maximum of numbers using common CRCW PRAM processors. 
15.6-4 Let be a sequence containing keys. Design an algorithm to determine the rank of any key using CREW PRAM processors and time. 
15.6-5 Design an algorithm having running time, which decides using common CRCW PRAM processors, whether element 5 is contained by a given array , and if is contained, then gives the largest index , for which holds. 
15.6-6 Design algorithm to merge two sorted sequence of length in time, using CREW PRAM processors. 
15.6-7 Determine the running time, speedup, work, and efficiency of all algorithms, discussed in this section. 
15.7 Mesh algorithms
To illustrate another model of computation we present two algorithms solving the prefix problem on meshes.
15.7.1 Prefix on chain
Let suppose that processor of the chain stores element in its local memory, and after the parallel computations the prefix will be stored in the local memory of . At first we introduce a naive algorithm. Its input is the sequence of elements , and its output is the sequence , containing the prefixes.
Chain-Prefix()
  1   sends  to  
  2   IN PARALLEL FOR  TO  
  3  FOR  TO  
  4    DO gets  from , then computes and stores  
             stores , and sends  to 
  5   gets  from , then computes and stores  
Saying the truth, this is not a real parallel algorithm.
Theorem 15.18 Algorithm Chain-Prefix determines the prefixes of p elements using a chain in time.
Proof. The cycle in lines 2–5 requires time, line 1 and line 6 requires time.
Since the prefixes can be determined in time using a sequential processor, and , so CHAIN-Prefix is not work-effective.
15.7.2 Prefix on square
An algorithm, similar to Chain-Prefix, can be developed for a square too.
Let us consider a square of size . We need an indexing of the processors. There are many different indexing schemes, but for the next algorithm Square-Prefix sufficient is the one of the simplest solutions, the row-major indexing scheme, where processor gets the index .
The input and the output are the same, as in the case of Chain-Prefix.
The processors form the processor row and the processors form the processor column . The input stored by the processors of row is denoted by , and the similar output is denoted by .
The algorithm works in 3 rounds. In the first round (lines 1–8) processor rows compute the row-local prefixes (working as processors of Chain-Prefix). In the second round (lines 9–17) the column computes the prefixes using the results of the first round, and the processors of this column send the computed prefix to the neighbour . Finally in the third round the rows determine the final prefixes. 
Square-Prefix()
  1   IN PARALLEL FOR  TO  
  2    DO sends  to  
  3   IN PARALLEL FOR  TO  
  4    FOR  TO  
  5       DO gets  from , then computes and 
  6          stores , and sends  to  
  7   IN PARALLEL FOR  TO  
  8    DO gets  from , then computes and stores  
  9   sends  to  
 10   IN PARALLEL FOR  TO  
 11    FOR  TO  
 12       DO gets  from , then computes and stores 
                stores , and sends  to 
 13   gets  from , then computes and stores  
 14   IN PARALLEL FOR  TO  
 15    DO send  to  
 16   IN PARALLEL FOR  TO  
 17    DO sends  to  
 18   IN PARALLEL FOR  DOWNTO  
 19    FOR  TO  
 20       DO gets  from , then computes and 
 21          stores , and sends  to  
 22   IN PARALLEL FOR  TO  
 23    DO gets  from , then computes and stores  
Theorem 15.19 Algorithm Square-Prefix solves the prefix problem using a square of size , major row indexing in time. 
Proof. In the first round lines 1–2 contain 1 parallel operation, lines 3–6 require operations, and line 8 again 1 operation, that is all together operations. In a similar way in the third round lines 18–23 require time units, and in round 2 lines 9–17 require time units. The sum of the necessary time units is .
Example 15.8 Prefix computation on square of size Figure Figure 15.23, “ Prefix computation on square.”(a) shows 16 input elements. In the first round Square-Prefix computes the row-local prefixes, part (b) of the figure show the results. Then in the second round only the processors of the fourth column work, and determine the column-local prefixes – results are in part (c) of the figure. Finally in the third round algorithm determines the final results shown in part (d) of the figure.
Figure 15.23.  Prefix computation on square.

 	CHAPTER NOTES	 	
Basic sources of this chapter are for architectures and models the book of Leopold [68], and the book of Sima, Fountaine and Kacsuk [95], for parallel programming the book due to Kumar et al. [38] and [68], for parallel algorithms the books of Berman and Paul, [15] Cormen, Leiserson and Rivest [22], the book written by Horowitz, Sahni and Rajasekaran [47] and the book [51], and the recent book due to Casanova, Legrand and Robert [20]. 
The website [101] contains the Top 500 list, a regularly updated survey of the most powerful computers worldwide [101]. It contains 42% clusters.
Described classifications of computers are proposed by Flynn [30], and Leopold [68]. The Figures Figure 15.1, “ SIMD architecture.”, Figure 15.2, “ Bus-based SMP architecture.”, Figure 15.3, “ ccNUMA architecture.”, Figure 15.4, “ Ideal, typical, and super-linear speedup curves.”, Figure 15.5, “ Locality optimisation by data transformation.”, Figure 15.7, “ Structure of an OpenMP program.” are taken from the book of Leopold [68], the program Figure 15.6, “ A simple MPI program.” from the book written by Gropp et al. [40].
The clusters are characterised using the book of Pfister [89], grids are presented on the base of the book and manuscript of Foster and Kellerman [32], [33].
With the problems of shared memory deal the book written by Hwang and Xu [48], the book due to Kleiman, Shah, and Smaalders [57], and the textbook of Tanenbaum and van Steen [98].
Details on concepts as tasks, processes and threads can be found in many textbook, e.g. in [94], [97]. Decomposition of the tasks into smaller parts is analysed by Tanenbaum and van Steen [98].
The laws concerning the speedup were described by Amdahl [3], Gustafson-Barsis [42] and Brent [16]. Kandemir, Ramanujam and Choudray review the different methods of the improvement of locality [55]. Wolfe [108] analyses in details the connection between the transformation of the data and the program code. In connection with code optimisation the book published by Kennedy and Allen [56] is a useful source.
The MPI programming model is presented according to Gropp, Snir, Nitzberg, and Lusk [40], while the base of the description of the OpenMP model is the paper due to Chandra, Dragum, Kohr, Dror, McDonald and Menon [21], further a review found on the internet [81]. 
Lewis and Berg [69] discuss pthreads, while Oaks and Wong [80] the Java threads in details. Description of High Performance Fortran can be found in the book Koelbel et al. [59]. Among others Wolfe [108] studied the parallelising compilers.
The concept of PRAM is due to Fortune and Wyllie and is known since 1978 [31]. BSP was proposed in 1990 by Valiant [102]. LogP has been suggested as an alternative of BSP by Culler et al. in 1993 [24]. QSM was introduced in 1999 by Gibbons, Matias and Ramachandran [37].
The majority of the pseudocode conventions used in Section the section called “15.6 PRAM algorithms” and the description of crossover points and comparison of different methods of matrix multiplication can be found in [22].
The Readers interested in further programming models, as skeletons, parallel functional programming, languages of coordination and parallel mobile agents, can find a detailed description in [68]. Further problems and parallel algorithms are analysed in the books of Leighton [65], [66] and in the chapter Memory Management of this book [9] and in the book of Horowitz, Sahni and Rajasekaran [47]. A model of scheduling of parallel processes is discussed in [36], [50], [107].
Cost-optimal parallel merge is analysed by Wu and Olariu in [109]. New ideas (as the application of multiple comparisons to get a constant time sorting algoritm) of parallel sorting can be found in the paper of Gararch, Golub and Kruskal [34]. 
Bibliography
[1] Author. Title. Kiadó. Volume. (year). pp.  To delete this reference
[2] R. G. Addie, M. Zukerman, T. Neame. Broadband Traffic Modeling: Simple Solutions to Hard Problems. IEEE Communications Magazine. 36(8). 1998. 88–95. 
[3] G. M. Amdahl. Validity of the single-processor approach to achieving large-scale computer capabilities, In AFIPS Conference Proceedings. 30. 1967. 483–485. 
[4] D. Anick, D. Mitra, M. Sondhi. Stochastic theory of a data handling system with multiple sources. The Bell System Technical Journal. 61. 1982. 1871–1894. 
[5] J. Aspnes , C. Busch , S. Dolev , F. Panagiota , C. Georgiou , A. Shvartsman , P. Spirakis , R. Wattenhofer . Eight open problems in distributed computing. Bulletin of European Association of Theoretical Computer Science of EATCS . 90. 2006. 109–126. 
[6] H. Attiya, C. Dwork , N. A. , L. J. Stockmeyer. Bounds on the time to reach agreement in the presence of timing uncertainty. Journal of the ACM . 41. 1994. 122–142. 
[7] H. Attiya, J. Welch . Distributed Computing, Fundamentals, Simulations and Advanced Topics. McGraw -Hill. 1998. 
[8] B. Awerbuch. Complexity of network synchronization. Journal of the ACM . 32(4). 1985. 804–823. 
[9] Á. Balogh, A. Iványi . Memory Management, In A. Iványi (Ed.) Algorithms of Informatics. Mondat Kiadó. 2. 2007. 799–850. 
[10] J. Banks , J. Carson, B. Nelson. Discrete-Event Simulation. Discrete-Event Simulation. 1996. 
[11] Trusoft Intl Inc. . Benoit 1.1. Trusoft Intl Inc.. 2007. 
[12] J. Beran. Statistics for Long-Memory Processes, Monographs on Statistics and Applied Probability. Chapman & Hall. 1986. 
[13] J. Beran, R. Sherman, M. Taqqu , W. Willinger. Long-range dependence in variable-bit-rate video traffic. IEEE Transactions on Communications . 43. 1995. 1566–1579. 
[14] P. Berman, J. Garay. Cloture votes: $n/4$-resilient distributed consensus in $t+1$ rounds. Mathematical Systems Theory. 26(1). 1993. 3–19. 
[15] K. A. Berman , J. L. Paul. Sequential and Parallel Algorithms. PWS Publishing Company. 1996. 
[16] R. P. Brent. The Parallel Evaluation of General Arithmetic Expressions. Journal of the ACM . 21. 1974. 201–206. 
[17] J. E. Burns. A formal model for message passing systems. Indiana University. (91. 1980. 
[18] J. E. Burns, N. A. Lynch . Bounds on shared memory for mutual exclusion. Information and Computation. 107(2). 1993. 171–184. 
[19] CACI . COMNET III. CACI Products Co.. 1997. 
[20] H. Casanova , A. Legrand , Y. Robert . Parallel Algorithms. Chapman & Hall. 2009. 
[21] R. Chandra, L. Dagum, D. Kohr, D. Maydan, J. McDonald, R. Menon. Parallel Programming in OpenMP. Morgan Kaufmann Publishers. 2000. 
[22] T. H. Cormen , C. E. Leiserson , R. L. Rivest . Introduction to Algorithms. The MIT Press/McGraw -Hill. 1990. 
[23] M. Crovella, A. Bestavros. Self-Similarity in World Wide Web Traffic: Evidence and Possible Causes. IEEE/ACM Transactions on Networking . 5(6). 1997. 835–846. 
[24] D. E. Culler , R. M. Karp , D. Patterson, A. Sahay, E. E. Santos, K. E. Schauser, R. Subramonian, T. von Eicken. LogP: A practical model of parallel computation. Communication of the ACM . 39(11). 1996. 78–85. 
[25] D. Dolev, R. Strong. Authenticated algorithms for Byzantine agreement. SIAM Journal on Computing. 12(4). 1983. 656–666. 
[26] N. Duffield, N. Oconnell. Large deviation and overflow probabilities for the general single-server queue, with applications. Mathematical Proceedings of the Cambridge Philosophical Society. 118. 1995. 363–374. 
[27] D. Duffy, A. McIntosh, M. Rosenstein, W. Willinger. Statistical analysis of CCSN/SS7 traffic data from working CCS subnetworks. IEEE Journal on Selected Areas Communications. 12. 1994. 544–551. 
[28] A. Erramilli, O. Narayan, W. Willinger. Experimental queueing analysis with long-range dependent packet-traffic. IEEE/ACM Transactions on Networking . 4(2). 1996. 209–223. 
[29] M. J. Fischer, N. A. Lynch , M. S. Paterson. Impossibility of distributed consensus with one faulty proces. Journal of the ACM . 32(2). 1985. 374–382. 
[30] M. J. Flynn . Very high-speed computer systems. Proceedings of the IEEE . 5(6). 1966. 1901–1909. 
[31] S. Fortune, J. Wyllie. Parallelism in Random Access Machines, In Proceedings of the Tenth Annual ACM Symposium on Theory of Computing. 1978. 114–118. 
[32] I. Foster , C. Kesselman . The Grid: Blueprint for a New Computing Infrastructure. Morgan Kaufman Publisher. 2004 (2nd edition). 
[33] I. Foster , C. Kesselman , J. M. Nick, S. Tuecke. The Physiology of the Grid: An Open Grid Services Architecture for Distributed Systems Integration. Available at www.globus.org/research/papers/ogsa.pdf . 2002. 
[34] W. Gararch, E. Evan , C. Kruskal . Proofs that yield nothing but their validity or all languages in NP. Journal of the ACM. 38(3). 1991. 691–729. 
[35] H. Garcia-Molina , J. Seiferas. Elections in a distributed computing systems. IEEE Transactions on Computers. C-31(1). 1982. 47–59. 
[36] P. Gács . Compatible sequences and a slow Winkler percolation. Combinatorics Probability and Computing. 13(6). 2004. 815–856. 
[37] P. B. Gibbons, Y. Matias , V. Ramachandran . Can a shared-memory model serve as a bridging model for parallel computation. Theory of Computing Systems. 32(3). 1999. 327–359. 
[38] A. Grama , A. Gupta , G. Karypis , V. Kumar . Introduction to Parallel Computing. Addison -Wesley. 2003 (2nd edition). 
[39] J. N. Gray. Notes on database operating system, In R. Bayer, R. M. Graham, G. Seegmuller (Eds.) Operating Systems: An Advanced Course, Lecture Notes in Computer Science . Springer -Verlag. 60. 1978. 393–481. 
[40] W. Gropp, M. Snir, B. Nitzberg, E. Lusk. MPI: The Complete Reference, Scientific and Engineering Computation Series. The MIT Press. 1998. 
[41] R. Gusella. A measurement study of diskless workstation traffic on an Ethernet. IEEE Transactions on Communications . 38. 1990. 1557–1568. 
[42] J. L. Gustafson. Reevaluating Amdahl's law. Communications of ACM . 28(1). 1988. 532–535. 
[43] T. Gyires . Simulation of the Harmful Consequences of Self-Similar Network Traffic. The Journal of Computer Information Systems. 42(4). 2002. 94–111. 
[44] T. Gyires . Extension of multiprotocol label switching for long-range dependent traffic: QoS routing and performance in IP networks. Computer Standards and Interfaces. 27. 2005. 117–132. 
[45] H. Hefles, D. Lucantoni. A Markov modulated characterization of packetized voice and data traffic and related statistical multiplexer performance. IEEE Journal on Selected Areas in Communication. 4. 1986. 856–868. 
[46] F. Heylighen . Principia Cybernetica Project. http://pespmc1.vub.ac.be/HEYL.html . 2004. 
[47] E. Horowitz , S. Sahni , S. Rajasekaran . Computer Algorithms. Computer Science Press. 1998. 
[48] K. Hwang, Z. Xu. Scalable Parallel Computing. McGraw-Hill . 1998. 
[49] Interoute . http://wwww.interoute.com/glossary.html . 2004. 
[50] A. Iványi . Density of safe matrices. Acta Universitatis Sapientiae . 1(2). 2009. 121–142. 
[51] A. Iványi . Párhuzamos algoritmusok (Parallel Algorithms). ELTE Eötvös Kiadó. 2003. 
[52] R. Jain, S. Routhier. Packet trains: Measurements and a new model for computer network traffic. IEEE Journal on Selected Areas in Communication. 4. 1986. 986–995. 
[53] K. Jones. Consultant's Guide to COMNET III. CACI Product Company. 1997. 
[54] G. A. Jones, J. Jones. Information and Coding Theory. Springer -Verlag. 2000. 
[55] M. Kandemir , J. Ramanujam , A. Choudhary . Compiler algorithms for optimizing locality and parallelism on shared and distributed-memory machines. Journal of Parallel and Distributed Computing. 60. 2000. 924–965. 
[56] K. Kennedy , R. Allen. Optimizing Compilers for Modern Architectures. Morgan Kaufman Publishers. 2001. 
[57] S. Kleiman, D. Shah, B. Smaalders. Programming with Threads. Prentice Hall. 1996. 
[58] L. Kleinrock . Queueing Systems. John Wiley & Sons. 1975. 
[59] C. H. Koelbel, D. B. Loveman, R. S. Schreiber, G. Steele Jr., M. E. Zosel. The High Performance Fortran Handbook. The MIT Press. 1994. 
[60] A. D. Kshemkalyani , M. Singhal . Distributed Computing. Cambridge University Press. 2008. 
[61] L. Lamport . A new solution of Dijkstra's concurrent programming problem. Communications of the ACM. 18(8). 1974. 453–455. 
[62] L. Lamport . A fast mutual exclusion algorithm. ACM Transactions on Computers. 5(1). 1987. 1–11. 
[63] L. Lamport , R. Shostak, M. Pease. The Byzantine generals problem. ACM Transactions on Programming Languages and Systems. 4(3). 1982. 382–401. 
[64] A. Law , W. Kelton . Simulation Modeling and Analysis. 3rd edition, McGraw-Hill Higher Education. 1976. 
[65] F. T. Leighton . Introduction to Parallel Algorithms and Architectures: Arrays-Trees-Hypercubes. Morgan Kaufman Publishers. 1992. 
[66] F. T. Leighton . Introduction to Parallel Algorithms and Architectures: Algorithms and VSLI. Morgan Kaufman Publishers. 2001. 
[67] W. Leland, M. Taqqu , D. Wilson. On the self-similar nature of Ethernet traffic. Computer Communication Reviews. 23. 1993. 183–193. 
[68] C. Leopold . Parallel and Distributed Computing, Wiley Series on Parallel and Distributed Computing. John Wiley & Sons. 2001. 
[69] B. Lewis, D. J. Berg. Multithreaded Programming with Phtreads. Prentice Hall. 1998. 
[70] M. Listanti, V. Eramo, R. Sabella. Architectural and Technological Issues for Future Optical Internet Networks. IEEE Communications Magazine. 8(9). 2000. 82–92. 
[71] N. A. Lynch . Distributed Algorithms. Morgan Kaufman Publisher. 2001 (5th edition). 
[72] N. A. Lynch , M. J. Fischer. On describing the behavior and implementation of distributed systems. Theoretical Computer Science. 13(1). 1981. 17–43. 
[73] B. Mandelbrot. The Fractal Geometry of Nature. W. H. Freeman . 1982. 
[74] B. Mandelbrot, J. W. Van Ness . Fractional Brownian Motions, Fractional Noises and Applications. SIAM Review. 10. 1968. 422–437. 
[75] McAfee . Sniffer Technologies. http://www.nai.com/us/index.asp . 2004. 
[76] B. Melamed . An overview of TES processes and modeling methodology, Lecture Notes in Computer Science, In L. Donatiello and A. R. Nelson (Eds.) Models and Techniques for Performance Evaluation of Computer and Communications Systems. Springer -Verlag. 1993. 359–393. 
[77] S. Molnár , A. Vidács , A. Nilsson. Bottlenecks on the way towards characterization of network traffic: Estimation and interpretation of the Hurst parameter. http://hsnlab.ttt.bme.hu/~molnar (Conference papers). 1997. 
[78] M. F. Neuts. A versatile Markovian point process. Journal of Applied Probability . 18. 1979. 764–779. 
[79] M. F. Neuts . Structured Stochastic Matrices of M/G/1 Type and Their Applications. Marcel Dekker . 1989. 
[80] S. Oaks, H. Wong. Java Threads. O'Reilly. 1999. 
[81] OpenMP Website. http://www.openmp.org . 2007. 
[82] OPNET . OPNET Modeler Documentation. www.opnet.com . 2007. 
[83] S. Owicki, D. Gries . An axiomatic proof technique for parallel programs I.. Acta Informatica . 6(4). 1976. 319–340. 
[84] S. Owicki, L. Lamport . Proving liveness properties of concurrent programs. ACM Transactions on Programming Languages and Systems. 4(3). 1982. 455–495. 
[85] C. Partridge. The end of simple traffic models. IEEE Network . 7(9). 1993. 
[86] V. Paxson, S. Floyd. Wide-area traffic: The failure of Poisson modeling. IEEE/ACM Transactions on Networking . 3. 1995. 226–244. 
[87] M. Pease, R. Shostak, L. Lamport . Reaching agreement in the presence of faults. Journal of the ACM . 27(2). 1980. 228–234. 
[88] G. Peterson, J. Fischer. Economical solutions for the critical section problem in distributed Systems, In Proceedings of the 9th ACM Symposium on Theory of Computing. IEEE Computer Society Press. 1977. 91–97. 
[89] G. F. Pfister. In Search of Clusters. Prentice Hall. 1998 (2nd edition). 
[90] A. Rényi. Probability Theory. Akadémiai Kiadó/North Holland Publ. House. 1970. 
[91] S. N. Ross. Simulation. Academic Press. 2006. 
[92] N. Santoro. Design and Analysis of Distributed Algorithms. John Wiley & Sons, New York. 2006. 
[93] A. Segall. Distributed network protocols. IEEE Transactions on Information Theory. IT-29(1). 1983. 23–35. 
[94] A. Silberschatz , P. Galvin , G. Gagne . Applied Operating System Concepts. John Wiley & Sons. 2000. 
[95] D. Sima , T. Fountain, P. Kacsuk . Advanced Computer Architectures: a Design Space Approach. Addison -Wesley Publishing Company. 1998 (2nd edition). 
[96] A. S. Tanenbaum . Computer Networks. Prentice Hall. 2004. 
[97] A. S. Tanenbaum . Modern Operating Systems. Prentice Hall. 2001. 
[98] A. S. Tanenbaum , M. van Steen . Distributed Systems. Principles and Paradigms. Prentice Hall. 2002. 
[99] M. Taqqu , W. Teverovsky, W. Willinger. Estimators for long-range dependence: an empirical study. Fractals . 3(4). 1995. 785–788. 
[100] G. Tel . Introduction to Distributed Algorithms. Cambridge University Press. 2000 (2nd edition). 
[101] Top 500 Supercomputer Sites. http://www.top500.org . 2007. 
[102] L. G. Valiant . A bridging model for parallel computation. Communications of the ACM . 33(8). 1990. 103–111. 
[103] Web Dictionary . http://pcp.wub.ac.be/ASC.html . 2007. 
[104] W. Willinger, V. Paxson. Discussion of {. The Annals of Statistics . 25(5). 1997. 1805–1869. 
[105] W. Willinger, M. Taqqu , R. Sherman, D. Wilson. Self-similarity through high-variability: statistical analysis of Ethernet LAN traffic at the source level. IEEE/ACM Transactions on Networking . 5. 1997. 71–86. 
[106] W. Willinger, D. Wilson, W. Leland, M. Taqqu . On traffic measurements that defy traffic models (and vice versa): self-similar traffic modeling for high-speed networks. Connections. 8(11). 1994. 14–24. 
[107] W. Winkler . Dependent percolation and colliding random walks. Random Structures & Algorithms. 16(1). 2000. 58–84. 
[108] M. J. Wolfe. High Performance Compilers for Parallel Computing. Addison Wesley Longman. 1996. 
[109] J. Wu, S.. On cost-optimal merge of two intransitive sorted sequences. International Journal of Foundations of Computer Science. 1(14). 2003. 99–106.