I Introduction
Traditional seismic acquisition systems employ cables to efficiently transmit data from the geophones (GPs) to the data center (DC). Telemetry cables, although efficient and reliable, are a burden in terms of cost and weight. With surveys growing larger in scale, costs of deployment and maintenance become proportionally higher. Furthermore, environments with complex terrains render deployment of cabled systems impractical. Hence, there is an increasing trend in oil exploration companies to shift to wireless technology in order to avoid these difficulties. Aspects for the transition to wireless acquisition include proper choice of wireless technology, the need for efficient communication protocols and reliable data transfer. Driven by the trend towards wireless seismic surveying, several authors studied cableless surveys and wireless architectures in the literature. Savazzi et al. [1] proposed to use ultrawide band wireless technologies in land seismic acquisition systems, where a hierarchical architecturebased, shortrange and longrange communication were developed for two stage data delivery. Freed [2] argued for nodal systems, where geophones receive seismic data continuously and store them for future retrieval. This eliminates the need for realtime data acquisition especially with a substantial number of channels. Ellis [3] proposed a hybrid system that combines cable and cableless acquisition systems. Reddy et al. [4] proposed a network architecture based on the IEEE 802.11af standard. The proposed architecture offers a high transmission data rate by operating in the TV white space. Most common wireless sensor network (WSN) topologies are tree, mesh and star topologies [5, 6, 7, 8, 9]. In a tree topology, the nodes are organized in a hierarchical structure. The root node serves intermediate nodes, that in turn serve other lowerlevel nodes. In a star topology, the central node is connected to all the other nodes. On the other hand, in a mesh topology, all the nodes are connected to each other through onehop or multihop distances. The structure adopted in this work is an enhanced version of the tree structure where the gateways (GWs) act as intermediate nodes between the GPs and the DC. In contrast to the conventional tree topology, any GP is assumed to be capable of transmitting its data to all GWs, rather than to a single GW only (as depicted in Fig. 1). This relaxes the constraints on the sumrate achieved by the GPs, which is essential in the formulation of the first problem related to sumrate maximization. To the best of our knowledge, this structure and the objective, i.e. sumrate maximization of the GPs have not been studied in a previous work. Allowing the GPs to send their data to all the GWs encourages random placement of the GWs. Hence, drones can be used as GWs to receive and store GPs’ data. This is an improvement compared to the work in [10], where GWs are assumed to be in fixed position. Furthermore, in contrast to the previous works, and building on the work done in [11], the focus of this paper is on the theoretical aspect of the geoseismic network, which was found to be underrepresented in the literature. Study of upper bounds on the network sumrate provides us with knowledge on the transmission capabilities of the geoseismic network. It also presents a valuable insight on how the decoding and selection processes at the GWs could affect the achievable individual rates of the GPs as well as the combined sumrate. For this purpose, comparisons between the case of optimizing the set of GPs to be decoded at each GW along with their decoding order, in order to maximize the achievable sumrate and the case of nooptimization are presented. In this work, the wireless system capability in terms of sumrate is optimized. Informationtheoretic bounds are used to gain insight into the achievable rates of the system. For this purpose, successive interference cancellation (SIC) is assumed at the (GWs). SIC is performed by successively decoding one of the GPs signals and subtracting it from the original signal, and then repeating the process for the rest of GPs, and thus reducing the interference after each decoding process. It is known that SIC achieves the capacity of the multiple access channel at each GW [12].
Furthermore, in this paper, the data transmission from the wireless GWs to the DC is considered. The DC uses SIC to decode the signals of the GWs and achieves the capacity of the multiple access channel (MAC). Although the problem in this work is formulated for seismic acquisition, it can be applied in general to other WSNs by taking the specifications into consideration. In summary, the contributions of the paper are as follows:

The problem of maximizing the sumrate of the GPs using SIC is formulated and three main metaheuristic algorithms are proposed to solve it.

The problem of data transmission from the GWs to the DC, for GWs with small buffers, is considered for two different objective functions and solved.

The problem of maximizing the weighted sum of the GWs rates is solved.
The rest of the paper is organized as follows. In Section II, the system model is presented. In Section III, the mathematical formulation of the proposed optimization problems is discussed. The proposed integer programming algorithms to solve the first optimization problem are presented in Section IV and in Section V findings and discussions are shown. Finally, in Section VI the concluding remarks are presented followed by an appendix.
Ii System Model and Problem Description
In this work, the deployment of wireless GPs in a gridlike network to cover the whole survey area is considered, where wireless GWs are deployed between the GP lines. The purpose of GWs is to collect data from the GPs, possibly store them for some time, and then forward the collected data to the DC for processing [13]. A typical land survey deploys 20,000 to 30,000 GPs over a large area on the order of 20 km^{2}. A GP equipped with a 24bit analogtodigital converter (ADC) generates sampled data at a rate of 144 kbps, with a sampling interval of 0.5 ms. The aggregate data rate for a typical survey can be about 4.32 Gbps, when all the GPs are active. A plane view of an orthogonal geometry survey is shown in Fig. 2, where two GWs are placed between the source and receiver lines to harvest the data from the GPs. Fig. 2 also depicts the first and second stages considered in this work, where the first stage is the data transmission from the GPs to the GWs and the second stage is the data delivery from the GWs to the DC. The GWs may have some buffering capability to store the collected data, then it is sent to the DC for processing.
In this paper, a system that consists of GPs and GWs is considered. The channels between the GPs and the GWs are assumed to follow a Rayleigh block fading model. This model is convenient for desert environment and seismic acquisition applications, where the GPs are placed in fixed positions and the scatterers are stationary. The channel matrix has dimensions. Some assumptions about the system model are listed below:

The channels between all the transmitters and receivers are constant for a fixed interval of time and change independently between time intervals.

The transmission period is a fixed time slot which is of equal width for all transmission sessions.

The channel state information (CSI) of all channels is available at all the GWs.

Each GP transmits with a fixed power, .
The GWs are using SIC, where they can choose to decode some of the GPs and consider others as interference to maximize the total sumrate of the network. Due to the use of SIC at the GWs, each GP will have different constraints on its rate, coming from the different GWs that decode its signal. The rate of each GP, , must satisfy all these constraints, which might lead to reducing its rate. Therefore, sometimes it is better not to decode the GP data at some GW as this will relax some of the constraints and help increase its rate. The Shannon normalized capacity is given as,
(1) 
where is in bits per second per Hz (bps/Hz), is the channel gain value taken from the channel matrix
. The additive noise is assumed to be Gaussian distributed with a variance
. If interference is present, it will appear in the denominator added to . Each GP in a multiple access channel can transmit with a rate. For multiple GPs transmitting simultaneously, a multiple access region is characterized by all the upper bounds on the rates. It is assumed that each GW will decode the GP signals in a decreasing order, starting from the signal with highest signaltonoise ratio (SNR) down to the signal with the minimum SNR. This will feature fairness to the weaker signals, while providing the maximum sumrate
[12]. To illustrate how decoding selection at the GWs can affect the total sumrate, consider a small representative network of GPs and GWs, , as shown in Fig. 3. The channel gain matrix is generated following Rayleigh random density function and is given as,In this case, by going through all possible decoding order combinations, the one that gives the maximum sumrate (denoted as maxsumrate) is illustrated in the first stage in Fig. 3, where GW_{1} decodes GP_{2} then GP_{1} (solid arrow), while treating GP_{3} as interference (dashed arrow), and GW_{2} decodes GP_{2} then GP_{3}, while treating GP_{1} as interference. To see the rate expressions for both cases please refer to Table I and Table II. The rate expressions can be verified from (2) by considering the individual constraints on the rates. By substituting the values, the total sumrate for the maxsumrate case is found to be bps/Hz. Now compare this with the decodeall case, where all the GWs will decode all the GPs without doing any kind of optimization. In this case, it can be found that bps/Hz. The individual rates for the two cases are shown in Table III. As can be seen, decoding a GP at multiple GWs puts more constraints on the rate of that GP, which might result in reducing its rate. However, this could present opportunities to increase the rates of other GPs, as is the case in maxsumrate where the rate of GP_{2} was decreased to 0.3668 bps/Hz, but it helped to increase the rate of GP_{1} to 3.011 bps/Hz. Therefore, for some GPs, it might be better not to decode them at all GWs, while for others it will improve the sumrate if they are decoded at all GWs. And this is the optimization problem that is considered in the first stage presented in the paper, i.e. how to optimize the set of GPs to be decoded at each GW so as to maximize the network sumrate.
GP  Maxsumrate 

1  
2  
3 
GP  Decodeall 

1  
2  
3 
Following that, the second stage is considered, where the delivery of data to the DC is studied. In the second stage, two optimization problems are considered based on the length of the GW buffer. For smallsize buffers, the goal is to transmit all the stored data in the transmission slot. In this case, two objective criteria are proposed. The first criterion considers the optimization of the total power expenditure at the GWs. Each GW is assumed to be constrained with a maximum transmission power which helps to save its battery. The second criterion considers power fairness between the GWs. In this case, only the maximum power is considered as the variable for optimization. It is proved that for both objective functions the problem is convex and, hence, convex optimization techniques can be used. The solution is obtained by finding the optimal values of power and the optimal decoding order at the DC. Furthermore, an optimization problem is formulated and proposed for largesize buffers at the GWs. In this case, the goal is to maximize the weightedsum rate of the GWs. Each GW might have a different length of stored data in its buffer and, hence, it might be necessary to give priority to those GWs with long queues. Giving priority can be mathematically done by increasing the weight of those GWs in the objective function of the optimization problem. The larger the weight of the GW, the more probable that it will have better rate. The problem here is a power allocation problem where the total power is constrained by a maximum value. In this case, the GWs’ rates are also constrained by all combinations of upper bounds on the achievable rates. If the survey area changes, the solution has to be updated by solving the relevant problem. As can be seen, this work can be generalized and applied to other WSNs, and not only to seismic acquisition systems.
GP  Rate (maxsumrate)  Rate (decodeall) 

1  
2  
3 
Iii Problems Formulation
The mathematical formulation of the proposed problems is presented in this section, where the first problem considers the sumrate maximization of the GPs in the first stage. The following two problems consider the second stage where the delivery of data to the DC is discussed.
Iiia Sumrate Maximization
The sumrate maximization problem is formulated as an integer constrained optimization problem as follows:
(2a)  
(2b) 
where in (2a) is the actual rate of the GP and in (2b) is an element in the set of all combinations , . Each link has an index that takes the value 1 if the GW decides to decode the signal of the GP, otherwise . The set of constraints in (2b) is a generalized form of the multiple access channel constraints (with multiple destinations), and it gives all the combinations of upper bounds on the achievable rates for all the GWs. The numerator inside the in (2b) accounts for the sum of GPs’ powers decoded at the GW, while the denominator represents the noise power added to the interference from all the GPs not decoded at this GW. Although allowing a GP signal to be decoded at multiple GWs can introduce redundant data, it is essential to maximize the sumrate of all GPs. This is because decoding a GP’s signal at multiple GWs could reduce the interference seen by other GPs’ signals and improve the overall sumrate. The optimization in (2) is performed on the links indices; ’s. This renders the problem an integer programming problem, where the search space is comprised of all decoding combinations and its size, , is given by the following formula:
(3) 
To give an idea about how rapidly the size of this search space increases, consider the network shown in Fig. 3, where and then from (3), combinations, while when the parameters are increased to and , there are
combinations. Going through all the combinations is a challenging task and consumes huge resources in terms of time and computations. Therefore, it is necessary to look for heuristic and metaheuristic methods to solve this problem. Heuristic algorithms are problemspecific lowcomplexity search protocols that provide a suboptimal solution by exploring the search space in a faster way. On the other hand, metaheuristics are highlevel algorithms that are applicable on a widerange of problems. They are generally used to solve complex problems where it is hard to develop a specific algorithm
[14]. In Section IV, three sets of metaheuristic algorithms are introduced to solve the problem in (2). Specifically, two variations of particle swarm optimization (PSO), two variations of ant colony optimization (ACO), and the simulated annealing (SA) method are discussed.
IiiB Mintotal Power
The mintotal power considers the case of data delivery from the GWs to the DC with minimal total transmission power. It is formulated as an optimization problem, where the variables of optimization are the individual transmission powers of the GWs, and the decoding order at the DC. The importance of the decoding order is that it affects the individual rate for each GW. The problem is formulated below,
(4) 
where and is the channel gain of the link between the GW and the DC. Here, refers to the normalized data rate in bits/s of the GW such that where is the amount of data stored in the GW, in bits, and is the transmission time, which is a fixed time slot for the GWs. Given the constraints in (4), the GWs powers are optimized, and from the optimal power values, the optimal decoding order can be constructed. The objective function, , is a linear function and is therefore a convex function. The first set of constraints on are linear constraints too. Thus, each one of them forms a convex set. The upper bounds on the normalized data are of the type where , and are constants. Since the function is a concave function, these set of constraints represent a convex set. Therefore, the problem is convex, and the global optimal solution can be efficiently obtained using convex optimization techniques such as Newton’s method. The optimal decoding order at the DC for the problem in (4), that will minimize the total power expenditure, is to decode the signals of the GWs based on the channel gain values in descending order. The proof is given in Appendix 1.
IiiC Minmax Power
In minmax power, the power fairness between the different GWs is considered. The objective of this problem is to find the decoding order at the DC that will minimize the maximum required transmission power at any of the GWs. The problem is formulated as a minimax problem, where the variables of optimization are the individual transmission powers of the GWs. In the case where the buffers are limited, the GWs must transmit all the stored data before the arrival of the new data. The set of all possible MAC achievable rate combinations form the constraints on the normalized data rate. For practical consideration, each GW is assumed to have a maximum power of . The problem can be formulated as:
(5) 
Similar to mintotal power, the convexity of the minmax power can directly be proved. Once the optimal powers are achieved, the decoding order can be found from the active constraints. A possible scenario is that the DC uses timesharing between two or more decoding orders in order to achieve the desired rates and at the same time minimize the maximum transmission power. For timesharing between two decoding order combinations, the sumrate will be the maximum value and its point will lie on the line connecting between the two corners in the capacity region characterized in (5). To achieve any point on the line, the DC will alternate between the two corners by giving a weight to each decoding order characterized by the corner. This may, sometimes, prove useful to achieve power fairness, by working at a point, where the maximum power of the GWs will be less than the value achieved by working at any of the individual decoding order combinations. A system of linear equations is solved to find the exact percentage of time for each decoding order at which the DC should employ successive decoding, to achieve the desired rate point on the MAC.
IiiD Max weightedsum
The objective of the max weightedsum problem is to maximize the weighted sum of the rates of the GWs. Each GW will enjoy a rate that is proportional to its buffer size. The longer the queue of data in the buffer, the larger is the weight in the objective function. A constraint of maximum power is held on the total power of the GWs. Hence, the problem is to allocate the powers for each GW such that the objective function is maximized. A set of constraints, that is comprised of all MAC achievable rate combinations, is also imposed on the objective function. The problem can be formulated as:
(6) 
where is the normalized weight of the buffer and . The total power , is the sum of all powers of the GWs. As can be seen in the objective function, the larger the queue of a buffer of a certain GW, the more weight it will have. In this problem, the queues of the buffers are assumed to be large. This is because the time slot is assumed to be equal in size for all transmission sessions. When the stored data in the buffers is small in size such that all the data can be sent within the time slot, a different problem will appear. The GW that has higher rate will finish sending its data faster and then it will have zero rate. When this happens, the upper bounds on the achievable rates will be altered as the rate of some of the GWs will be zero at some point during the transmission period. This renders the problem very complex to solve. Therefore, the focus here is on the case of long buffers filled with large queues of data, such that all upper bounds on the achievable rates hold throughout the transmission time.
Iv Optimization Algorithms
In this section, the algorithms to solve the optimization problem in (2) are adopted and conclusion is made based on performance and complexity.
Iva Particle Swarm Optimization
Particle Swarm Optimization (PSO) uses simple velocity and position equations to update the particles positions which correspond to the possible solution [15, 16, 17, 18, 19, 20]. Velocity here represents the rate of change of a particle’s position in the search space. Two variations of PSO are used in this work for the problem in (2), namely, anglemodulated particle swarm optimization (AMPSO) and discrete particle swarm optimization (DPSO). These variations are able to solve the binary problem at hand.
IvA1 Ampso
The AMPSO algorithm has proven efficiency, especially when the dimensionality of the binary vector is large
[21]. The function given in (7), which is derived from anglemodulation theory, is optimized using PSO. It is given as,(7) 
where . , , , and are horizontal shift, maximum frequencies of the sine and cosine functions, vertical shift and a single element from a set of evenly separated intervals based on the number of bits, respectively. The parameters in (7) are substitued into the velocity equation to update their values and obtain the positions of the new parameters. The velocity equation is given below:
(8) 
where is the inertia factor, and are the local and global best values, respectively. and are the weights given to the local best and global best values, respectively. and
are random numbers that are drawn from a uniform distribution between 0 and 1,
. is the current optimized parameter, i.e. , , or . The position update equation is given simply as: .IvA2 Dpso
In DPSO, the indices ’s in (2) are considered as a matrix of bits, where each row represents the links between the GPs and a specific GW, and each bit (index) is optimized separately. The position of the optimized variable is restricted to 0 or 1. Velocity here is the probability that a bit is 1, where represents the index of the particle and represents the index of the dimension. The dimensionality of the bits matrix is . The velocity formula is given as [22]:
(9) 
where is the position of the particle in the dimension, is the bit position of the best performance the particle has achieved so far. is the bit position for the best performance of all the particles so far (global solution). Obviously,
is a realvalued number and to translate it to a probability, a logistical transformation called Sigmoid function is given as:
(10) 
To determine , a random value is generated from a uniform distribution in the range . Then, it goes through the following decision rule:
(11) 
A bound is usually set on velocity so that the algorithm keeps exploring the search space. The flowchart in Fig. 4 illustrates the detailed steps of the AMPSO and DPSO algorithms. It shows the variables initialization stage, followed by the process for each iteration. At a given iteration, each particle updates its potential solution based on the corresponding velocity equations and then updates the position. The local and global best solutions are then updated by comparison with the new solutions offered by the particles.
IvB Ant Colony Optimization
Ant Colony Optimization (ACO) is based on observations of the behavior of real ants[23, 24, 25]. In order to remember the path towards home and guide other ants, while searching randomly for food, they release a chemical substance called pheromone. The pheromone evaporates with time, which allows the ants to explore other parts of the search space and avoid being trapped in a local optimum. If there is more than one path that is taken by ants, the shorter path will hold pheromone longer since ants take less time walking through it. Thus, eventually all ants use this shorter path. The detailed steps of the general ACO algorithm can be followed in Fig. 5
. Firstly, an initialization stage is run by introducing the pheromone and simulation variables. The probability for each binary value is estimated by each ant for all possible paths. Then, ants construct the solution by generating random values and comparing them with the estimated probabilities. Each solution is used to find the sumrate and update the pheromone. Varieties of ACO algorithms have been developed in the literature.
IvB1 Ant System (AS)
In AS algorithm, all ants that successfully constructed a solution update the pheromone values at each iteration.
(12) 
where is the evaporation coefficient, is number of ants, is the pheromone value, is the index of possible outcomes, , is the dimensionality of the system, , is the new value of pheromone for a link with possible outcome and is the pheromone deposit, made by the ant , based on its sumrate value. When constructing a solution, the probability of ant choosing for the bit is given by:
(13) 
where is the set of feasible components, i.e. and is the heuristic information, which depends on the problem at hand. It is usually linked to a priori information about the problem. In our problem, is linked to the channels’ gains, as will be shown in Section V. and represent the weight given to the pheromone and to the heuristic information, respectively.
IvB2 MaxMin Ant System (MMAS)
In MMAS, only the best ant updates the pheromone values. Moreover, the pheromone update is bounded as follows,
(14) 
where is the maximum possible value of the pheromone and is the minimum possible value of pheromone.
IvC Simulated Annealing(SA)
Simulated annealing is based on the metal annealing in material science [26]. It uses restarts to avoid being trapped in a local maximum. It always accepts moves that improve the value of the objective function. If this move does not improve the objective function, then the algorithm may accept the move with some probability. This probability assumes its highest value at the beginning of the algorithm. Therefore, the algorithm is more willing to accept moves that do not improve the objective function at the beginning to explore the search space. As the number of iteration increases, that probability decreases. The willingness of the algorithm to explore new parts of the search space is signified by a temperature degree. Therefore, as it decreases gradually, the algorithm becomes less willing to trade a worse position by its current position. Moreover, the temperature value is equivalent to the number of iterations. This algorithm is illustrated in Fig. 6 where is the temperature value, is the difference in sumrate between the best solution so far and the next solution.
In Section V, the algorithms are tested and a conclusion is derived about them based on complexity and performance.
V Results and Discussion
To solve the first optimization problem, a limited computational budget is given to the proposed algorithms in the simulations. The computational budget is represented as the maximum temperature () for SA, number of particles/ants () and iterations () for MMAS/AS/DPSO/AMPSO. The constants used in aforementioned algorithms, (), , for AMPSO, for DPSO, (), () are set to , , , , and , respectively. These values are obtained using extensive simulation and based on [27]. In AS, equals the channel gain when and equals the average of the other channel gains for the same GW when . However, as will be discussed shortly, must be adapted when the problem involves a very large search space. In MMAS, a moderate range of is found to give acceptable results. Furthermore, a big range might cause convergence stagnation, as the algorithm will keep exploring without benefitting from favourable findings. The mean square error (MSE) for the algorithms is calculated with respect to the optimal solution. The optimal solution is obtained via exhaustive search (ES) which tries all possible decoding order combinations and finds the one that gives the maximum sumrate. Fig. 7 shows the MSE versus iterations for GPs and GWs. The size of the search space that the ES looks through is combinations. For the algorithms here, a computational budget of and is used. Furthermore, the results are averaged over 100 random channel realizations. As can be seen in Fig. 7, AS and SA achieve the best convergence results. MMAS and DPSO achieve comparable results. AMPSO, however, has a slower convergence and scores the worst results for the final MSE value. It should be noted that in all results, SA algorithm’s iterations are modified to match the computational complexity of other algorithms since it does not entail ants/particles. For example, in Fig. 7, SA result is modified by multiplying the number of iterations by , which is 30. Table IV shows a comparison between the performance of exhaustive search (ES) and the proposed algorithms for a small network of 8 GPs and 2 GWs. Here, three different levels of computational budget of the proposed algorithms are shown. Comparable results between all the algorithms can be seen, with DPSO achieving higher results for the first two cases and AS achieving better results for the third case where the computational budget is lower. To describe realistic surveys, number of GPs and GWs must be larger. Therefore, an example of and is considered. The number of decoding order combinations here is beyond the capabilities of any practical computer; combinations. Therefore, the use of the proposed algorithms with a limited computational budget becomes very crucial in such cases. Two different scenarios for how the GPs transmit data to the GWs are considered. In Scenario 1, the GP is assumed to have a small buffer. In this case, the GW can hold data only up to the next transmission session. On the other hand, Scenario 2 assumes that the GP is equipped with a large buffer to store data for future time transmissions. This gives the GPs the choice to transmit at a given transmission session or not based on its channels’ gains. However, sending the stored data later may affect the accuracy of the seismic information that reaches the DC, especially if the channels of those GPs are in a deep fade for several continuous transmission periods.
Scenario 1) (Realtime) where all the GPs always send data:
In this scenario the powers of all GPs are fixed. Here, it is more beneficial to decode all the signals sent by the GPs due to the existence of interference all the time. Fig. 8 shows the normalized data rate in bps/Hz versus iterations. The results are averaged over 50 random channel realizations. The channels follow Rayleigh fading with . Here, the used parameters are , and . As can be seen, AS achieves the best result followed by MMAS, DPSO, AMPSO and finally SA.
ES  MMAS  AS  DPSO  SA  AMPSO  Computational budget  
5.49  5.34  5.41  5.49  5.32  5.31 


5.49  5.11  5.17  5.21  4.77  4.76 


5.49  4.34  4.68  4.32  3.94  4.03 

Scenario 2) BufferAided GPs:
All the GPs are assumed to have buffers where they can store the sampled data before transmission. In this scenario the decoding order of the GWs can be optimized so that the sumrate is even improved over Scenario 1. This is achieved through deactivation of some of the GPs that generally have low channel quality. Fig. 9 shows the normalized rate versus iterations and the results are averaged over 50 random channel realizations. Here, the used parameters are , and . The best result is found to be achieved by AS, followed by MMAS, DPSO and then AMPSO, while SA displays the worst performance. Table V shows an example comparing the performance of nooptimization and the proposed algorithms for the two scenarios simulating a network of 100 GPs and 8 GWs. The nooptimization case is the special case where the GWs neither cooperate nor share information about the decoding order. It can be seen that Scenario 2 results give comparative advantage in sumrate over Scenario 1 for most algorithms. Compared to the small network of and , it can be seen that SA achieves a lower sumrate value compared to the nooptimization case.
Algorithms  Scenario 1  Scenario 2 

Nooptimization  4.33  4.33 
DPSO  6.32  7.51 
AMPSO  3.36  4.33 
SA  1.38  1.10 
AS (adapted)  7.77  7.87 
AS (not adapted)  3.62  4.66 
MMAS (adapted)  6.92  7.17 
MMAS (not adapted)  1.87  2.93 
Also, it is noticed that AS and MMAS achieve lower sumrate values when their heuristic information is not adapted. This is an interesting observation which indicates the effect of the convergence speed of metaheuristic algorithms on their capability to explore the search space in a limited time. Although AS and MMAS belong to ACO algorithms, which have been proved to always converge, the convergence speed is a problem that is still not fully explored [28, 29]. To overcome the issue of the algorithm’s stagnation, the effects of various parameters such as are studied. Except for , it is found through simulation that the effects of all parameters on convergence speed are minimal. However, for , the choice of the heuristic information will greatly influence the results. This is because biases the search towards regions that are expected to be promising. For Scenario 1, a weight is given to each GW that is proportional to the average of the channels gains associated between all the GPs and this GW. The higher the average for a certain GW, it is more likely to decode signals from GPs. This helps increase the sumrate because when a GW decodes GPs signals, it will successively eliminate the interference, which increases the total sumrate. For Scenario 2, it would be useful to turn off GPs which have all their links as weak channels. The reason is that those GPs have insignificant rates, and they cause interference for other users. Therefore, removing them reduces the unnecessarily limiting interference present in the network, which will help increase the rate of other GPs, while only losing the small rate that belongs to those GPs. To do this, the probability of having an index of 0 for the links associated with these GPs is increased. This makes it more probable that GWs will not decode GPs that have low channel gains. This measure is combined with the Scenario 1 measure, where each GW will have a weight to its links that is proportional to the average of the channel gains. The combined effect directs the search towards regions that have a higher sumrate and thus increases the algorithm performance with limited budget. Table VI shows the time taken by each algorithm to search for optimal solution for the same scenario considered in Table V. It can be shown that the algorithms have a complexity on order of . By comparing the time, it is found that SA achieves the best results by taking the least time. Following it, AS and MMAS give comparable results, while DPSO and AMPSO achieve the worst results. However, taking performance and complexity into consideration, it can be concluded that AS achieves the best results.
Importantly, our problem formulation cna find the total number of GWs that are required to serve a typical seismic survey. GP’s rate(144 Kbps[10]) is used for this purpose. AS and DPSO are used to find the average GP’s rate for various number of GPs in Fig. 10 and 11, respectively. It is noted that with AS, for example, 16 GWs can support up to around 360 GPs while maintaining the required average GP rate. However, when using DPSO, only up to 250 GPs are served.
In the following, two buffer sizes are considered at the GWs, namely, smallsize and largesize buffers.
Algorithm  Time in seconds 

DPSO  40 
AMPSO  69.95 
SA  6.67 
AS  13.16 
MMAS  14.5 
i. Smallsize Buffer
The two variations of this problem are solved using the fmincon algorithm in MATLAB. The fmincon algorithm utilizes a number of techniques to solve convex optimization problems. These techniques include interiorpoint method and sequential quadratic programming (SQP). These methods are based on approximation of the nonlinear constraints with linear constraints. A case study is considered for eight GWs. Random sets of stored data , and channel gain values are generated. The optimal power values for each GW as well as the optimal decoding order at the DC are found via simulation. The total power is also calculated. The same channels and stored data are considered for the two parts of the problem for ease of comparison.
Case Study
The normalized stored data rates in the buffer of each GW and the channel gain values between the GWs and the DC are given below,
where is in bps/Hz. By solving the mintotal power, the individual powers for the GWs are given as,
and the total power, . The optimal decoding order at the DC is found to be as follows:
GW_{3} GW_{7} GW_{6} GW_{5} GW_{1} GW_{4} GW_{8} GW_{2}.
However, for the minmax power, the solution is found by setting all of the powers of the GWs to be,
and the total power is, . In this case, the DC uses timesharing to achieve power fairness. All power values are equal and this is the best possible solution as any decrease of a power value prompts an increase in another to satisfy the achievable rate constraints. The only upper bound constraint that holds with equality is the sumrate constraint, and this shows that the DC uses timesharing between a number of decoding order combinations. A system of eight equations and eight variables is solved to find the exact percentage for each decoding order that the DC uses. The decoding order combinations can be chosen randomly but if one or more of the variables is negative, this indicates that the sumrate point is not in the targeted area. In this case, the decoding order associated with that variable has to be flipped. The decoding order combinations at the DC are found to be as shown in Table VII. Each decoding order combination is shown with the associated percentage of the time. This will minimize the maximum power used for transmission by the GWs.
Decoding order combination  Time percentage 

1 2 3 4 5 6 7 8  22.38% 
7 6 5 4 3 2 1 8  3.11% 
7 8 1 2 3 4 5 6  5.94% 
6 7 8 1 2 3 4 5  15.63% 
5 6 7 8 1 2 3 4  14.77% 
4 5 6 7 8 1 2 3  20.93% 
3 4 5 6 7 8 1 2  1.14% 
2 3 4 5 6 7 8 1  16.37% 
ii. Largesize Buffer
To consider the largesize buffer scenario, an example of a network of eight GWs is considered. The normalized stored data rate in the buffer of each GW and the channel gain value between the GWs and the DC is given below,
where is in bps/Hz. is generated following the uniform distribution between 30 and 100, i.e. . The problem is solved using the fmincon algorithm in MATLAB, where SQP is used. The maximum total power is set to be . From , the weights for each buffer are calculated to be
where . The individual powers for the GWs that maximize the objective function are found to be,
and the total power, . The optimal decoding order at the DC is found to be as follows:
GW_{7} GW_{3} GW_{1} GW_{8} GW_{4}, while GW_{2}, GW_{5} and GW_{6} are turned off.
The normalized data rate for each of the GWs is found to be,
In this case, GW_{2}, which has the lowest channel gain and lowest stored normalized data rate, is turned off. Moreover, GW_{5} and GW_{6}, which have relatively low channel gains and stored normalized data rate, are turned off.
Vi Conclusion
In this paper, the optimization of both stages of a geoseismic wireless acquisition system were considered. In the first stage, GPs send their sampled data to a few GWs. In this stage, the problem of maximizing the information theoretic bounds on the sumrate was formulated, and several metaheurisitic algorithms were proposed to solve it. The proposed algorithms are based on concepts from swarm intelligence and material science. These algorithms are anglemodulated particle swarm optimization, discrete particle swarm optimization, ant system, maxmin ant system and simulated annealing. Among them, the ant system algorithm proved to be the best in solving the sumrate maximization problem, considering the performance and computational complexity. Furthermore, it was shown that when the search space gets very large, convergence speed for the ant colony optimization algorithms can be increased by adapting the heuristic information accordingly. In the second stage, which is the data delivery from the GWs to the DC, two problems were explored based on the length of the buffer at the GW. For smallsize buffers, an optimization problem to minimize the total power of the GWs and another problem to minimize the maximum power of all the GWs were proposed and solved. For largesize buffers, a problem to maximize the weighted sum of the rates was proposed. Convexity was discussed for the three problems, and simulation results and discussion were provided.
Appendix 1
Consider the last constraint on the sum of rates in (4), i.e. . It is obvious that for minimum total power the constraint should hold with equality. Let , where is a constant. Assume the channel gains are ordered in a descending order, and assume is the minimum possible value for , such that,
(15) 
Now, subtract from and add to such that the equality still holds, where,
(16) 
where . Therefore, . The total power in this case is,
(17) 
since is a positive term. Hence, the minimum value for is used by decoding it last. By repeating the same procedure for up to , it is found that the decoding order should follow the decreasing order of the channel gain values.
Acknowledgment
This work is supported by the Center for Energy and GeoProcessing (CeGP) at King Fahd University of Petroleum and Minerals (KFUPM) and Georgia Institute of Technology, under grant number GTEC1601.
References
 [1] S. Savazzi, U. Spagnolini, L. Goratti, D. Molteni, M. LatvaAho, and M. Nicoli, “Ultrawide band sensor networks in oil and gas explorations,” IEEE Communications Magazine, vol. 51, no. 4, pp. 150–160, 2013.
 [2] D. Freed, “Cablefree nodes : The next generation land seismic system,” The Leading Edge, vol. 27, pp. 878–881, jul 2008.
 [3] R. Ellis, “Current cabled and cablefree seismic acquisition systems each have their own advantages and disadvantages  Is it possible to combine the two?,” First Break, vol. 32, no. 1, pp. 91–96, 2014.
 [4] V. A. Reddy, G. L. Stuber, and S. I. AlDharrab, “Energy Efficient Network Architecture for Seismic Data Acquisition via Wireless Geophones,” in 2018 IEEE International Conference on Communications (ICC), pp. 1–5, 2018.
 [5] V. J. Hodge, S. O. Keefe, M. Weeks, and A. Moulds, “Wireless Sensor Networks for Condition Monitoring in the Railway Industry : A Survey,” IEEE Transactions on Intelligent Transportation Systems, vol. 16, no. 3, pp. 1088–1106, 2015.
 [6] X. Yang, D. Deng, and M. Liu, “An overview of routing protocols on Wireless Sensor Network,” in 4th International Conference on Computer Science and Network Technology, pp. 1000–1003, 2015.
 [7] P. Park, S. C. Ergen, C. Fischione, C. Lu, and K. H. Johansson, “Wireless Network Design for Control Systems : A Survey,” IEEE Communications Surveys & Tutorials, vol. 20, no. 2, pp. 978–1013, 2018.
 [8] D. Sharma, S. Verma, and K. Sharma, “Network Topologies in Wireless Sensor Networks : A Review,” International Journal of Electronics & Communication Technology, vol. 4, no. 3, 2013.
 [9] D. G. Reina, S. L. Toral, F. Barrero, N. Bessis, and E. Asimakopoulou, “The Role of Ad Hoc Networks in the Internet of Things : A Case Scenario for The Role of Ad Hoc Networks in the Internet of Things : A Case Scenario for Smart,” in Internet of Things and Intercooperative Computational Technologies for Collective Intelligence, pp. 89–113, Springer, 2013.

[10]
N. Iqbal, S. AlDharrab, A. Muqaibel, W. Mesbah, and G. Stüber, “Analysis of Wireless Seismic Data Acquisition Networks using Markov Chain Models,” in
2018 IEEE 29th Annual International Symposium on Personal, Indoor and Mobile Radio Communications (PIMRC), (Bologna), pp. 1–5, 2018.  [11] A. Othman, W. Mesbah, N. Iqbal, S. AlDharrab, A. Muqaibel, and G. Stüber, Sumrate maximization for wireless seismic data acquisition systems, pp. 181–185. Society of Exploration Geophysicists, 2018.
 [12] A. El Gamal and Y.H. Kim., Network Information Theory. Cambridge University Press, first ed., 2011.
 [13] G. J. O. Vermeer, 3D Seismic Survey Design. Society of Exploration Geophysicists, second ed., 2012.
 [14] I. Boussaïd, J. Lepagnot, and P. Siarry, “A survey on optimization metaheuristics,” Information Sciences, vol. 237, pp. 82–117, 2013.
 [15] N. Iqbal, A. Zerguine, and N. AlDhahir, “Decision Feedback Equalization using Particle Swarm Optimization,” Signal Processing, vol. 108, pp. 1–12, 2015.
 [16] F. Marini and B. Walczak, “Particle swarm optimization (PSO). A tutorial,” Chemometrics and Intelligent Laboratory Systems, vol. 149, pp. 153–165, dec 2015.
 [17] A. Rezaee Jordehi and J. Jasni, “Particle swarm optimisation for discrete optimisation problems: a review,” Artificial Intelligence Review, vol. 43, no. 2, pp. 243–258, 2012.
 [18] T. Kaur and D. Kumar, “Particle Swarm OptimizationBased Unequal and Fault Tolerant Clustering Protocol for,” IEEE Sensors Journal, vol. 18, no. 11, pp. 4614–4622, 2018.

[19]
L. Yan, R. Zhu, Y. Liu, and N. Mo, “TrAdaBoost Based on Improved Particle Swarm Optimization for CrossDomain Scene Classification With Limited Samples,”
IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, vol. 11, no. 9, pp. 3235–3251, 2018.  [20] N. Iqbal, A. Zerguine, and N. AlDhahir, “Adaptive equalisation using particle swarm optimisation for uplink SCFDMA,” Electronic Letters, vol. 50, no. 6, pp. 469–471, 2014.

[21]
G. Yavuz and D. Aydin, “Angle Modulated Artificial Bee Colony Algorithms for Feature Selection,”
Applied Computational Intelligence and Soft Computing, vol. 2016, pp. 1–6, 2016.  [22] Z. L. Yang, A. Wu, and H. Q. Min, “Deployment of wireless sensor networks for oilfield monitoring by multiobjective discrete binary particle swarm optimization,” Journal of Sensors, 2016.
 [23] M. Dorigo and T. Stützle, Ant Colony Optimization: Overview and Recent Advances, pp. 311–351. Cham: Springer International Publishing, 2019.
 [24] M. Marzband, E. Yousefnejad, A. Sumper, and J. L. DomínguezGarcía, “Real time experimental implementation of optimum energy management system in standalone Microgrid by using multilayer ant colony optimization,” International Journal of Electrical Power & Energy Systems, vol. 75, pp. 265–274, feb 2016.
 [25] J. Zhang, H. Huang, T. Liu, and C. Zhang, “Multiscale Fractures Characterization Based on Ant Colony Optimization and TwoDimensional Variational Mode Decomposition,” IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, vol. 11, no. 7, pp. 2562–2570, 2018.
 [26] M. M. Mafarja and S. Mirjalili, “Hybrid Whale Optimization Algorithm with simulated annealing for feature selection,” Neurocomputing, vol. 260, pp. 302–312, oct 2017.

[27]
G. Pampara, N. Franken, and A. P. Engelbrecht, “Combining Particle Swarm
Optimisation with Angle Modulation to Solve Binary Problems,”
Congress on Evolutionary Computation
, pp. 89–96, 2005.  [28] K. Krynicki, M. E. Houle, and J. Jaen, “A Nonhybrid Ant Colony Optimization Heuristic for Convergence Quality,” Proceedings  2015 IEEE International Conference on Systems, Man, and Cybernetics, SMC 2015, pp. 1706–1713, 2016.
 [29] Q. Zhu and L. Wang, “The analysis of the convergence of ant colony optimization algorithm,” Frontiers of Electrical and Electronic Engineering in China, vol. 2, no. 3, pp. 268–272, 2007.
Comments
There are no comments yet.