Maximizing Computation Rate for Sustainable Wireless Powered MEC network: An Efficient Dynamic Task Offloading Algorithm with User Assistance

1. Introduction

In the era of the Internet of Things (IoT) [1,2,3,4], the surge in the number of mobile devices and network traffic has created a high demand for data processing capabilities and response speeds. Mobile Edge Computing (MEC) technology effectively enhances the operational efficiency of smart devices by providing powerful computing resources at the network’s edge, particularly in delay-sensitive applications such as autonomous driving and AR/VR [5]. By offloading complex computational tasks to MEC servers [6,7,8,9,10,11], resource-constrained mobile devices (MDs) can significantly alleviate the pressure experienced when running high-demand applications and achieve a notable leap in performance. However, the limited battery capacity of mobile devices has become a bottleneck that restricts their further development [12]. This limitation is especially evident for devices that require continuous operation over long periods and are not easily recharged, such as in remote areas or emergency situations, where insufficient battery life can severely impact the functionality and reliability of the devices. Therefore, despite the significant advantages of MEC in enhancing network performance, the limitation of battery endurance remains a key issue that urgently needs to be addressed.

In addition to battery constraints, reducing the energy consumption of IoT devices during the data offloading process is equally important. Trillions of tiny smart sensors make up the Internet of Things, facing significant limitations in computational capabilities and energy supplied by batteries. Advancements in wireless energy harvesting technologies, including renewable energy harvesting and wireless power transfer (WPT) [13], can alleviate the challenges previously posed by battery capacity limitations. Renewable energy sources like solar, wind, and ocean energy can provide power to some extent, but they are significantly influenced by natural conditions such as weather and climate [14]. To address this issue, green wireless charging technology has emerged. This technology can offers stable energy to devices through radio frequency signals and store it in the batteries of IoT nodes for future use, extending battery life [15,16]. To ensure that nodes do not fail due to energy depletion, green wireless charging adheres to the principle of energy neutrality operation [17] state, ensuring that the energy consumed in any operation will never exceed the energy collected. Green Wireless Powered Mobile Edge Computing (WPMEC) combines the strengths of WPT and MEC, enhancing devices’ computational capabilities and energy self-sufficiency. In the upcoming 6G networks, green WPMEC will provide IoT devices with quick response and real-time experiences [18], while reducing operational costs and extending the lifespan of devices.

However, WPMEC networks face the challenge of the dual far-near effect [19] caused by positional differences, which has prompted the development of edge collaborative networks [20,21,22] to optimize application offloading performance. By introducing a user cooperation (UC) mechanism, where nearby users assist distant users in accelerating the processing of computational tasks while effectively offloading their own tasks, this collaborative approach leverages the superior channel conditions of nearby users to gain more energy during the WPT phase. This not only addresses the unfairness caused by geographical location but also enhances the efficiency of energy utilization. The dense deployment of smart IoT devices further facilitates the opportunity to utilize the unused computing resources of idle devices and wireless energy harvesting. These devices, by assisting in completing the computational tasks of distant users, contribute to improving the overall computational performance of the WPMEC network.

To demonstrate the effectiveness of UC, recent studies, such as references [23,24], have effectively addressed the dual far-near effect in WPMEC networks through UC. The D2D communication in reference [25] and the incentive mechanism in reference [26] are both designed to facilitate resource sharing and collaborative offloading. In references [20,24,27,28], authors have focused on studying the basic three-node WPMEC model, which involves a Far User (FU) being allowed to offload computational input data to a Near User (NU). In reference [29], researchers designed a Non-Orthogonal Multiple Access (NOMA) based computation offloading scheme aimed at enhancing the performance of multi-user MEC systems. Google has also developed federated learning technology, which enables multiple devices to collaborate on machine learning tasks. Despite this, these studies are often based on the assumption of determinism or predictability of future information, failing to fully integrate the dynamic changes of the network environment, which may affect the efficiency and success rate of task offloading and processing.

This paper investigates the basic three-node green WPMEC network shown in Figure 1, focusing on the use of collaborative communication technology to accomplish the computation-intensive and delay-sensitive tasks powered by the HAP. Our goal is to maximize the network’s data processing capability in a real-time dynamic offloading system, taking into account the randomness of data generation and the high dynamics of wireless channel conditions. The challenges we face in addressing this issue include the unpredictability of task arrivals and the dynamics of channel conditions, as well as the coupling of variables in resource allocation, which makes traditional convex optimization methods inapplicable. To tackle these challenges, we have designed an efficient dynamic task offloading algorithm, the User-Assisted Dynamic Resource Allocation Algorithm (UADRA). This algorithm employs Lyapunov optimization techniques to transform the problem into a simplified form that relies on current information, and performs dynamic resource allocation in each time slot to enhance the network’s data processing capability. Our primary contributions are summarized as follows:

• We propose a model for maximizing the long-term weighted computation rate in a green sustainable WPT-MEC network, subject to energy consumption constraints. Our model takes into account the randomness of task arrivals and the dynamic variations in wireless channel states. By extending the models in [28,30], we have effectively mitigated the double near-far effect and fostered collaboration among users. This is achieved by considering the different data weights of near and far nodes, which enhances data transmission efficiency, system flexibility, and alignment with real-world applications.

• We design a low-complexity dynamic control algorithm called UADRA, based on Lyapunov optimization theory, while providing rigorous mathematical theory analysis to demonstrate its performance. A virtual queue is introduced to transforms the time-average energy constraint to a queue stable requirement. Leveraging the drift-plus-penalty technique, we decouple the original multi-stage stochastic problem into a non-convex deterministic sub-problem for each time slot. Furthermore, by using variable substitution and convex optimization theory, we convert these sub-problems into convex problems with a minimal number of decisions , that can be solved efficiently. Our algorithm can work efficiently without depending on the prior knowledge of the system information.

• We have conducted extensive simulations to evaluate the effectiveness and practicality of our proposed algorithm, particularly examining how the control parameter V, network bandwidth, energy constraint, task arrival rate, and geographical distance affect the average computation rate and network stability. Experimental results demonstrate that the proposed algorithm not only shows superior performance over benchmark methods but also achieves this improvement while simultaneously ensuring the stability of the system queues, with an enhancement of up to 4% in overall performance. Moreover, our algorithm effectively achieves a trade-off between computation rate and stability, adhering to the theoretical bounds of $\left[O\right(1/V),O(V\left)\right]$.

The remainder of this paper is organized as follows. Section II we presents the system model of the user-assistance green WPMEC network and sets forth a MSSO problem. In Section III, we utilize the Lyapunov optimization approach to tackle the problem, putting forward an effective dynamic offloading algorithm with an accompanying theoretical analysis of its performance. Section IV evaluate the efficacy of the suggested algorithm via simulation outcomes. Conclusively, Section V concludes the paper.

1.1. Related Work

The integration of WPT technology with MEC networks provides an effective solution for IoT devices, enhancing their energy and computing capabilities with controllable power supply and low-latency services. Recent research has extensively explored the potential of these wirelessly powered MEC networks. For instance, in [31], researchers optimized charging time and data offloading rates for WPT-MEC IoT sensor networks to improve computational rates in scenarios. Furthermore, authors in [32] investigated a NOMA-assisted WPT-MEC network with a nonlinear EH model, successfully enhancing the system’s Computational Energy Efficiency (CEE) by fine-tuning key parameters within the network. Specifically, to meet the energy consumption requirements of devices, the authors in [33] proposed a Particle Swarm Optimization (PSO)-based algorithm. The goal was to reduce the latency of devices processing computational data streams by jointly optimizing charging and offloading strategies. Additionally, in [34], the authors focused on the computational latency issue in WPT-MEC networks. They found suitable offloading ratio strategies to achieve synchronized latency for all WDs, effectively reducing the duration of the overall computational task.

To tackle the dual far-near effect issue, researchers have begun to focus on user-assisted WPMEC networks and have confirmed their effectiveness in enhancing the computing performance of distant users. Specifically, in [35], the study analyzed a three-node system composed of distant users, nearby users, and the base station within a user-assisted MEC-NOMA network model, addressing the optimization problem of joint transmission time and power allocation for users. Furthermore, [36,37] respectively explored joint computing and communication collaboration schemes and the application of Device-to-Device (D2D) communication in MEC. The method proposed in [36] aims to maximize the total computing rate of the network with the assistance of nearby users, while [37] focuses on minimizing the overall network response delay and energy consumption through joint multi-user collaborative partial offloading, transmission scheduling, and computing allocation. [38] extended this research, expanding from a single collaboration pair to multiple collaboration pairs, proposing a scheme to achieve the minimization of the total energy consumption of the AP.

In user assistance networks, the online collaborative offloading method, which is highly adaptable and can promptly respond to changes in task arrivals, has garnered significant attention from the research community. For instance, in [39], to address the randomness of energy and data arrivals, a Lyapunov optimization-based method was proposed to maximize the long-term system throughput. Furthermore, in [40], the authors studied and proposed a Lyapunov-based Profit Maximization (LBPM) task offloading algorithm in the context of the Internet of Vehicles (IoV), which aims at maximizing the time-averaged profit as the optimization goal. Additionally, in [41], within the application of MEC in the industrial IoT, a Lyapunov-based privacy-aware framework was introduced, which not only addressed privacy and security issues but also achieved optimization effects by reducing energy consumption. In [42], focusing on a multi-device single MEC system, the energy-saving task offloading problem was formulated as a time-averaged energy minimization problem considering queue length and resource constraints.

Different from prior studies, this paper is dedicated to addressing the challenges of dynamic task offloading in a green sustainable WPMEC networks with user assistance. We take into account the total system’s energy consumption constraint, the dynamically arriving tasks in real-time scenarios, and the high dynamics of wireless channel conditions. Moreover, the temporal coupling between WPT and user collaborative communication, as well as the coupling of data offloading timing and transmission power in collaborative communication, impose significant challenge to the problem.

2. System Model

2.1. Communication Model

As depicted in Figure 1, we consider a WPMEC system that comprises two MDs and one HAP. The HAP is equipped with both an RF energy transmitter and an MEC server, offering wireless energy and computation offloading services to MDs within the base station’s coverage area. Both MDs work on the same frequency band and are tasked with completing computation-intensive, delay-sensitive operations. One MD, known as the Far User (FU), is situated at a considerable distance from the HAP. In contrast, the other MD, called the Near User (NU), is in closer proximity to the HAP and serves as a facilitator, aiding the FU by offloading a portion of its tasks to the HAP.

To minimize interference, each MD utilizes a Time Division Duplex (TDD) approach to switch between communication and energy harvesting phases. We implement a Time Division Multiple Access (TDMA) framework with each time slot lasting T seconds. At the start of each time slot, both MDs capture RF signals emitted by the HAP for energy harvesting. Given the suboptimal channel conditions between the FU and the HAP, along with the compounded near-far effect, the FU is capable of transferring some of its computational data to the NU. The NU then forwards this data to the HAP. Additionally, the NU is scheduled to offload its own computational tasks to the HAP. The primary symbols and definitions used are enumerated in Table 1.

2.2. Task Queue Model

Both FU and NU maintain a queue to buffer the randomly arriving tasks for awaiting processing. Let $Q_f\left(t\right)$ and $Q_n\left(t\right)$ denote the queue lengths at the FU and NU at slot t, respectively, which cab be observed at the beginning of time slot t. The backlog of task queue update according to the following equations:

$$Q_f\left(t+1\right)=max\left[Q_f\left(t\right)-\left(d_f^{\mathrm{loc}}\left(t\right)+d_f^{\mathrm{off}}\left(t\right)\right),0\right]+A_f\left(t\right)$$

(1)

$$Q_n\left(t+1\right)=max\left[Q_n\left(t\right)-\left(d_n^{\mathrm{loc}}\left(t\right)+d_{n,a2}^{\mathrm{off}}\left(t\right)\right),0\right]+A_n\left(t\right)$$

(2)

where $A_f\left(t\right)$ and $A_n\left(t\right)$ represent the raw task data arriving at the FU and NU data queues during time slot t, respectively.

In reality, the data queue lengths $Q_n$ and $Q_f$ are always positive. Additionally, since the data processing amount by the MDs cannot exceed the current queue length, we have the following constraints:

$$d_f^{\mathrm{loc}}\left(t\right)+d_{f,n}^{\mathrm{off}}\left(t\right)\le Q_f\left(t\right)$$

(3)

$$d_n^{\mathrm{loc}}\left(t\right)+d_{n,a2}^{\mathrm{off}}\left(t\right)\le Q_n\left(t\right)$$

(4)

Additionally, if the local computation of the MD is capable of processing the entire current queue, then local computation is given priority , which is also in line with practical requirements.

2.3. Computation Model

(1) Local Computation Mode: MDs are capable of performing continuous local computations. Let $f_f$ and $f_n$ denote the local CPU frequencies of the FU and NU, respectively. Additionally, let ${\varphi }_f$ and ${\varphi }_n$ represent the CPU cycles required to process one bit of data at the FU and NU, respectively. Consequently, the amount of raw data (in bits) processed locally by MDs within a time slot t and the corresponding energy expenditure over that period can be determined by the following threshold rules.

For the FU, the local data processing in bits is given by

$$d_f^{\mathrm{loc}}\left(t\right)=\left\{\begin{array}{cc}{\displaystyle \frac{f_{f}T}{{\varphi }_f}},\hfill & \mathrm{if}{\displaystyle \frac{f_{f}T}{{\varphi }_f}}<Q_f\left(t\right)\hfill \\ Q_f\left(t\right),\hfill & \mathrm{others}\hfill \end{array}\right.$$

(5)

At the begin of each time slot t, we obtain the queue length of the FU. If the queue length exceeds its local processing capacity, the FU will utilize a user collaboration mechanism to offload partial tasks; if the queue length is within the computing capacity, the FU will perform complete local computation (that is, ${\tau }_1=0,{\tau }_2=0$), and the corresponding energy consumed during this period, as combined with $d\left(t\right)=ft/\varphi$ (where $d\left(t\right)$ denotes the amount of data processed locally at slot t) [30], is:

$$e_f^{\mathrm{loc}}\left(t\right)=\left\{\begin{array}{cc}\kappa {\left(f_f\right)}^{3}T,\hfill & \mathrm{if}{\displaystyle \frac{f_{f}T}{{\varphi }_f}}<Q_f\left(t\right)\hfill \\ \kappa {\left(f_f\right)}^2{Q}_f\left(t\right){\varphi }_f,\hfill & \mathrm{others}\hfill \end{array}\right.$$

(6)

Similarly, for the NU, the local data processing is

$$d_n^{\mathrm{loc}}\left(t\right)=\left\{\begin{array}{cc}{\displaystyle \frac{f_{n}T}{{\varphi }_n}},\hfill & \mathrm{if}{\displaystyle \frac{f_{n}T}{{\varphi }_n}}<Q_n\left(t\right)\hfill \\ Q_n\left(t\right),\hfill & \mathrm{others}\hfill \end{array}\right.$$

(7)

When the NU engages in complete local computation, ${\tau }_3=0$, with the energy consumed being

$$e_n^{\mathrm{loc}}\left(t\right)=\left\{\begin{array}{cc}\kappa {\left(f_n\right)}^{3}T,\hfill & \mathrm{if}{\displaystyle \frac{f_{n}T}{{\varphi }_n}}<Q_n\left(t\right)\hfill \\ \kappa {\left(f_n\right)}^2{Q}_n\left(t\right){\varphi }_n,\hfill & \mathrm{others}\hfill \end{array}\right.$$

(8)

Here, $\kappa >0$ represents the computational energy efficiency coefficient of the processing chip [43].

(2) Offloading Mode: We adopt a partial task offloading strategy for both FU and NU. As illustrated in Figure 2. Initially, the FU transfers a portion of its task data to the NU during time ${\tau }_1^t$, using a transmission power $P_1^t$. Subsequently, the NU forwards the FU’s data to the HAP during time ${\tau }_2^t$ with a transmission power $P_2^t$. Furthermore, the NU offloads its own data to the HAP during time ${\tau }_3^t$ using transmission power $P_3^t$. Throughout this process, the NU consume the energy it has harvested. It is essential that the amount of data offloaded by the FU for remote computation at the HAP corresponds to the data initially offloaded to the NU.

In practical applications, the HAP possesses significantly greater computational capacity and transmission power compared to the MDs. As a result, the time required for data computation and feedback at the HAP is negligible. Therefore, the total time allocated for energy harvesting and task offloading by MDs must not exceed the time slot duration T, as represented by the following inequality:

$${\tau }_0^t+{\tau }_1^t+{\tau }_2^t+{\tau }_3^t\le T$$

(9)

We focus solely on the energy consumption of the MDs during the offloading process. The transmission power of the FU $P_1^t$ is constrained by $P_1^t\le P_f^{max}$. According to Shannon’s theorem, the amount of data offloaded from the FU to the NU is given by

$$d_{f,n}^{\mathrm{off}}\left(t\right)={\tau }_1^{t}W{log}_2\left(1+{\displaystyle \frac{P_1^t{g}_f^t}{{\sigma }^2}}\right)$$

(10)

The corresponding energy consumption for this process is

$$e_{f,n}^{\mathrm{off}}\left(t\right)=P_1^t{\tau }_1^t$$

(11)

where W is the channel bandwidth, $g_f^t$ is the channel gain from the FU to the NU during time slot t, and ${\sigma }^2$ represents the power of additive white Gaussian noise.

After the initial time $\left({\tau }_0^t+{\tau }_1^t\right)$, the NU has received the offloaded tasks from the FU. The NU must then determine the transmission powers $P_2^t$ for relaying the FU’s data to the HAP and $P_3^t$ for offloading its own data. The quantity of data relayed by the NU from the FU to the HAP is expressed as

$$d_{n,a1}^{\mathrm{off}}\left(t\right)={\tau }_2^{t}W{log}_2\left(1+{\displaystyle \frac{P_2^t{g}_n^t}{{\sigma }^2}}\right)$$

(12)

with the associated energy consumption being

$$e_{n,a1}^{\mathrm{off}}\left(t\right)=P_2^t{\tau }_2^t$$

(13)

Similarly, the amount of data offloaded by the NU to the HAP is

$$d_{n,a2}^{\mathrm{off}}\left(t\right)={\tau }_3^{t}W{log}_2\left(1+{\displaystyle \frac{P_3^t{g}_n^t}{{\sigma }^2}}\right)$$

(14)

and the corresponding energy consumption is

$$e_{n,a2}^{\mathrm{off}}\left(t\right)=P_3^t{\tau }_3^t$$

(15)

Here, $P_2^t,P_3^t\le P_n^{max}$ denote the transmission power of the NU, and $g_n^t$ represents signifies the channel gain from the NU to the HAP during the respective time slot.

During each time slot t, the NU is required to process all offloading tasks received from the FU, adhering to the constraint that the amount of data offloaded by the FU does not exceed that of the NU’s capability to relay it

$$d_{f,n}^{\mathrm{off}}\left(t\right)\le d_{n,a1}^{\mathrm{off}}\left(t\right)$$

(16)

2.4. Energy Consumption Model

As illustrated in Figure 2, during the initial phase ${\tau }_0^t$, the HAP disseminates radio frequency energy to the MDs. The energy harvested by each MD can be represented as

$$e_f^{\mathrm{eh}}\left(t\right)={\alpha }_1{\tau }_0^t{P}_0^t$$

(17)

$$e_n^{\mathrm{eh}}\left(t\right)={\alpha }_2{\tau }_0^t{P}_0^t$$

(18)

where $e_f^{\mathrm{eh}}\left(t\right)$ and $e_n^{\mathrm{eh}}\left(t\right)$ denote the energy harvested by the FU and NU from the HAP during time slot t, respectively. The coefficients ${\alpha }_1$ and ${\alpha }_2$ are defined as ${\alpha }_1=\mu h_f^t$, ${\alpha }_2=\mu h_n^t$, where $0<\mu <1$, represents the energy conversion efficiency. $P_0^t$ signifies the RF energy transmission power of the HAP, and the channel power gains from the HAP to the FU and NU are represented by $h_f^t$ and $h_n^t$, respectively. Assuming block fading, these gains remain consistent within a time slot but fluctuate independently between frames. It is important to note that the WPT from the HAP is the sole energy source for performing computational tasks.

In each time slot t, the energy consumption of the FU and NU must adhere to the following constraints:

$$e_f^{\mathrm{loc}}\left(t\right)+e_{f,n}^{\mathrm{off}}\left(t\right)\le {\alpha }_1{P}_0^t{\tau }_0^t$$

(19)

$$e_n^{\mathrm{loc}}\left(t\right)+e_{n,a1}^{\mathrm{off}}\left(t\right)+e_{n,a2}^{\mathrm{off}}\left(t\right)\le {\alpha }_2{P}_0^t{\tau }_0^t$$

(20)

For a sustainable wireless charging IoT network, all energy for the wireless nodes originates from the RF signals emitted remotely by the HAP. Therefore, controlling the energy consumption of the HAP is of significant importance to the development of the entire IoT network. Consequently, we assume that the time-average energy consumption of the HAP is subject to the following constraints.

$$\overline{P_0^t{\tau }_0^t}\le \gamma$$

(21)

where $\overline{x}\triangleq \underset{K\to \infty }{\lim \sup }{\displaystyle \frac{1}{K}}{\displaystyle \sum _{t=0}^{K-1}}\mathbb{E}\left\{x\right\}$. $P_0^t{\tau }_0^t$ is the energy consumption for HAP at the tth time slot, and $\gamma$ is the energy threshold.

2.5. Problem Formulation

For a dynamically changing WPT-MEC network system, maintaining system stability is crucial due to the stochastic arrival of tasks and the dynamic changes in the channel environment [44] . Therefore, we first provide the definition of system queue stability as follows.

Definition 1. (Queue Stability): The task data queue is strong stable if it satisfies [44]

$$\underset{K\to \infty }{\lim }{\displaystyle \frac{1}{K}}\stackrel{K}{{\displaystyle \sum }_{t=1}}\mathbb{E}\left\{Q_f\left(t\right)\right\}<\infty$$

,

$$\underset{K\to \infty }{\lim }{\displaystyle \frac{1}{K}}\stackrel{K}{{\displaystyle \sum }_{t=1}}\mathbb{E}\left\{Q_n\left(t\right)\right\}<\infty$$

In this paper, we aim to design a dynamic offloading algorithm to maximize the long-term average weighted sum computation rate for all MDs, by making decision of time allocation ${\tau }^t=\left[{\tau }_0^t,{\tau }_1^t,{\tau }_2^t,{\tau }_3^t\right]$, and the power allocation ${\mathit{P}}^t=\left[P_0^t,P_1^t,P_2^t,P_3^t\right]$ at each time slot t under the long term power consumption constraint of IoT network. Our proposed algorithm should ensure the stability of data queues and work without prior knowledge of random channel conditions and data arrival patterns, optimizing within each time slot. By denoting $\tau ={\left\{{\tau }^t\right\}}_{t=1}^K$, $\mathit{P}={\left\{{\mathit{P}}^t\right\}}_{t=1}^K$, the maximization of weighted computation rate problem of WPT-MEC can be formulated as a MSSO problem (P1) as:

$$\begin{array}{cc}\hfill P1:\underset{\tau ,\mathit{P}}{max}& \underset{K\to \infty }{\lim }{\displaystyle \frac{1}{K}}\stackrel{K}{{\displaystyle \sum }_{t=1}}{D}_{tot}\left(t\right)\hfill \end{array}$$

(22a)

$$\begin{array}{cc}\hfill s.t.& {\tau }_0^t+{\tau }_1^t+{\tau }_2^t+{\tau }_3^t\le \tau ,\hfill \end{array}$$

(22b)

$$\begin{array}{cc}\hfill & e_f^{\mathrm{loc}}\left(t\right)+e_{f,n}^{\mathrm{off}}\left(t\right)\le {\alpha }_1{P}_0^t{\tau }_0^t,\hfill \end{array}$$

(22c)

$$\begin{array}{cc}\hfill & e_n^{\mathrm{loc}}\left(t\right)+e_{n,a1}^{\mathrm{off}}\left(t\right)+e_{n,a2}^{\mathrm{off}}\left(t\right)\le {\alpha }_2{P}_0^t{\tau }_0^t,\hfill \end{array}$$

(22d)

$$\begin{array}{cc}\hfill & \overline{P_0^t{\tau }_0^t}\le \gamma ,\hfill \end{array}$$

(22e)

$$\begin{array}{cc}\hfill & \underset{K\to \infty }{\lim }{\displaystyle \frac{1}{K}}\stackrel{K}{{\displaystyle \sum }_{t=1}}\mathbb{E}\left\{Q_f\left(t\right)\right\}<\infty ,\underset{K\to \infty }{\lim }{\displaystyle \frac{1}{K}}\stackrel{K}{{\displaystyle \sum }_{t=1}}\mathbb{E}\left\{Q_n\left(t\right)\right\}<\infty ,\hfill \end{array}$$

(22f)

$$\begin{array}{cc}\hfill & d_f^{\mathrm{loc}}\left(t\right)+d_{f,n}^{\mathrm{off}}\left(t\right)\le Q_f\left(t\right),\hfill \end{array}$$

(22g)

$$\begin{array}{cc}\hfill & d_n^{\mathrm{loc}}\left(t\right)+d_{n,a2}^{\mathrm{off}}\left(t\right)\le Q_n\left(t\right),\hfill \end{array}$$

(22h)

$$\begin{array}{cc}\hfill & d_{f,n}^{\mathrm{off}}\left(t\right)\le d_{n,a1}^{\mathrm{off}}\left(t\right),\hfill \end{array}$$

(22i)

$$\begin{array}{cc}\hfill & {\tau }_0^t,{\tau }_1^t,{\tau }_2^t,{\tau }_3^t,\ge 0,\hfill \end{array}$$

(22j)

$$\begin{array}{cc}\hfill & 0\le p_0^t\le P_a^{\max },0\le p_1^t\le P_f^{\max },0\le p_2^t,p_3^t\le P_n^{\max },\hfill \end{array}$$

(22k)

where $D_{tot}\left(t\right)={\omega }_1\left(d_f^{\mathrm{loc}}\left(t\right)+d_{f,n}^{\mathrm{off}}\left(t\right)\right)+{\omega }_2\left(d_n^{\mathrm{loc}}\left(t\right)+d_{n,a2}^{\mathrm{off}}\left(t\right)\right)$, ${\omega }_1$ and ${\omega }_2$ denotes the fixed weight of the FU and NU, respectively. Constraint (22b) ensures the total time allocation does not exceed the time slot. Constraints (22c) and (22d) represent the energy constraints for the FU and Near NU, respectively. Constraint (22e) captures the average power constraints for the system. Constraints (22f) guarantee the stability of data queues. Constraints (22g) and (22h) denote the maximum data processing limits for the FU and NU within time slot t. Constraint (22i) ensures that the data offloaded from the FU can be processed within the same time slot. Problem (P1) poses a complex stochastic optimization challenge that conventional convex optimization methods cannot easily solve. To address this, we introduce a Lyapunov-based optimization algorithm that transforms the stochastic problem into a deterministic one within each time slot.

3. Algorithm Design

3.1. Lyapunov Optimization

To address the average power constraints, we introduce a virtual energy queue $\mathrm{Y}\left(t+1\right)={\left[\mathrm{Y}\left(\mathrm{t}\right)-\gamma +P_0^t{\tau }_0^t\right]}^{+}$, where ${\left[x\right]}^{+}\triangleq max\left\{0,x\right\}$[44]. Here $\mathrm{Y}\left(\mathrm{t}\right)$ can be seen as a queue with random "energy arrivals" $P_0^t{\tau }_0^t$ and a fixed "service rate" $\gamma$. Thus, we derive the following Lemma 1.

Lemma 1.

The long-term average power constraints will be met if the virtual queue $\mathrm{Y}\left(\mathrm{t}\right)$ satisfies average rate stability

$$\underset{K\to \infty }{\lim }{\displaystyle \frac{1}{K}}\mathbb{E}\left\{\mathrm{Y}\left(\mathrm{t}\right)\right\}=0$$

Proof. 

According to the above formula $\mathrm{Y}\left(\mathrm{t}+1\right)={[\mathrm{Y}\left(\mathrm{t}\right)-\gamma +P_0^t{\tau }_0^t]}^{+}$ , we can conclude that

$$\mathrm{Y}\left(\mathrm{t}+1\right)\ge \mathrm{Y}\left(\mathrm{t}\right)-\gamma +P_0^t{\tau }_0^t$$

(23)

we can expand all the terms of $t\in \{0,1,2,...,K-1\}$, sum them up, and then take the average, yielding

$${\displaystyle \frac{\mathrm{Y}\left(K\right)-\mathrm{Y}\left(0\right)}{K}}\ge {\displaystyle \frac{1}{K}}\stackrel{K-1}{{\displaystyle \sum }_{t=0}}{P}_0^t{\tau }_0^t-\gamma$$

(24)

Next, we simultaneously take the expectation on both sides and set K, obtaining:

$$\underset{K\to \infty }{\lim \sup }{\displaystyle \frac{\mathbb{E}\left\{\mathrm{Y}\left(\mathrm{t}\right)\right\}}{K}}\ge \underset{K\to \infty }{\lim \sup }{\displaystyle \frac{1}{K}}\stackrel{K-1}{{\displaystyle \sum }_{t=0}}\mathbb{E}\left\{P_0^t{\tau }_0^t\right\}-\gamma$$

(25)

By our assumption$\underset{K\to \infty }{\lim }{\displaystyle \frac{1}{K}}\mathbb{E}\left\{\mathrm{Y}\left(\mathrm{t}\right)\right\}=0$, it follows that $\underset{K\to \infty }{\lim \sup }{\displaystyle \frac{1}{K}}{\displaystyle {\displaystyle \sum }_{t=0}^{K-1}}\mathbb{E}\left\{P_0^t{\tau }_0^t\right\}-\gamma \le 0$, i.e., $\overline{P_0^t{\tau }_0^t}\le \gamma$ , the lemma is proven.    □

By defining a network queue vector $\mathsf{\Theta }\left(t\right)=\left\{Q_n\left(t\right),Q_f\left(t\right),\mathrm{Y}\left(\mathrm{t}\right)\right\}$, which encompasses the queue lengths for the NU, FU, and the virtual queue, respectively, we can obtain the associated Lyapunov function $L\left(\mathsf{\Theta }\left(t\right)\right)$ and the Lyapunov drift $\Delta \left(\mathsf{\Theta }\left(t\right)\right)$ as

$$L\left(\mathsf{\Theta }\left(t\right)\right)\triangleq {\displaystyle \frac{1}{2}}\left(Q_n{\left(t\right)}^2+Q_f{\left(t\right)}^2+\mathrm{Y}{\left(\mathrm{t}\right)}^2\right)$$

(26)

$$\Delta \left(\mathsf{\Theta }\left(t\right)\right)\triangleq \mathbb{E}\left\{\mathit{L}\left(\mathit{t}+\mathit{1}\right)-\mathit{L}\left(\mathit{t}\right)\mid \mathsf{\Theta }\left(\mathit{t}\right)\right\}$$

(27)

Employing the Lyapunov optimization theory, we drive the drift-plus-penalty as

$${\Delta }_V\left(\mathsf{\Theta }\left(t\right)\right)=\Delta \left(\mathsf{\Theta }\left(t\right)\right)-V\mathbb{E}\left\{D_{\mathrm{tot}}\left(t\right)\mid \mathsf{\Theta }\left(t\right)\right\}$$

(28)

Here,$V>0$ signifies the penalty’s importance weight. The Lyapunov optimization method aims to minimize the upper bound of the drift plus penalty, thereby maximizing the objective function while ensuring queue stability. Optimizing the objective function across each time slot leads to long-term optimality. It’s important to note that a higher value of V prioritizes objective maximization, whereas a lower value emphasizes queue stability. To obtain the upper bound of ${\Delta }_V\left(\mathsf{\Theta }\left(t\right)\right)$, we derive the following Lemma 2.

Lemma 2.

At each time slot t, for any control strategy $\left\{{\tau }^t,{\mathit{P}}^t\right\}$, the one-slot Lyapunov drift-plus-penalty${\Delta }_V\left(\mathsf{\Theta }\left(t\right)\right)$is bounded as per the following inequality

$$\begin{array}{cc}\hfill {\Delta }_V\left(\mathsf{\Theta }\left(t\right)\right)& \le B-V\mathbb{E}\left\{D_{\mathrm{tot}}\left(t\right)\mid \mathsf{\Theta }\left(t\right)\right\}+Q_f\left(t\right)\mathbb{E}\left\{A_f\left(t\right)-d_f^{\mathrm{loc}}\left(t\right)-d_{f,n}^{\mathrm{off}}\left(t\right)\mid \mathsf{\Theta }\left(t\right)\right\}\hfill \\ \hfill & +Q_n\left(t\right)\mathbb{E}\left\{A_n\left(t\right)-d_n^{\mathrm{loc}}\left(t\right)-d_{n,a2}^{\mathrm{off}}\left(t\right)\mid \mathsf{\Theta }\left(t\right)\right\}+\mathrm{Y}\left(\mathrm{t}\right)\mathbb{E}\left\{P_0^t{\tau }_0^t-\gamma \mid \mathsf{\Theta }\left(t\right)\right\}\hfill \end{array}$$

(29)

where B is a constant that satisfies the following $\forall t$

$$\begin{array}{cc}\hfill B& \ge {\displaystyle \frac{1}{2}}\mathbb{E}\left\{{\left(d_f^{\mathrm{loc}}\left(t\right)+d_{f,n}^{\mathrm{off}}\left(t\right)\right)}^2+A_f{\left(t\right)}^2\mid \mathsf{\Theta }\left(t\right)\right\}\hfill \\ \hfill & +{\displaystyle \frac{1}{2}}\mathbb{E}\left\{{\left(d_n^{\mathrm{loc}}\left(t\right)+d_{n,a2}^{\mathrm{off}}\left(t\right)\right)}^2+A_n{\left(t\right)}^2\mid \mathsf{\Theta }\left(t\right)\right\}\hfill \\ \hfill & +{\displaystyle \frac{1}{2}}\mathbb{E}\left\{{\gamma }^2+{\left(P_0^t{\tau }_0^t\right)}^2\mid \mathsf{\Theta }\left(t\right)\right\}\hfill \end{array}$$

(30)

Proof. 

For all $\forall a,b,c\ge 0$, we have the inequality ${\left(max[a-b,0]+c\right)}^2\le a^2+b^2+c^2+2a\left(c-b\right)$. By using the inequality, we have

$$\begin{array}{cc}\hfill \Delta Q_f\left(t\right)& ={\displaystyle \frac{1}{2}}\left(Q_f{\left(t+1\right)}^2-Q_f{\left(t\right)}^2\right)\hfill \end{array}$$

$$\begin{array}{cc}\hfill & \le {\displaystyle \frac{{\left(d_f^{\mathrm{loc}}\left(t\right)+d_{f,n}^{\mathrm{off}}\left(t\right)\right)}^2+A_f{\left(t\right)}^2}{2}}+Q_f\left(t\right)\left(A_f\left(t\right)-d_f^{\mathrm{loc}}\left(t\right)-d_{f,n}^{\mathrm{off}}\left(t\right)\right)\hfill \\ \hfill \Delta Q_n\left(t\right)& ={\displaystyle \frac{1}{2}}\left(Q_n{\left(t+1\right)}^2-Q_n{\left(t\right)}^2\right)\hfill \end{array}$$

(31)

$$\begin{array}{cc}\hfill & \le {\displaystyle \frac{{\left(d_n^{\mathrm{loc}}\left(t\right)+d_{n,a2}^{\mathrm{off}}\left(t\right)\right)}^2+A_n{\left(t\right)}^2}{2}}+Q_n\left(t\right)\left(A_n\left(t\right)-d_n^{\mathrm{loc}}\left(t\right)-d_{n,a2}^{\mathrm{off}}\left(t\right)\right)\hfill \\ \hfill \Delta \mathrm{Y}\left(\mathrm{t}\right)& ={\displaystyle \frac{1}{2}}\left(\mathrm{Y}{\left(\mathrm{t}+1\right)}^2-\mathrm{Y}{\left(\mathrm{t}\right)}^2\right)\hfill \end{array}$$

(32)

$$\begin{array}{cc}\hfill & \le {\displaystyle \frac{{\left(P_0^t{\tau }_0^t\right)}^2+{\gamma }^2}{2}}+Y\left(t\right)\left(P_0^t{\tau }_0^t-\gamma \right)\hfill \end{array}$$

(33)

Upon combining inequalities (31), (32), and (33), the resulting expression yields the upper bound of the Lyapunov drift-plus-penalty.    □

By applying the drift-plus-penalty minimization technique, and eliminating the constant term observable at the start of time slot t in (29), we can obtain the optimal solution to problem (P1) by solving the following problem in each individual time slot.

$$\begin{array}{cc}\hfill P2:\underset{{\tau }^t,{\mathit{P}}^t}{min}& \mathrm{Y}\left(\mathrm{t}\right)P_0^t{\tau }_0^t-\left[Q_f\left(t\right)+{\omega }_{1}V\right]\left\{{\tau }_1^{t}W{log}_2\left(1+{\displaystyle \frac{P_1^t{g}_f^t}{{\sigma }^2}}\right)+d_f^{\mathrm{loc}}\left(t\right)\right\}\hfill \end{array}$$

$$\begin{array}{cc}\hfill & -\left[Q_n\left(t\right)+{\omega }_{2}V\right]\left\{{\tau }_3^{t}W{log}_2\left(1+{\displaystyle \frac{P_3^t{g}_n^t}{{\sigma }^2}}\right)+d_n^{\mathrm{loc}}\left(t\right)\right\}\hfill \end{array}$$

(34a)

$$\begin{array}{cc}\hfill s.t.& \left(22b\right)-\left(22d\right),\left(22g\right)-\left(22k\right)\hfill \end{array}$$

(34b)

Due to the the non-convex constraints (34b), (P2) remains a non-convex problem. We introduce auxiliary variables ${\phi }_0=p_0^t{\tau }_0^t,{\phi }_1=p_1^t{\tau }_1^t,{\phi }_2=p_2^t{\tau }_2^t,{\phi }_3=p_3^t{\tau }_3^t$ here. According to equation (22j), we have $0\le {\phi }_0\le {\tau }_0{P}_a^{\max },0\le {\phi }_1\le {\tau }_1{P}_f^{\max },0\le {\phi }_2\le {\tau }_2{P}_n^{\max },0\le {\phi }_3\le {\tau }_3{P}_n^{\max }$. We denote $\phi =\left[{\phi }_0,{\phi }_1,{\phi }_2,{\phi }_3\right]$. To simplify the mathematical expressions, we omit t here. So (P2) can be equivalently substituted as

$$\begin{array}{cc}\hfill P3:\underset{\tau ,\phi }{max}& \mathrm{Y}{\phi }_0+C_1\left\{{\tau }_{1}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{1}A}{{\tau }_1}}\right)+d_f^{\mathrm{loc}}\right\}+C_2\left\{{\tau }_{3}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{3}B}{{\tau }_3}}\right)+d_n^{\mathrm{loc}}\right\}\hfill \end{array}$$

(35a)

$$\begin{array}{cc}\hfill s.t.& {\tau }_0+{\tau }_1+{\tau }_2+{\tau }_3\le \tau ,\hfill \end{array}$$

(35b)

$$\begin{array}{cc}\hfill & {\phi }_1+k{\left(f_f\right)}^{3}T\le {\alpha }_1{\phi }_0,\hfill \end{array}$$

(35c)

$$\begin{array}{cc}\hfill & {\phi }_2+{\phi }_3+k{\left(f_n\right)}^{3}T\le {\alpha }_2{\phi }_0,\hfill \end{array}$$

(35d)

$$\begin{array}{cc}\hfill & {\tau }_{1}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{1}A}{{\tau }_1}}\right)\le {\tau }_{2}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{2}B}{{\tau }_2}}\right),\hfill \end{array}$$

(35e)

$$\begin{array}{cc}\hfill & {\displaystyle \frac{f_{f}T}{{\varphi }_f}}+{\tau }_{1}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{1}A}{{\tau }_1}}\right)\le Q_f,\hfill \end{array}$$

(35f)

$$\begin{array}{cc}\hfill & {\displaystyle \frac{f_{n}T}{{\varphi }_n}}+{\tau }_{3}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{3}B}{{\tau }_3}}\right)\le Q_n,,\hfill \\ \hfill & 0\le {\phi }_0\le {\tau }_0{P}_a^{\max },0\le {\phi }_1\le {\tau }_1{P}_f^{\max },\hfill \end{array}$$

(35g)

$$\begin{array}{cc}\hfill & 0\le {\phi }_2\le {\tau }_2{P}_n^{\max },0\le {\phi }_3\le {\tau }_3{P}_n^{\max },\hfill \end{array}$$

(35h)

$$\left(22j\right)$$

(35i)

The coefficients are defined as $C_1=-\left[Q_f+{\omega }_{1}V\right]$ and $C_2=-\left[Q_n+{\omega }_{2}V\right]$, where $A={\displaystyle \frac{g_f}{{\sigma }^2}}$, and $B={\displaystyle \frac{g_n}{{\sigma }^2}}$. Owing to the non-convex constraint (35e), problem (P3) remains non-convex. To address this, we introduce auxiliary variables ${\psi }_1$ and ${\psi }_2$, defined such that ${\psi }_1\le {\tau }_{1}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{1}A}{{\tau }_1}}\right),{\psi }_2\le {\tau }_{3}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{3}B}{{\tau }_3}}\right)$. By making these definitions, problem (P3) can be reformulated in terms of these auxiliary variables as follows:

$$\begin{array}{cc}\hfill P4:\underset{\tau ,\phi }{max}& \mathrm{Y}{\phi }_0+C_1\left\{{\tau }_{1}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{1}A}{{\tau }_1}}\right)+d_f^{\mathrm{loc}}\right\}+C_2\left\{{\tau }_{3}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{3}B}{{\tau }_3}}\right)+d_n^{\mathrm{loc}}\right\}\hfill \end{array}$$

(36a)

$$\begin{array}{cc}\hfill s.t.& {\tau }_0+{\tau }_1+{\tau }_2+{\tau }_3\le \tau ,\hfill \end{array}$$

(36b)

$$\begin{array}{cc}\hfill & {\phi }_1+k{\left(f_f\right)}^{3}T\le {\alpha }_1{\phi }_0,\hfill \end{array}$$

(36c)

$$\begin{array}{cc}\hfill & {\phi }_2+{\phi }_3+k{\left(f_n\right)}^{3}T\le {\alpha }_2{\phi }_0,\hfill \end{array}$$

(36d)

$$\begin{array}{cc}\hfill & {\psi }_1\le {\tau }_{2}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{2}B}{{\tau }_2}}\right),\hfill \end{array}$$

(36e)

$$\begin{array}{cc}\hfill & {\psi }_1\le {\tau }_{1}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{1}A}{{\tau }_1}}\right),\hfill \end{array}$$

(36f)

$$\begin{array}{cc}\hfill & {\psi }_1\le Q_f-{\displaystyle \frac{f_{f}T}{{\varphi }_f}},\hfill \end{array}$$

(36g)

$$\begin{array}{cc}\hfill & {\psi }_2\le {\tau }_{3}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{3}B}{{\tau }_3}}\right),\hfill \end{array}$$

(36h)

$$\begin{array}{cc}\hfill & {\psi }_2\le Q_n-{\displaystyle \frac{f_{n}T}{{\varphi }_n}},,\hfill \\ \hfill & 0\le {\phi }_0\le {\tau }_0{P}_a^{\max },0\le {\phi }_1\le {\tau }_1{P}_f^{\max },\hfill \end{array}$$

(36i)

$$\begin{array}{cc}\hfill & 0\le {\phi }_2\le {\tau }_2{P}_n^{\max },0\le {\phi }_3\le {\tau }_3{P}_n^{\max },\hfill \end{array}$$

(36j)

$$\left(22j\right),$$

(36k)

${\psi }_1={\tau }_{1}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{1}A}{{\tau }_1}}\right),{\psi }_2={\tau }_{3}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{3}B}{{\tau }_3}}\right)$ , when the problem (P4) reaches the optimal solution, it aligns with the original problem (P3). Specifically, for constraint (36e), ${\tau }_{2}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{2}B}{{\tau }_2}}\right)$ is a concave function of ${\tau }_2$ since the perspective operation applied to $W{log}_2\left(1+{\phi }_{2}B\right)$ preserves convexity [45] . Thus, constraint (36e) is convex, making problem (P4) a convex optimization problem. To further analyze the solution, we will employ convex optimization tools, such as CVX [46].

We introduce an efficient dynamic task offloading algorithm with user assistance to tackle Problem 4. Additionally, we apply the Lagrange method to gain meaningful insights into the optimal solution’s properties.

Theorem 1.

Given non-negative Lagrange multipliers ${\lambda }_i$, $i=1,2,\dots ,12$, the optimal power allocation $\mathit{P}=\left[P_0,P_1,P_2,P_3\right]$ must fulfill certain conditions

$$\begin{array}{cc}\hfill & P_0=\left\{\begin{array}{cc}0,\hfill & \mathit{if}{\tau }_0=0\hfill \\ {\left[P_a^{\max }\right]}^{+},\hfill & \mathit{others}\hfill \end{array}\right.\hfill \end{array}$$

(37a)

$$\begin{array}{cc}\hfill & P_1=\left\{\begin{array}{cc}0,\hfill & \mathit{if}{\tau }_1=0\hfill \\ {\left[{\displaystyle \frac{W\left(Q_f+V+{\lambda }_5\right)}{ln2\left({\lambda }_2+{\lambda }_{10}\right)}}-{\displaystyle \frac{{\sigma }^2}{g_f}}\right]}^{+},\hfill & \mathit{others}\hfill \end{array}\right.\hfill \end{array}$$

(37b)

$$\begin{array}{cc}\hfill & P_2=\left\{\begin{array}{cc}0,\hfill & \mathit{if}{\tau }_2=0\hfill \\ {\left[{\displaystyle \frac{{\lambda }_{4}W}{ln2\left({\lambda }_3+{\lambda }_{11}\right)}}-{\displaystyle \frac{{\sigma }^2}{g_n}}\right]}^{+},\hfill & \mathit{others}\hfill \end{array}\right.\hfill \end{array}$$

(37c)

$$\begin{array}{cc}\hfill & P_3=\left\{\begin{array}{cc}0,\hfill & \mathit{if}{\tau }_3=0\hfill \\ {\left[{\displaystyle \frac{W\left(Q_n+V+{\lambda }_7\right)}{ln2\left({\lambda }_3+{\lambda }_{12}\right)}}-{\displaystyle \frac{{\sigma }^2}{g_n}}\right]}^{+},\hfill & \mathit{others}\hfill \end{array}\right.\hfill \end{array}$$

(37d)

Proof. 

Let ${\lambda }_i\ge 0\left(i=1,2,\dots ,12\right)$ for denote the Lagrange multipliers corresponding to the constraints. The Lagrangian function for problem (P4), constructed based on these multipliers, is as follow

$$\begin{array}{cc}\hfill L\left(\tau ,\phi ,\lambda \right)& =\mathrm{Y}\left(\mathrm{t}\right){\phi }_0+C_1\left\{{\tau }_{1}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{1}A}{{\tau }_1}}\right)+{\displaystyle \frac{f_{f}T}{{\varphi }_f}}\right\}+C_2\left\{{\tau }_{3}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{3}B}{{\tau }_3}}\right)+{\displaystyle \frac{f_{n}T}{{\varphi }_n}}\right\}\hfill \\ \hfill & +{\lambda }_1\left[{\tau }_0^t+{\tau }_1^t+{\tau }_2^t+{\tau }_3^t-\tau \right]\hfill \\ \hfill & +{\lambda }_2\left[{\phi }_1+k{\left(f_f\right)}^{3}T-{\alpha }_1{\phi }_0\right]\hfill \\ \hfill & +{\lambda }_3\left[{\phi }_2+{\phi }_3+k{\left(f_n\right)}^{3}T-{\alpha }_2{\phi }_0\right]\hfill \\ \hfill & +{\lambda }_4\left[{\psi }_1-{\tau }_{2}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{2}B}{{\tau }_2}}\right)\right]\hfill \\ \hfill & +{\lambda }_5\left[{\psi }_1-{\tau }_{1}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{1}A}{{\tau }_1}}\right)\right]\hfill \\ \hfill & +{\lambda }_6\left[{\psi }_1-Q_f-{\displaystyle \frac{f_{f}T}{{\varphi }_f}}\right]\hfill \\ \hfill & +{\lambda }_7\left[{\psi }_2-{\tau }_{3}W{log}_2\left(1+{\displaystyle \frac{{\phi }_{3}B}{{\tau }_3}}\right)\right]\hfill \\ \hfill & +{\lambda }_8\left[{\psi }_2-Q_n-{\displaystyle \frac{f_{n}T}{{\varphi }_n}}\right]\hfill \\ \hfill & +{\lambda }_9\left[{\phi }_0-{\tau }_0{P}_a^{\max }\right]\hfill \\ \hfill & +{\lambda }_{10}\left[{\phi }_1-{\tau }_1{P}_f^{\max }\right]\hfill \\ \hfill & +{\lambda }_{11}\left[{\phi }_2-{\tau }_2{P}_n^{\max }\right]\hfill \\ \hfill & +{\lambda }_{12}\left[{\phi }_3-{\tau }_3{P}_n^{\max }\right]\hfill \end{array}$$

(38a)

We can use the first-order optimality conditions. Taking the derivative of the Lagrangian function yields

$$\begin{array}{cc}\hfill & {\phi }_0={\left[{\tau }_0{P}_a^{\max }\right]}^{+},\hfill \end{array}$$

(39a)

$$\begin{array}{cc}\hfill & {\phi }_1={\left[{\displaystyle \frac{{\tau }_{1}W\left(Q_f+V+{\lambda }_5\right)}{ln2\left({\lambda }_2+{\lambda }_{10}\right)}}-{\displaystyle \frac{{\tau }_1{\sigma }^2}{g_f}}\right]}^{+},\hfill \end{array}$$

(39b)

$$\begin{array}{cc}\hfill & {\phi }_2={\left[{\displaystyle \frac{{\tau }_2{\lambda }_{4}W}{ln2\left({\lambda }_3+{\lambda }_{11}\right)}}-{\displaystyle \frac{{\tau }_2{\sigma }^2}{g_n}}\right]}^{+},\hfill \end{array}$$

(39c)

$$\begin{array}{cc}\hfill & {\phi }_3={\left[{\displaystyle \frac{{\tau }_{3}W\left(Q_n+V+{\lambda }_7\right)}{ln2\left({\lambda }_3+{\lambda }_{12}\right)}}-{\displaystyle \frac{{\tau }_3{\sigma }^2}{g_n}}\right]}^2,\hfill \end{array}$$

(39d)

By examining the first-order derivatives, we can establish the necessary conditions for optimality. The relationship between the auxiliary variables and the original variables is leveraged to derive the theorem.    □

According to this theorem, we can infer that during the process of radio frequency energy transfer, higher power leads to better results. As the value of W increases, both FU and NU entities are more motivated to offload data, which results in a reduction of the computational tasks performed locally. Furthermore, an increase in V prompts FU and NU to offload a larger portion of their tasks. This shift towards offloading is driven by the desire to enhance the computation rate, which in turn leads the MDS to increase the volume of data offloaded to meet its objectives.

The process of solving the original MSSO problem, denoted as (P1), is encapsulated within Algorithm 1.

Algorithm 1: User-Assisted Dynamic Resource Allocation Algorithm

3.2. Algorithm Complexity Analysis

We have employed the Lyapunov algorithm, which not only ensures the stability of the system but also allows us to effectively decompose the complex overall problem into multiple sub-problems (P2). By implementing Algorithm 1, we are able to solve (P2) within each time slot. Therefore, the solution complexity at each time slot is crucial in determining the overall performance and responsiveness of the algorithm. When solving the problem (P4), we have utilized the interior-point method, which possesses a computational complexity of approximately $O(n^{3.5}log(1/ϵ))$.

In this context, n denotes the total of decision variables, and $ϵ$ reflects the degree of precision in computation. In our algorithm, with $n=8$, the number of decision variables is sufficiently small. This design not only ensures the efficiency of the algorithm but also enables us to complete the optimization of resource allocation within a reasonable time, thus meeting the performance requirements in practical applications.

3.3. Algorithm Performance Analysis

In this section, we demonstrate that the proposed scheme can achieve an optimal long-term time-average solution. First, we establish a key assumption as follows

$$\underset{K\to \infty }{\lim }{\displaystyle \frac{1}{K}}\stackrel{K-1}{{\displaystyle \sum }_{t=0}}{D}_{\mathrm{tot}}\left(t\right)=\overline{D\left(t\right)}$$

(40)

Subsequently, we deduce that the expected value will also converge to the same set of solutions

$$\underset{K\to \infty }{\lim }{\displaystyle \frac{1}{K}}\stackrel{K-1}{{\displaystyle \sum }_{t=0}}\mathbb{E}\left\{[D_{\mathrm{tot}}\left(t\right)\right\}=\overline{D\left(t\right)}$$

(41)

Furthermore, we establish the existence of an optimal solution founded on the existing conditions of the queue as follows.

Lemma 3.

Should problem (P1) be solvable, there is a set of decisions ${\{{\tau }^t,{\mathit{P}}^t\}}^{*}$ that adhere to the conditions $\forall t,\epsilon >0$

$$\mathbb{E}\left\{A_f\left(t\right)-d_f^{\mathrm{loc},*}\left(t\right)-d_f^{\mathrm{off},*}\left(t\right)\right\}\le \epsilon$$

(42)

$$\mathbb{E}\left\{A_n\left(t\right)-d_n^{\mathrm{loc},*}\left(t\right)-d_{n,a2}^{\mathrm{off},*}\left(t\right)\right\}\le \epsilon$$

(43)

$$\mathbb{E}\left\{P_0^{t,*}{\tau }_0^{t,*}-\gamma \right\}\le \epsilon$$

(44)

Here, the asterisk * denotes the value associated with the optimal solution.

Proof. 

Here we omit the proof details for brevity. See the parts 4 and 5 of [44]. □

Theorem 2.

The optimal long-term average weighted computation rate derived from (P1) is bounded below by a lower bound, which is independent of time and space. The algorithm is capable of achieving the following solutions:

• $\overline{D\left(t\right)}\ge \mathbb{E}\left\{D_{tot}^{*}\left(t\right)\right\}-B/V,$
• All queues $Q_f\left(t\right)$, $Q_n\left(t\right)$, $\mathrm{Y}\left(t\right)$ are mean rate stable, thereby satisfying the constraints.

Proof. 

For any $\epsilon >0$, we consider the policy and queue state as defined in equations (40). Given that the values $P_0^{t,*}{\tau }_0^{t,*},d_f^{\mathrm{loc},*}\left(t\right),d_n^{\mathrm{loc},*}\left(t\right),d_{n,a2}^{\mathrm{off},*}\left(t\right),d_f^{\mathrm{off},*}\left(t\right)$, are independent of the queue status $\mathsf{\Theta }\left(t\right)$, we can deduce that

$$\begin{array}{cc}\hfill & \mathbb{E}\left\{A_f\left(t\right)-d_f^{\mathrm{loc},*}\left(t\right)-d_f^{\mathrm{off},*}\left(t\right)\mid \mathsf{\Theta }\left(t\right)\right\}=\mathbb{E}\left\{A_f\left(t\right)\right\}-\mathbb{E}\left\{d_f^{\mathrm{loc},*}\left(t\right)+d_f^{\mathrm{off},*}\left(t\right)\right\}\le \epsilon \hfill \end{array}$$

(45a)

$$\begin{array}{cc}\hfill & \mathbb{E}\left\{A_n\left(t\right)-d_n^{\mathrm{loc},*}\left(t\right)-d_{n,a2}^{\mathrm{off},*}\left(t\right)\mid \mathsf{\Theta }\left(t\right)\right\}=\mathbb{E}\left\{A_n\left(t\right)\right\}-\mathbb{E}\left\{d_n^{\mathrm{loc},*}\left(t\right)+d_{n,a2}^{\mathrm{off},*}\left(t\right)\right\}\le \epsilon \hfill \end{array}$$

(45b)

$$\begin{array}{cc}\hfill & \mathbb{E}\left\{P_0^{t,*}{\tau }_0^{t,*}-\gamma \mid \mathsf{\Theta }\left(t\right)\right\}=\mathbb{E}\left\{P_0^{t,*}{\tau }_0^{t,*}\right\}-\mathbb{E}\left\{\gamma \right\}\le \epsilon \hfill \end{array}$$

(45c)

By integrating these results into equation (29) and taking $\epsilon \to 0$, we obtain

$$\begin{array}{cc}\hfill {\Delta }_V\left(\mathsf{\Theta }\left(t\right)\right)& \le B-V\mathbb{E}\left\{D_{\mathrm{tot}}^{*}\left(t\right)\right\}+Q_f\left(t\right)\mathbb{E}\left\{A_f\left(t\right)-d_f^{\mathrm{loc},*}\left(t\right)-d_{f,n}^{\mathrm{off},*}\left(t\right)\left(t\right)\right\}\hfill \\ \hfill & +Q_n\left(t\right)\mathbb{E}\left\{A_n\left(t\right)-d_n^{\mathrm{loc},*}\left(t\right)-d_{n,a2}^{\mathrm{off},*}\left(t\right)\right\}+\mathrm{Y}\left(\mathrm{t}\right)\mathbb{E}\left\{P_0^{t,*}{\tau }_0^{t,*}-\gamma \right\}\hfill \\ \hfill & \le B-V\mathbb{E}\left\{D_{\mathrm{tot}}^{*}\left(t\right)\right\}\hfill \end{array}$$

(46)

Utilizing the iterated expectation, and summing the aforementioned inequality over the time horizon $t\in \{0,1,...,T-1\}$, we derive the following result

$$\begin{array}{cc}\hfill \mathbb{E}\left\{L\left(\mathsf{\Theta }\left(K\right)\right)\right\}& -\mathbb{E}\left\{L\left(\mathsf{\Theta }\left(0\right)\right)\right\}-V\stackrel{K-1}{{\displaystyle \sum }_{t=0}}\mathbb{E}\left\{D_{tot}\left(t\right)\right\}\hfill \\ \hfill & \le K\left[B-V\mathbb{E}\left\{D_{tot}^{*}\left(t\right)\right\}\right]\hfill \end{array}$$

(47a)

By dividing both sides of equation (47a) by $VK$, applying Jensen’s inequality, and considering that $\mathbb{E}\left\{L\left(\mathsf{\Theta }\left(t\right)\right)\right\}\ge 0$, we obtain

$$\begin{array}{c}\hfill -{\displaystyle \frac{1}{K}}\stackrel{K-1}{{\displaystyle \sum }_{t=0}}\mathbb{E}\left\{D_{tot}\left(t\right)\right\}\le {\displaystyle \frac{B}{V}}-\mathbb{E}\left\{D_{tot}^{*}\left(t\right)\right\}\end{array}$$

(48)

Furthermore, letting $K\to \infty$, we have

$$\begin{array}{ccc}\hfill -\underset{K\to \infty }{\lim }{\displaystyle \frac{1}{K}}\stackrel{K-1}{{\displaystyle \sum }_{t=0}}& \mathbb{E}\left\{D_{tot}\left(t\right)\right\}\le {\displaystyle \frac{B}{V}}\hfill & \hfill -\mathbb{E}\left\{D_{tot}^{*}\left(t\right)\right\}\end{array}$$

(49)

From Equation (41), we have

$$\begin{array}{c}\hfill -\overline{D\left(t\right)}\le {\displaystyle \frac{B}{V}}-\mathbb{E}\left\{D_{tot}^{*}\left(t\right)\right\}\end{array}$$

(50)

Furthermore we obtain

$$\begin{array}{c}\hfill \overline{D\left(t\right)}\ge \mathbb{E}\left\{D_{tot}^{*}\left(t\right)\right\}-{\displaystyle \frac{B}{V}}\end{array}$$

(51)

□

Theorem 3.

A defined upper limit confines the time-averaged sum of queue lengths

$$Q^{\sum }\le {\displaystyle \frac{B}{\epsilon }}+V{\displaystyle \frac{\overline{D\left(t\right)}}{\epsilon }}$$

(52)

Proof. 

By employing the method of iterated expectations and applying telescoping sums iteratively for each $t\in \{0,1,...,T-1\}$, we can derive

$$\begin{array}{cc}\hfill \mathbb{E}\left\{L\left(\mathsf{\Theta }\left(K\right)\right)\right\}& -\mathbb{E}\left\{L\left(\mathsf{\Theta }\left(0\right)\right)\right\}-V\mathbb{E}\left\{D_{tot}\left(t\right)\mid \mathsf{\Theta }\left(t\right)\right\}\hfill \\ \hfill & \le KB-\epsilon \stackrel{K-1}{{\displaystyle \sum }_{t=0}}\mathbb{E}\left\{Q_n\left(t\right)+Q_f\left(t\right)+\mathrm{Y}\left(t\right)\right\}\hfill \end{array}$$

(53)

Dividing both sides of (53) by $K\epsilon$, and considering $T\to \infty$, we rearrange the terms to obtain the desired result

$$\begin{array}{c}\hfill {\displaystyle \frac{B}{\epsilon }}-\underset{K\to \infty }{\lim }{\displaystyle \frac{1}{K}}\stackrel{K-1}{{\displaystyle \sum }_{t=0}}\mathbb{E}\left\{Q_n\left(t\right)+Q_f\left(t\right)+\mathrm{Y}\left(t\right)\right\}+V{\displaystyle \frac{\overline{D\left(t\right)}}{\epsilon }}\ge 0\end{array}$$

(54)

Specifically

$$\begin{array}{c}\hfill Q^{\sum }\le {\displaystyle \frac{B}{\epsilon }}+V{\displaystyle \frac{\overline{D\left(t\right)}}{\epsilon }}\end{array}$$

(55)

□

Theorems 2 and 3 underpin our proposed algorithm by establishing that as V increases, the computation rate $\overline{D\left(t\right)}$ improves at a rate of $O\left(1/V\right)$, whereas the queue length grows at a rate of $O\left(V\right)$. This suggests that by choosing an appropriately large value for V, we can achieve an optimal $\overline{D\left(t\right)}$. Furthermore, the time-averaged queue length, denoted as $Q^{\sum }$, is shown to increase linearly with V. This linear relationship implies a trade-off $\left[O\left(1/V\right),O\left(V\right)\right]$ between the optimized computation rate and the queue length, . In line with Little’s Law [47], which posits that delay is proportional to the time-averaged length of the data queue. In other words, our proposed algorithm enables a trade-off between $\overline{D\left(t\right)}$ and the average network latency.

4. Simulation Results

In this section, plenty of numerical simulations are performed to assess the efficiency of our proposed algorithm. Our experiments were conducted on a computational platform with an Intel(R) Xeon(R) Gold 6148 CPU 2.40 GHz, 20 cores and four GeForce RTX 3070 GPUs. In our simulations, we employed a free-space path-loss model to depict the wireless channel characteristics [48]. The averaged channel gain $\overline{h}$ here is denoted as

$$\overline{h}=A_d{\left({\displaystyle \frac{3\times {10}^8}{4\pi f_c{d}_i}}\right)}^{d_e},$$

(56)

where $A_d$ denotes the antenna gain, $f_c$ denotes the carrier frequency, $d_e$ denotes the path loss exponent, and $d_i$ in meters denotes the distance between two nodes. The time-varying WPT and task offloading channel gains are represented by the vector ${\mathbf{h}}^{\mathbf{t}}=\left[a_1^t{h}_m^t,a_2^t{h}_h^t,a_3^t{g}_m^t,a_4^t{g}_h^t\right]$, adhering to the Rayleigh fading channel model. In this model, the random channel fading factors $a_i$ are characterized by an exponential distribution with a unit mean, capturing the variability inherent in wireless communication channels. For the sake of simplicity, we assume that the vector of fading factors $\left[a_1^t,a_2^t,a_3^t,a_4^t\right]$ is constant and equal to $\left[1.0,1.0,1.0,1.0\right]$ for any given time slot, implying that the channel gains are static within that slot. The interval between task arrivals at FU and NU both follow an exponential distribution with constant average rates ${\lambda }_f$ and ${\lambda }_n$ respectively. The parameters are all listed in Table 2.

4.1. Impact of System Parameters on Algorithm Performance

Figure 3 illustrates trend for the average task computation rate $D_{tot}$ and the average task queue lengths of FU and NU over a 5000 time slots period. The the task arrival rates for FU and NU are set at 1.75 Mbps and 2 Mbps, respectively. Initially, $D_{tot}$ is low, but it rapidly increases and eventually stabilizes as time progresses. This initial surge in $D_{tot}$ is due to the system’s initial adjustment to the initial task queue fluctuations, which demands more intensive processing for FU tasks, resulting in increased energy consumption and a temporary reduction in the overall computation rate. Moreover, the average queue length decreases and stabilizes, reflecting the system’s ability to self-regulate and reach a steady state.

Figure 4 demonstrates the average task computation rate of our proposes algorithm under different control parameter $V=100,300,500,700$.The results show that the average task computation rates converge similarly across different V. Notably, as V increases, there is a corresponding increase in the average task computation rates. This trend is attributed to the fact that a larger V compels the algorithm to prioritize computation rates over queue stability, which is consistent with theoretical analysis, which is corresponding to Theorems 2. Here, the parameter V serves as a balancing factor between the task computation rate and queue stability, reflecting a trade-off that is consistent with our theoretical predictions.

Figure 5 shows the trend of average task queue lengths of FU and NU after with different V. As V increases from 100 to 900, the task queue lengths of FU and HD declines from $3.1\times {10}^8$ bits to $2.9\times {10}^8$ bits and from $3.11\times {10}^8$ bits to $2.84\times {10}^8$ bits respectively. In the user-assisted offloading paradigm, processing a task from NU involves only a single offloading step, which is markedly more efficient than the two-step offloading data from FU. Consequently, with the increasing of V, indicating a shift in the algorithm’s focus towards prioritizing computation rates over queue stability, the algorithm will naturally prefer to offload tasks from NU, resulting in a more rapid decline in the NU’s task queue length. By tuning the value of V, an optimal balance in task processing between FU and HU can be effectively achieved, reflecting the algorithm’s adaptability to different operational priorities.

Figure 6 evaluates the impact of the average energy constraint $\gamma$ on system performance with V fixed at 500. As $\gamma$ increases from $1.75$ joule to $2.25$ joule, the average task computation rate rises from $1.79\times {10}^6$ bits/s to $1.86\times {10}^6$ bits/s, while the average task queue length of FU and NU decreases from $4\times {10}^8$ bits to $2.5\times {10}^8$ bits. This reduction in queue length and increase in computation rate are attributed to the higher energy availability for WPT, enabling FU and NU to offload tasks more effectively. Notably, after the average energy constraint $\gamma$ reaching 2.1 joule, the variation of task computation rate and task queue is reduced. This observation suggests that there is an upper bound to the energy consumption of our algorithm. Beyond this threshold, additional energy has minimal impact on performance enhancement. Hence, energy constraint $\gamma$ as a critical parameter, significantly influences both the data processing rates and the stability of task queues within the system. The findings underscore the importance of energy management in optimizing system performance.

Figure 7 presents the total system energy consumption of our proposed algorithm under different values of the parameter V, specifically $V=300$ and $V=700$. Initially, the total energy consumption exhibits substantial fluctuations. However, with the progression of time, the system’s energy consumption stabilizes and hovers around the average energy constraint $\gamma$ after approximately 2500 time slots. Notably, an elevated V value is correlated with higher average energy consumption, as the algorithm pay more attention to the system’s computation rate, consequently incurring greater energy costs. Figure 7 highlight the algorithm’s efficacy in managing average energy consumption, a critical feature for the sustainability of IoT networks.

In Figure 8, we evaluate the offloading power across varying bandwidths W, ranging from $0.96\times {10}^6$ Hz to $1.05\times {10}^6$ Hz. It is observed that all the offloading powers increases as bandwidth W escalates. Consistent with the results of the analysis in Theorem 1, the increase in bandwidth makes the system more inclined to perform task offloading, reflected in the increase in offloading power.

4.2. Comparison with Baseline Algorithms

To evaluate the performance of our propose algorithm, we choose the the following three representative benchmarks as the baseline algorithms.

(1) All offloading scheme: Both FU and NU do not perform local computing and consumes all the energy for task offloading.

(2) No cooperation scheme: FU offloads tasks directly to HAP without soliciting assistance from NU, similar to the method in [49].

(3) No lyapunov scheme: Disregarding the dynamics of task queues and the energy queue, this scheme focuses solely on maximizing the average task computation rate, similar to the Myopic method in [30]. To ensure a fair comparison, we constrain the energy consumption of each time slot by solving.

$$\begin{array}{cc}\hfill max& D_{tot}\left(t\right)\hfill \\ \hfill s.t.& \left(22b\right)-\left(22d\right),\left(22g\right)-\left(22k\right)\hfill \\ \hfill & P_0^t{\tau }_0^t\le \gamma \hfill \end{array}$$

(57)

Figure 9 shows the average task computation rates under four schemes over a 5000 time-slots period, with the control parameter V set to 500. All schemes converge after 1000 time slots. Our proposed algorithm achieves the best task computation rate after convergence, outperforming the other three by 0.8%, 3.9%, and 4.1% respectively. Our algorithm’s key strengths lie not only in achieving the highest data processing rate but also in ensuring the stability of the system queues. This prevents excessively long queues that could lead to prolonged response times and a negative user experience. The no-Lyapunov scheme, while achieving the second-highest computation rate, neglects queue stability in its pursuit of maximizing computation speed. This oversight can lead to system instability, prolonged user service times, and potential system failure. The all-offloading scheme, relying solely on edge computing, consumes more energy and thus underperforms in energy-limited scenarios. In the no-cooperation scheme, the system initially benefits from a high computation rate due to the NU’s lack of assistance to the FU. But as the NU’s tasks are completed and its resources are no longer available, the average computation rate falls sharply. The FU’s communication with the HAP is further impeded by the dual proximity effect, causing a notable decline in the system’s long-term computation performance.

Figure 10 show the impact of varying network bandwidth W from $0.96\times {10}^6$ Hz to $1.05\times {10}^6$ Hz on the performance of different schemes. As W increase, task computation rates for all schemes rise, reflecting improved transmission efficiency for wireless power transfer and task offloading. This allows HAP to handle more offloaded tasks, highlighting the critical role of bandwidth in system performance. Notably, our proposed scheme consistently outperforms others across all bandwidth levels, showcasing its adaptability and robustness in varying network conditions.

Figure 11 evaluates how the distance between FU and NU affects the performance of all four schemes, with distances varying from 120 m to 160 m. We observe that as the distance increases, the computation rates for both our proposed scheme and the all-offloading scheme decrease. This suggests that proximity plays a crucial role in task offloading efficiency. In contrast, the no-cooperation scheme shows a stable computation rate, consistent with its design that excludes task offloading between FU and NU. Interestingly, the no-Lyapunov scheme performs best at a distance of about 140 meters. However, its performance drops as the distance decreases, contrary to the expectation that a shorter distance would enhance task offloading from FU to NU. This unexpected trend is likely due to instances where the FU’s task queue depletes faster than new tasks arrive, leading to lower computation rates for the no-Lyapunov scheme. This highlights the importance of balancing task computation rates with queue stability in system design.

In Figure 12, we evaluate the performance of four schemes as the task arrival rate of NU varies from $1.85\times {10}^6$ bits/s to $2.35\times {10}^6$ bits/s. Our proposed scheme’s task computation rate demonstrates a modest increase and maintains the highest computation rate as tasks arrive more rapidly. This trend underscores the scheme’s robustness across diverse scenarios. Correspondingly, the no cooperation scheme exhibits a more pronounced increase in task computation rate. This is attributable to its vigorous task processing capacity at the NU, which allows it to capitalize on the higher task arrival rates effectively.

5. Conclusions

The joint optimization of computation offloading and resource allocation in WPMEC systems poses a significant challenge due to the time-varying network environments and the time-coupling nature of energy consumption constraint. In this study, we aim at maximizing the average task computation rate in an MEC system with WPT and user collaboration. A task computation maximization problem was formulated considering uncertainty of load dynamics alongside energy consumption constraint. The complexity of this problem is amplified by the interplay of multiple coupled parameters. To surmount these difficulties, we introduce an innovative online algorithm, named UAORA. This algorithm leverages Lyapunov optimization theory to manage long-term task queue lengths and energy consumption, effectively transforming the sequential decision-making dilemma into a series of deterministic optimization problems for each time slot. Extensive simulation results substantiate the effectiveness of the proposed UAORA algorithm, demonstrating a significant enhancement in average task computation performance when compared with benchmark methods. Simulations also underscores the advantages of jointly consider the task computation rate and the task queues in our algorithm. An enticing avenue for future research is to account for dynamically changing local computation rates. Incorporating such variability could potentially augment the system’s flexibility and applicability within the complexities of real-world scenarios.