With the development of increasing node numbers in High Performance Computing (HPC), i.e. strongscaling, we get the problem that many HPC applications become communication-bound. Thus, there is a need for optimizing communication. Apart from improving the communication protocols, one can work on problem-specific communication patterns. A commonly used, eficient communication pattern for the simulation of homogeneously distributed particle systems is the Shift communication. Aiming to use the Shift for our own molecular dynamics simulation library, FMSolvr, we modeled it with timed Petri nets in order to get an estimation of the communication time. For the time analysis we use MaxPlus Algebra. We put our model's timing predictions up against an implementation of the Shift and found that our model is correct. However, our approach is over-modeling the problem.
small cubic boxes (spatial decomposition with an octree2). With distributing inter-dependent sub-regions of the simulation space over diferent nodes, the question of eficient data exchange arises.
At first, the time cost for computation in High Performance Computing (HPC) was much
higher than for communication. Simulation was computation-bound. But processors became
faster, calculation cheaper, and HPC workstations grew into HPC clusters, until today we
have the problem that strong-scaling leads to communication-boundness regarding inter-node
communication. Thus, communication is a topic of optimization again. There are developments
in the field of communication protocols such as LCI [
Since the Shift appears to be dominant in this field of research, we wanted to examine if it is
suitable for FMSolvr. Given certain system properties and constraints, we modeled the Shift upon
those properties using Timed Petri Nets (TdPN)3. Petri nets appear to be a good choice because
other than for example CSP (Communicating Sequential Processes) they can easily be evaluated
regarding timing behavior with the help of Max-Plus Algebra (MPA, also (max,+) algebra).
Having a formal model of the Shift will allow us to compare the Shift to other communication
schemes. One such comparison that is of special interest to us, is the Team Shift developed by
Driscoll, Geroganas, and Koanantakool [
As for the structure of this work: Section 2 will explain the mechanics behind MDS and the meaning of the terms used in the introduction such as near-field, Section 3 will introduce the Shift communication, and Section 4 will walk the reader through the necessary theory of TdPNs and MPA. Then, Section 5 explains how we modeled the Shift based on our system properties and constraints, and Section 6 will compare the model predictions with the measurements from the implementation. Finally, Section 7 will give a discussion of the methods and results as well as an outlook to our future work.
This section will introduce the reader to the basics of MDS: how to distribute the particles over
the computation nodes in Section 2.1, and the specifics of near-field interactions, which we are
focusing on with this work, in Section 2.2.
2Octree means that we split the simulation space into eight big boxes, then we split those big boxes into eight
smaller boxes, which we split again into eight smaller boxes, and so on, until we have our required granularity of
splitting up the simulation space.
3We choose the abbreviation TdPN because TPN stands for Time Petri Nets, which are a diferent type of
timedependent nets [
In MDS there are three ways to distribute particles over nodes: atom, force, and spatial
decomposition [
Spatial decomposition means that we take our simulation space, divide it into small boxes,
and distribute these boxes over our computation nodes. Of course, there are diferent ways to
divide the space. When one handles sparse, clustered particle systems, one would like to use
methods similar to Barnes-Hut [
When a particle is moving among other particles, it is subject to and the cause of Coulomb forces (aka. electrostatic forces). The Coulomb force is strong at a short range but its influence decreases with growing distance (∼ 1/2), which means that when two particles move apart, the individual influence that one particle has on another slowly melts into the background noise of other particles’ force fields. The distance at which the individual particle can no longer be distinguished from the background is the end of the near-field, this distance is called the cut-of radius cut.
In order to supply for each particle the data of particles within its cut-of radius, a node must import the data of particles, which surround its box within the distance of cut. Figure 1 shows the import area in relation to the cut-of radius. The node holding the dark-blue box must communicate with every node that holds a (partly) light-blue box. This means that there is no diference in the number of communication partners between Figure 1a and 1b. Only with a radius bigger than the box length, the number of communication partners changes.
In order to determine the number of communication partners of a node, we use = ⌈cut/⌉
to calculate the number of boxes covered by the import area in one direction. For Figures 1a
and 1b is one, and for Figure 1c is two. Then the total number of communication partners
of one node in a three-dimensional simulation is defined by
= (2 + 1)3 − 1
(1)
4For more details on atom and force decomposition as well as their communication complexity, please refer to [
With a growing , the number of communication partners increases quickly. In order to
ensure an eficient exchange of data, one needs an eficient communication algorithm. Here
is where Shift communication enters the game. A predecessor of the Shift is well-known in
parallel programming [8, Ch.4.2], the MDS version was introduced by [
Figure 2a serves to give the reader an understanding of the surroundings of a box. The whole simulation space looks like a grid of boxes where this 3 × 3 grid is a sub-grid of the whole. To describe the Shift, we want to focus on the blue and yellow box, so the rest of the grid is made invisible (Fig. 2b). We will use = 1 for the description.
In the first step of the Shift the node exchanges its data with its direct neighbors Left and Right (Fig. 2c). At the end of this step the node holds all the row data within = 1. For the second step it combines the Left and Right data as well as its own into one message and exchanges it for the package of its Front and Back neighbors (Fig. 2d). Now it hold all the knowledge about the plane. In the third step of the Shift the accumulated data of the plane is exchanged for the plane data of the Up and Down neighbors (Fig. 2e). So, in order to get all the required data for = 1, a node must only exchange six messages with its direct six neighbors.
In case of > 1, the neighbors help daisy-chaining the data further down the row. This means, the Right neighbor would hand the node the data of the second neighbor to the right and the Left neighbor of the second left neighbor. The data of the third neighbor to the right would pass through the second neighbor to the right and the Right neighbor before it arrives at the target node. So, a node must exchange a total of 6 packages to receive all the required data (this means 12 messages in total). However, the reader might have noticed that a node must wait for specific data to arrive before it can either pass it on or start the next phase. This is where time dependency between nodes comes into play.
In order to capture the time dependencies between nodes, Petri nets with timed transitions can be used. A timed transition here represents a message passing process. This section will do a brief recapture of Timed Petri Nets (TdPNs) and explain how to transfer TdPNs into the mathematically evaluable form of a max-plus matrix.
A marked Petri net is characterized by a five-tuple = (, , , , 0) where and are
places and transitions respectively, are the edges connecting a place with a transition or a
transition with a place, are the multiplicity (weights) assigned to edges, and 0 is the initial
marking. If we assign each transition ∈ a duration ∈ , we gain a TdPN, which is a
six-tuple of = (, , , , 0, ) where is called the duration function6. [
To make a practical example we shall look at a node that wants to send data to another node. The sending process can be captured by a timed transition, while the places will represent certain internal states of the node. In this case, we have two states: ‘before sending’ and ‘done sending’. Figure 3a depicts the net with its initial marking, and Figure 3b shows the marking after the node has finished sending. The time that passed between state one and state two is 6In our case is not a function but simple time values for each transaction. 2 the time assigned to the transition , denoted by the number in angle brackets.
The Max-Plus Algebra (MPA) comprises of two operations: the maximum operation ⊕ and the addition ⊗ . They are defined as follows: ⊕ = ⊕ = max{, }
⊗ = ⊗ = + .
The addition (⊗ ) has priority over the maximum (⊕ ) if not ordered otherwise by brackets. MPA also knows a neutral element for the maximum, the zero element = −∞ , and a neutral element for the addition, the identity element = 0, for which holds:
Equivalent to the sum ∑︀ and the product ∏︀ in linear algebra, MPA defines an accumulative maximum and sum: ⊕ = ⊗ = .
⨁︁ = max{0, 1, ..., } =0 ⨂︁ = 0 + 1 + ... + . =0
⊕ = [ ⊕ ] . ⊗ = [︃ ]︃ ⨁︁ ⊗ . =1 And just like in linear algebra, MPA can be applied to matrices. Taking two matrices = [ ] and = [ ] of size × , the maximum is an element-wise maximum defined as Taking a matrix = [ ] of size × and a matrix = [ ] of size × , the MPA matrix plus operation is defined as We also find a neutral element for the MPA matrix plus operation, the identity matrix , which is defined as Petri nets. For a matrix A of size × and a vector ⃗ of length , it is defined as In regard to vectors, the MPA matrix-vector addition (⊗ ) has an important role in analyzing ⎡
= ⎢⎢⎢. . .
⎤ . .. ⎥ . ⎥⎥ . . .. ⎥
⎦
⊗ ⃗ =
⨁︁ ⊗ .
=1
]︃
For more detailed information or as a reference for the next section, please find [
MPA is a tool to transform TdPN into a mathematically evaluable matrix form. For the
transformation we use a matrix of size | | × | | (| | is the total number of places in the net), and if
there is a transition that connects place with place , the matrix element gets the time
value that the transition takes, otherwise it gets the value = −∞
[
Let us extend the example of a node sending a message. Now, we have two nodes that want to send each other a message. Let us assume they cannot send and receive at the same time, then one must send after the other. If one node wants to send, it needs the “ok” of the other node, so in order for the transition to fire, both nodes must provide a token: the sending node to signal that it has data ready, the receiving node to signal that it is ready to receive that data. The respective Petri net is shown in Figure 4a. Place 0 is the internal state of node 0 that means ‘data ready’, and 1 is 0’s state ‘done sending’. Place 2 is the internal state of node 1 that means ‘ready to receive’ while 3 is 1’s state ‘done receiving’. Equivalently, the sending process of node 1 is shown in Figure 4b. When we combine the two nets into one TdPN, we obtain Figure 4c, the net of a classical ping-pong message exchange. One can see that the two nodes have timing dependencies.
Our final Petri net has six places, so our connection matrix is of size 6 × 6. When we fill in the connections between places like described in the first paragraph of this section and assume that = 1, we get the connection matrix = ⎢⎢⎢⎢ ⎡ ⎢0 ⎢⎣0 4 4
⎤ 4 0 0 ⎡ ⎤ ⎢
1 1 ⎥ 4 ⎥⎥⎥⎥ = ⎢⎢⎢⎢ 1 1 ⎥⎥⎥⎥ . ⎢ ⎣1 1 ⎦⎥
1 1 ⎥ ⎦⎥ 1 3 1
What we would like to find now is the longest path. This can be achieved by MPA adding ( ⊗ ) the connection matrix with itself. is able to return the longest path of length (with ∈ N) and is defined as = ⊗ ⊗ ... ⊗ = ⨂︁ .
=1 * = ⊕ ⊕ 2 ⊕ ... ⊕ But we want the overall longest path, which is returned by the star matrix * defined as where ≥ | | − 1, because this guaranties that all possible connections are included even in the case of a node ring. Now we can MPA add (⊗ ) our input vector to the star matrix and get as the result the times when the places of the net are reached. Let us assume that node 0 starts sending immediately and node 1 is immediately ready to receive, then our input vector is a zero vector ⃗ = ⃗0. By performing MPA matrix-vector addition * ⊗ ⃗ we get an output vector ⟨0⟩ 0 4 1 (c) ⃗ that tells us what state is reached after what time7: Vector ⃗ shows that the states represented by place 0 and 3 are immediately reached (due to our assumption that both nodes are ready at time 0), that the states that are represented by place 1 and 4 (‘0 done sending’ and ‘1 done receiving’) are reached after one time unit, and that the states represented by place 2 and 5 (‘0 done receiving’ and ‘1 done sending’) are reached after two time units. If we compare this to our TdPN in Figure 4c, the reader will ifnd that the numbers are congruent with the paths through the net.
In this section we will give a detailed description of how we modeled the Shift on the basis of our system properties. Therefore, we will summarize our system properties and additional constraints and then show the Petri net model for the first dimension (row) and how to extend it to the second and third dimension (column and tower).
The communication on the JURECA-DC cluster is location invariant8. This means that it does not matter if two communication nodes are placed in the same rack or on diferent racks. The physical distance of nodes has a neglectable impact on the latency of the exchanged message.
For our implementation, we constrain the system by demanding that sender and receiver
must have a handshake for the transfer. For small messages below 256 kB9, MPI does not
demand a handshake but uses the eager protocol (see OpenMPI eager [
First, we will look at the Petri net for a one-dimensional Shift communication. As we already introduced with the examples in Section 4, each node is represented by a number of places that stand for diferent states which the node reaches with progressing time (see example: ready to send → done sending = ready to receive → done receiving). In order to send a message, a sending node needs the agreement of the receiving node that the latter is ready to receive data. Together they can enable the transition that represents the sending process. The firing of the transition takes a certain time and produces a token for the sending node, which means that the sending is finished, and a token for the receiving node, which means that the data is received.
In order to avoid non-determinism (for now), we dictate a communication pattern: the nodes can either start with all transactions to the right and then switch to finish the transactions to the left or vice versa. Figure 5 shows the first pattern for = 1. Places, which are in a horizontal line, represent the states of one node . The transactions are indexed with the number of the preceding place of the sending node; the durations have been omitted since we have location invariance, i.e. all transitions take the same time. For the nodes 2 and 3 the presented section of the whole net is suficient for a MPA analysis of the time dependencies for = 1. The number of nodes required for the analysis is increasing quickly with a growing . It can be put into a formula in the following way10:
net = 8 − 2 .
For the other two phases (column and towers), the nets look identical to the first dimension, only the transition times change due to higher message load. One net of the column dimension requires net neighboring nets from the row dimension. Equivalently, one net of the tower dimension will need just as many neighboring column TdPNs. The diferent phases can be handled independently with the output of the previous phase being the input for the next one. For this work we limited our analysis to the first phase to see whether the model can reflect the actual time behavior.
In this section we will present the results of the theoretical analysis with MPA and the actual measurements from the implementation.
For the MPA analysis we used Python to generate and analyze the matrices. The input for the
transition durations is the average message passing time between nodes of the cluster because
we wanted to see whether measuring the average latency and feeding it into the matrix would
be enough to gain an accurate prediction for the communication time. With an OpenMPI
pingpong message (using MPI_Ssend) we measured the average message passing time between
nodes using 1,000 ping-pongs. The results of the ping-pong measurements for diferent message
sizes are listed in Table 4 in the Appendix.
10The formula results from analyzing the nets for ∈ [
The implementation of the Shift was done with OpenMPI 4.1.2 (GCC 11.2.0) on the JURECADC cluster (InfiniBand network). We used MPI_Ssend for all send routines in order to force a rendezvous for each message. The amount of part-taking nodes appears to be neglectable (see Tables 5-7 in Appendix), so we will focus on the measurements performed with 64 nodes to get the full spectrum of -values.
The results from the MPA analysis compared to the implementation measurements can be found in Tables 1-3. The predictions are within the standard deviation of the measurements, so the model can be considered correct for the range of -values that we used. However, we 16 6 27 22 17 12 7 28 23 18 13 8 28 18 8 29 24 19 14 9 can also see that the prediction is more accurate for message sizes of 10 and 1,000 bytes than for 100,000 bytes.
This section is the conclusion of our work and discusses our approach to modeling, our model, our measurements, and our result. Last, it gives an outlook on what we want to do next.
7.1.1. Measurements We noticed that the predictions are less accurate for the message size of 100,000 bytes. Our thought is that for messages of this size we can see network saturation efects which we do not have when we perform a single ping pong measurement. A way to improve the model here would be to create network noise while measuring the ping pong.
Another thing we noticed it that if we limit the ping pong measurements and the Shift measurements to the 20% smallest values (the best cases), the prediction becomes stunningly accurate (diferences of a few hundred nanoseconds). However, we cannot name a reason for the 20% threshold, so this needs further evaluation. 7.1.2. Model Our model can be considered correct. However, the MPA matrices are growing exponentially with the -value, which makes them time consuming to evaluate. Moreover, we are overmodeling our system. Our model incorporates more detail than what is necessary for the analysis. The ability to assign diferent times to all transitions in order to model slower and faster nodes is not necessary due to location invariance of the cluster. Here, we reduced the model’s complexity to a simple calculation of average message latency multiplied by the number of message exchanges/number of rendezvous.
By using a duration function that is dependent on the bandwidth and trafic in the network instead of duration values, a model like ours could be helpful to analyze the influence of network congestion. It would also be more suitable for analyzing the timing behavior that results from the lack of location invariance, i.e. diferent duration values per node pair.
We chose Petri nets in the first place because of their ability to include time easily and handle non-determinism, which we will need in the future once we switch to an implementation with MPI_Isend. However, TdPNs do not seem like a suitable approach for our problem or our cluster due to the reasons mentioned above.
For the presented analysis, we forced MPI not to use the eager protocol but the rendezvous protocol, which is increasing the latency. A future model will incorporate the eager protocol. Also, non-determinism in the form of e.g. MPI_Isend will be included since it might be able to reduce communication time by introducing flexibility in the case where a communication partner is not instantly available.
We thank Andreas Beckmann for assistance with the implementation of the one-dimensional Shift. This work is supported by DFG project FMhub, project Nr. 443189148.
= 64 [ns] = 32 [ns] = 16 [ns] = 8 [ns] = 4 [ns] 1 2 3 4 5 6 7 8 9 10 0 5 10 50 100 250 260 500