CN106714301A - Carrier optimization method in wireless positioning network - Google Patents

Abstract

The invention discloses a MILP (Mixed Integer Linear Programming) solution method of a carrier resource optimization problem. When a MIP (Mixed Integer Programming) model is solved by a YALMIP tool, firstly, a MIP problem is converted into MILP, and thus, a target function needs to be subjected to linear approximate processing. A global SPEB (Squared Position Error Bound) is used as a target function; under the conditions that lower bounds and upper bounds of transmitting powder, signal bandwidths and carrier frequencies of nodes are limited, total available power in an integral system is limited, a total bandwidth in an integral wireless network is limited and signal carrier frequencies and bandwidths among the nodes do not interfere mutually, allocation resources obtained by optimization comprise transmitting powder, carrier frequency and a signal bandwidth; and then a destination node acquires position information of the destination node by carrying out one-way TOA (Time of Advent) ranging with anchor nodes and cooperative destination nodes. The method disclosed by the invention can rapidly and accurately achieve the joint power and spectrum resource optimization of non-cooperative and cooperative wireless positioning networks.

Description

A kind of carrier wave optimization method in wireless positioning network

Technical field

The present invention relates to the carrier wave optimization side in wireless positioning network technical field, more particularly to a kind of wireless positioning network Method.

Background technology

With the development of applications of wireless technology, the service based on positioning is increasingly paid close attention to by people.It is modern at some In the application scenarios of radio communication service and equipment, such as the path in indoor positioning and navigation, medical aid work, car networking Planning, tracking of logistics etc., the acquisition of high precision position information play more and more important effect.GLONASS (Global Navigation Satellite System, GNSS) is the most well known technological means for providing location aware, but It is that the technology can fail in some extreme environment performances, such as：High building room block under modern metropolitan cities, indoor scene, underground Environment etc..It is a kind of using can effectively make up the deficiency of GNSS using the positioning network technology of wireless signal, but wireless location Network can there are problems that resource using being restricted, anchor node covers and becomes study hotspot in recent years.

As shown in Figure 1, two class nodes are generally comprised in wireless positioning network：Anchor node known to position (anchor) With it needs to be determined that the destination node (agent) of position.Destination node try to the positional information of oneself by with anchor node or The distance between person's cooperation agent nodes measurement such as arrival time (TOA), angle of arrival (AOA), received signal strength (Received Signal Strength Intensity, RSSI) etc. is completed.In existing research contents, one kind is based on etc. The form of valency Fisher's information matrix (Equivalent Fisher Information Matrix, EFIM) gives wireless fixed Under the network condition of position, square position error lower bound with compact form (Squared Position Error Bound, SPEB) definition.SPEB can be as the index for weighing broadband signal positioning precision, and its essence is CramerRao lower limit CRLB A kind of expression-form.On this basis, can be drawn the following conclusions according to existing research, the resource in positioning network is as sent out Penetrating the introducing cooperated between power, signal bandwidth, the reasonable distribution of carrier frequency and suitable agent nodes can be largely The upper performance for improving positioning.

Non-cooperating positions co-allocation (the Joint Power and Spectrum of network power and frequency spectrum Allocation, JPSA) propose to optimize the transmission power of network node simultaneously in research, signal bandwidth, carrier frequency is obtaining High-precision positioning performance, studying us according to this can show that carrier wave is important based on during TOA network positions Factor, however, because in solution carrier wave optimization process, the non-convex situation on carrier frequency constraints is exhaustive using one kind Method (Full Enumeration, FE) solves the problems, such as above-mentioned JPSA.The so-called method of exhaustion is all possible carrier waves of node of traversal Combination of frequency situation, although it has preferable solving precision, but its solution ability is limited by nodes quantity, If anchor node number be on the increase or co-positioned in cooperative node introducing, the complication of network will cause to calculate and solve Complexity with the factorial form order of magnitude improve, it is difficult to obtain optimal solution in relative time.

The content of the invention

The purpose of the present invention is to propose to the carrier wave optimization method in a kind of wireless positioning network, asking for carrier wave optimization is being solved Optimal solution can be more rapidly drawn during topic, to adapt to the increase of node high quantity, is omited while proposing co-positioned network JPSA and surveying Middle solution mode, highlights advantage of the co-positioned in wireless positioning network.

It is that, up to above-mentioned purpose, the present invention is achieved through the following technical solutions：

A kind of carrier wave optimization method in wireless positioning network, carrier wave is excellent in being applied to non-cooperating positioning network Change, the set of destination node and anchor node in the positioning network is expressed asWithThe positional representation of node k is p_k=[x_k,y_k]^T,Using square Error lower bound SPEB weighs the positioning precision of whole network, and it is derived by equivalent Fei Sheer information matrixs EFIM；λ_kjNode k Ranging information intensity RII and j between is the inverse of range error CRLB,

Wherein β_jIt is the effective bandwidth of node j, ξ_kjIt is the related parameter of signal propagation channel, is one for representing transmission The positive number of the characteristic of channel, P_jFor representing the transmission power of node j,It is the path attenuation between node k and j；Non-cooperating is determined The SPEB of position location estimation of the network on destination node is usedIt is right to representOne is carried out in a zonule R Rank Taylor linear launches, and makes θ=[P^T,B^T,f^T]^TTo simplify the expression of symbol, then

And θ^(m-1)It is the solution of (m-1) secondary iterative, equally it can be applied in next iteration expansion；Non-cooperating is positioned Network power is with the co-allocation model of frequency spectrum：

The MILP Mathematical Modelings ζ of non-cooperating positioning network JPSA₃:

s.t.P_min≤P_j≤P_max

B_min≤B_j≤B_max

f_min≤f_j≤f_max

||θ-θ^(m-1)| |=| | Δ θ | |≤R；

Wherein M is sufficiently large positive number, and u is binary variable, u ∈ { 0,1 },

φ_kjIt is the angle between node k and node j, φ_kiIt is the angle between node k and node i；

Non-cooperating positioning network JPSA strategy MILP method for solving idiographic flows are as follows：

S1：Selection initial value, θ=θ^(m), m=0；

S2：MILP problems ζ is solved with YALMIP instruments₃, output Δ θ；

S3：Update θ^(m+1)=θ^(m)+ Δ θ, m=m+1；

S4：Whether≤R is restrained for | | the Δ θ | | that judges the condition of convergence, and S5 is jumped to if setting up, and otherwise returns to S3；

S5：Export the result of SPEB and θ.

A kind of carrier wave optimization method in wireless positioning network, carrier wave is excellent in being applied to co-positioned network Change, the set of destination node and anchor node in the positioning network is expressed asWithThe positional representation of node k is p_k=[x_k,y_k]^T,Using square Error lower bound SPEB weighs the positioning precision of whole network, and it is derived by equivalent Fei Sheer information matrixs EFIM；λ_kjNode k Ranging information intensity RII and j between is the inverse of range error CRLB,

Wherein β_jIt is the effective bandwidth of node j, ξ_kjIt is the related parameter of signal propagation channel, is one for representing transmission The positive number of the characteristic of channel, P_jFor representing the transmission power of node j,It is the path attenuation between node k and j；On matrix Invert mark theorem：N rank symmetric positive definite matrixs A, μ₁=tr (A), μ₂=tr (A*A), a is the minimal eigenvalue of matrix A, then to square Battle array is inverted and asks the computing of mark to have following form again

The SPEB of the destination node k of collaborative network is defined as：Use in its upper boundCarry out table Show, wherein

WhereinThe EFIM expression matrix forms that destination node location is estimated in co-positioned network are represented, it and non-association Make network difference, because the introducing of communication distance measuring between destination node is not such as the simple diagonal formation of non-cooperating scene EFIM matrixes Formula；With the upper bound of SPEBUsed as the object function of model, the MIP Mathematical Modelings of co-positioned network JPSA are：

ζ₄:

s.t.B_min≤B_j≤B_max

f_min≤f_j≤f_max

Wherein M is sufficiently large positive number, and u is binary variable, u ∈ { 0,1 }；Above-mentioned model is still unsatisfactory for MILP solutions Model, therefore JPSA problems in collaborative network are solved using a kind of iterative linear algorithm IL, comprise the following steps that：

S1：Selection initial value, θ=θ^(m-1), m=1 must according to formula (*)With the value of a；

S2：Make θ^(m)=θ^(m-1)+ Δ θ, solves the upper dividing value of SPEB

WhereinIt is μ₁In a first order Taylor expansion value of zonule R, by object functionReplace withIt is right to representThe result after linear process is carried out, wherein, make θ=[P^T,B^T,f^T]^TAccorded with simplifying Number expression, and adds constraints on this basis | | θ-θ^(m-1)| |=| | Δ θ | |≤R, thus set up and complete a MILP mould Type problem ζ₄′；

S3：ζ is solved using YALMIP instruments₄', export θ^(m)=θ^(m-1)+ Δ θ, updates m=m+1；

S4：Whether≤R is restrained for | | the Δ θ | | that judges the condition of convergence, and S5 is jumped to if setting up, and otherwise returns to S2；

S5：Export the result of SPEB and θ.

The beneficial effects of the invention are as follows：The present invention proposes the MILP in a kind of carrier resource optimization problem Method (Mixed Integer Linear Programming, MILP), the proposition of the method can quickly, accurately solve non-association Make and JPSA optimization problems in co-positioned network.

Brief description of the drawings

Fig. 1 is wireless positioning network model schematic；

Fig. 2 is the SINC time domains and frequency-domain waveform figure of normalized energy；

Fig. 3 is FE and MILP solving results comparison diagram of the invention；

Fig. 4 is the SPEB result figures under different resource allocation strategy.

Specific embodiment

The present invention is described in further detail below by specific embodiment combination accompanying drawing.

There is setting N in a 2-D wireless positioning network_aIndividual destination node, N_bAnchor node node known to individual position.It is false If all nodes (agent and anchor) clock synchronization in network, the set of destination node and anchor node can be represented respectively ForWithThe position of node k can be expressed as p_k=[x_k,y_k ]^T,In measurement process, anchor node (may participate in the agent nodes of cooperation) broadcast transmission is by optimizing The distribution resource that goes out, transmission power, carrier frequency, signal bandwidth, then destination node by with anchor node and cooperate Agent nodes carry out unidirectional TOA range findings to obtain the positional information of oneself.

In synchronously positioning network, the unidirectional range error under multi-path environment can be created as Gaussian distribution model.Range finding Performance can represent that it is range finding with ranging information intensity (Ranging Information Intensity, RII) The inverse of error CRLB.RII between node k and j can be expressed as：

Wherein β_jIt is the effective bandwidth of node j, c is the light velocity, χ_kj∈ (0,1) is the overlapping coefficient in path between node k and j (Path Overlap Coefficient, POC), it is used for describing the influence of Multipath Transmission.It is in the letter for receiving Number r_jThe signal to noise ratio (SNR) of l paths in (t).By corresponding abbreviation, the related parameter ξ of signal propagation channel_kjCan use One positive number represents transmission channel characteristic, and P is used for representing the transmission power of node j, the path attenuation between node k and j can With withTo represent.

The present invention weighs whole net using mean square error lower bound (squared position error bound, SPEB) The positioning precision of network, its by equivalent Fei Sheer information matrixs (Equivalent Fisher Information Matrix, EFIM) derive, i.e., so-called CRLB.The SPEB of destination node k is defined as：

WhereinIt is p_kEstimated location.J_e(p_k) it is EFIM of the agent k nodes by measurement acquisition.

In co-positioned network, N_aThe 2N that the EFIM of individual destination node can be written as_a×2N_aMatrix form, should (k, j) individual element of EFIM is：

In (3),And C_kjThe ranging information between node k and j is represented,

Wherein q_kj=[cos (φ_kj),sin(φ_kj)]^TRepresent the angle information between node k and j.SPEB describes position The minimum lower limit of estimated accuracy is put, therefore can be as the good and bad parameter of the performance for weighing positioning network.

In formula (1), the definition of effective bandwidth β is：

If using the SINC waveforms of normalized energy, as shown in Figure 2, the frequency domain of the waveform can be expressed as：

So effective bandwidth β²Can be rewritten as：

Wherein B represents actual bandwidth, and f represents carrier frequency, and formula (8) is updated in formula (1).Between node k and j RII can be write as：

According to the result of study having had, JPSA solving-optimizings problem into following Mathematical Modeling can be set up：

ζ₁:

s.t.P_min≤P_j≤P_max (11)

B_min≤B_j≤B_max (12)

f_min≤f_j≤f_max (13)

Model ζ₁In,i∈N_b, j ≠ i.Overall situation SPEB is used as object function for the present invention, and (11)-(13) represent basis The design requirement of system each node is owned by oneself transmission power, signal bandwidth, the lower bound and upper limit of carrier frequency System, (14) represent available total power limitation in whole system, and (15) represent total bandwidth limitation in whole wireless network, and (16) are used for Expression carrier wave is frequent, the restrictive condition of signal bandwidth non-interference condition, it is ensured that do not overlapped in bandwidth between each node.

In model ζ₁In, most difficult problem is that the absolute value expression-form of constraints (16) is non-convex, therefore ζ₁ The problem, i.e. master mould can be processed with the method for exhaustion (Full Enumeration, FE) can be by solving K！Individual subproblem (K represents the ＜ K ＜ N of number of nodes 0 for participating in positioning distance measuring_a+N_b), one of them typical subproblem model is：

ζ₁ ^FE:

s.t.(11)-(15) (18)

All solutions for solving subproblem draw master mould ζ by contrast₁Optimal solution.This traversal method of exhaustion can be examined Consider each possibility, obtain relatively accurate result, but when anchor node number than it is larger when, computational complexity will very It is high.

In order to reduce computation complexity, constraints (16) can be modified to

Wherein M is sufficiently large positive number, and u is binary variable, u ∈ { 0,1 }.So ζ₁Can be changed to：

ζ₂:

s.t.(11)-(15),(21)-(22) (24)

ζ₂An also mixed integer programming problem for non-convex (Mixed Integer Programming, MIP), But be available with for example well known branch and bound method (Branch and Bound Methods) of the algorithm of comparative maturity while Existing solution instrument MOSEK, YALMIP can effectively realize the solution to algorithm.

The present invention solves ζ with YALMIP instruments₂, it is necessary to MIP problems are converted into linear programming (Linear during model Programming, LP), it is therefore desirable to propose that object function is converted into linear structure by the approximate evaluation method of linearisation.

Carrier wave optimization in non-cooperating positioning network

Proposition 1：In non-cooperating positioning network, the SPEB's of destination node k can be expressed as following approximate form：

According to proposition 1, can be rightFirst order Taylor linear expansion is carried out in a zonule R.Use θ=[P^T, B^T,f^T]^TTo simplify the expression of symbol, therefore nonlinear object function (23) can be rewritten into：

Wherein,

||θ-θ^(m-1)| |=| | Δ θ | |≤R (26)

And θ^(m-1)It is the solution of (m-1) secondary iterative, equally it can be applied in next iteration expansion.

Therefore ζ₁MILP models can be expressed as following form：

ζ₃:

s.t.(11)-(15),(21)-(22),(26) (28)

Algorithm 1：A kind of MILP method for solving of non-cooperating positioning network JPSA carrier wave optimizations is as follows：

S1：Selection initial value, θ=θ^(m), m=0.

S2：MILP problems ζ is solved with YALMIP instruments₃, output Δ θ

S3：Update θ^(m+1)=θ^(m)+ Δ θ, m=m+1.

S4：Whether≤R is restrained for | | the Δ θ | | that judges the condition of convergence, and S5 is jumped to if setting up, and otherwise returns to S3.

S5：Export the result of SPEB and θ.

Carrier wave optimization in co-positioned network

Determine network for cooperation, EFIM matrixes are no longer diagonal matrixs long, so the SPEB expression-forms under the collaborative network Will not the form as given by proposition 1 again, it is therefore desirable to propose another linearization approximate algorithm by collaborative network JPSA problems Object function carries out linear process.

Theorem 1:The upper bound of matrix inversion mark：

Assuming that n rank symmetric positive definite matrixs A, μ₁=tr (A), μ₂=tr (A*A), a are the minimal eigenvalues of matrix A, then

According to the statement of theorem 1, the upper bound of the SPEB expression-forms of collaborative network can be usedTo express, among these

WhereinRepresent the EFIM expression matrix forms in co-positioned network.

First with the upper bound of SPEBUsed as the object function of master mould, therefore original model can be rewritten as：

ζ₄:

s.t.(11)–(15),(21)–(22) (31)

Now, it is possible to find ζ₄In object functionIt is still non-linear, so the present invention is proposed on this basis A kind of iterative linear algorithm (Iterative Linearization, IL) carries out linear process.

Algorithm 2：A kind of IL method for solving of co-positioned network JPSA carrier waves optimization is as follows：

S1：Selection initial value, θ=θ^(m-1), m=1. must according to formula (29)With the value of a.

S2：Make θ^(m)=θ^(m-1)+ Δ θ, solves the upper dividing value of SPEB

WhereinIt is μ₁In a first order Taylor expansion value of zonule R, by object functionReplace withIt is right to representThe result after linear process is carried out, wherein making θ=[P^T,B^T,f^T]^TAccorded with simplifying Number expression, and on this basis addition constraints (26), thus set up complete a MILP model problems ζ₄′。

S3：ζ is solved using YALMIP instruments₄', export θ^(m)=θ^(m-1)+ Δ θ, updates m=m+1.

S4：Whether≤R is restrained for | | the Δ θ | | that judges the condition of convergence, and S5 is jumped to if setting up, and otherwise returns to S2.

S5：Export the result of SPEB and θ.

In order to prove the lifting of carrier wave optimization method proposed by the invention to JPSA problem solving overall performances, the present invention Verified by emulation.Under equivalent environment, compared for MILP methods and FE methods process carrier frequency frequency optimize when Solving speed.The contrast of the solving result precision of two methods is also given simultaneously.MILP methods proposed by the invention can With the JPSA in effectively solving the problems, such as under co-positioned network environment.Finally compared for successively under Different Strategies resource allocation Position error SPEB effects.In simulation process, the power of all destination nodes and anchor node, signal bandwidth, carrier frequency quilt Normalization, topological structure is arranged as having N in the square region of 10 × 10_aIndividual agent nodes and N_b=3 anchor nodes.

Table 1：FE and MILP solving speeds contrast (unit：Second)

Accompanying drawing 3 compared for being utilized respectively the solution JPSA treatment carrier wave optimizations of MILP and FE methods in co-positioned network As a result, although FE methods can effectively try to achieve optimal solution, but work as N_aNumber than it is larger when, solution procedure very take, such as accompanying drawing 3 It is shown, it is also possible to find two methods process JPSA problems when, it is as a result close, while it can be seen that when agent quantity more than 4 When using FE methods solve carrier wave optimization problem it is considerably complicated and also be difficult relative time obtain optimal solution, work as utilization During MILP methods, working as N_aWhen quantity is more than 6, still can be with normal work.Two methods complexity is relatively by phase With the run time that program is given under emulation platform.According to table 1, when the run time of MILP is only solved just with FE methods 1% or so, therefore MILP have solution efficiency higher.

Accompanying drawing 4 compares three kinds of different resource allocative decisions and obtains position error SPEB values, including, non-cooperating positioning net JPBA (Joint Power and Bandwidth Allocation, JPBA) in network, JPSA, in co-positioned network JPSA.Can draw the following conclusions, first no matter which kind of Resource Allocation Formula, position error all can be with the quantity of agent Increase and increase, this is relevant with the limitation of general power, secondly, the performance that two kinds of JPSA schemes all can be than JPBA is good, this with show Some work sutdy results are consistent, and this is primarily due to the phase information that JPSA resource allocation policies make use of carrier wave.Finally, Cooperating between agent and agent has very great help to the lifting of positioning precision, thus under co-positioned network condition JPSA is the optimal policy in three kinds of strategies.

Above content is to combine specific preferred embodiment further description made for the present invention, it is impossible to assert Specific implementation of the invention is confined to these explanations.For general technical staff of the technical field of the invention, On the premise of not departing from present inventive concept, some simple deduction or replace can also be made, should be all considered as belonging to of the invention Protection domain.

Claims (3)

1. the carrier wave optimization method in a kind of wireless positioning network, is applied to carrier wave optimization, its feature in non-cooperating positioning network It is：The set of destination node and anchor node in the positioning network is expressed asWithThe positional representation of node k is p_k=[x_k,y_k]^T,Using square Error lower bound SPEB weighs the positioning precision of whole network, and it is derived by equivalent Fei Sheer information matrixs EFIM；λ_kjNode k Ranging information intensity RII and j between is the inverse of range error CRLB,

${\mathrm{\λ}}_{kj}={\mathrm{\ξ}}_{kj}\frac{P_j{\mathrm{\β}}_j^2}{d_{kj}^2}$

Wherein β_jIt is the effective bandwidth of node j, ξ_kjIt is the related parameter of signal propagation channel, is one for representing transmission channel The positive number of characteristic, P_jFor representing the transmission power of node j,It is the path attenuation between node k and j；Non-cooperating positions net The SPEB of location estimation of the network on destination node is usedIt is right to representSingle order is carried out in a zonule R safe Linear expansion is strangled, θ=[P is made^T,B^T,f^T]^TTo simplify the expression of symbol, then

And θ^(m-1)It is the solution of (m-1) secondary iterative, equally it can be applied in next iteration expansion；Non-cooperating positions network Power is with the co-allocation model of frequency spectrum：

The MILP Mathematical Modelings ζ of non-cooperating positioning network JPSA₃:min.

s.t.P_min≤P_j≤P_max

B_min≤B_j≤B_max

f_min≤f_j≤f_max

$(f_j-f_i)+M*u\≥\frac{1}{2}(B_j+B_i)$

$-(f_j-f_i)+M*(1-u)\≥\frac{1}{2}(B_j+B_i)$

||θ-θ^(m-1)| |=| | Δ θ | |≤R；

Wherein M is sufficiently large positive number, and u is binary variable, u ∈ { 0,1 },

φ_kjIt is the angle between node k and node j, φ_kiIt is the angle between node k and node i；

Non-cooperating positioning network JPSA strategy MILP method for solving idiographic flows are as follows：

S1：Selection initial value, θ=θ^(m), m=0；

S2：MILP problems ζ is solved with YALMIP instruments₃, output Δ θ；

S3：Update θ^(m+1)=θ^(m)+ Δ θ, m=m+1；

S4：Whether≤R is restrained for | | the Δ θ | | that judges the condition of convergence, and S5 is jumped to if setting up, and otherwise returns to S3；

S5：Export the result of SPEB and θ.

2. method according to claim 1, it is characterised in that：The SPEB of destination node k is defined as：

WhereinIt is p_kEstimated location, J_e(p_k) it is EFIM of the agent k nodes by measurement acquisition.

3. the carrier wave optimization method in a kind of wireless positioning network, is applied to carrier wave optimization in co-positioned network, The set of destination node and anchor node in the positioning network is expressed asWithThe positional representation of node k is p_k=[x_k,y_k]^T,Using square Error lower bound SPEB weighs the positioning precision of whole network, and it is derived by equivalent Fei Sheer information matrixs EFIM；λ_kjNode k Ranging information intensity RII and j between is the inverse of range error CRLB,

${\mathrm{\λ}}_{kj}={\mathrm{\ξ}}_{kj}\frac{P_j{\mathrm{\β}}_j^2}{d_{kj}^2}$

Wherein β_jIt is the effective bandwidth of node j, ξ_kjIt is the related parameter of signal propagation channel, is one for representing transmission channel The positive number of characteristic, P_jFor representing the transmission power of node j,It is the path attenuation between node k and j；On matrix inversion Mark theorem：N rank symmetric positive definite matrixs A, μ₁=tr (A), μ₂=tr (A*A), a is the minimal eigenvalue of matrix A, then to Matrix Calculating There is following form against the computing of mark is asked again

$tr\left(A^{-1}\right)\≤\left[\begin{array}{cc}{\mathrm{\μ}}_1& n\end{array}\right]*{\left[\begin{array}{cc}{\mathrm{\μ}}_2& {\mathrm{\μ}}_1\\ a^2& a\end{array}\right]}^{-1}*\left[\begin{array}{c}n\\ 1\end{array}\right],$

The SPEB of the destination node k of collaborative network is defined asUse in its upper boundCarry out table Reach, wherein

$\begin{array}{cccc}{\mathrm{\μ}}_1=tr\left(J_e^{co}\right),& {\mathrm{\μ}}_2=tr(J_e^{co}*J_e^{co}),& a=\min eig\left(J_e^{co}\right)& (*)\end{array}$

The EFIM expression matrix forms of destination node location estimation in co-positioned network are represented, it is with un-coordinated net frk not Together, because the introducing of communication distance measuring between destination node is not such as the simple diagonal matrix form SPEB of non-cooperating scene EFIM matrixes The upper boundUsed as the object function of model, the MIP Mathematical Modelings of co-positioned network JPSA are：

ζ₄:min.

s.t.B_min≤B_j≤B_max

f_min≤f_j≤f_max

$(f_j-f_i)+M*u\≥\frac{1}{2}(B_j+B_i)$

$-(f_j-f_i)+M*(1-u)\≥\frac{1}{2}(B_j+B_i),$

Wherein M is sufficiently large positive number, and u is binary variable, u ∈ { 0,1 }；Above-mentioned model is still unsatisfactory for MILP and solves mould Type, therefore JPSA problems in collaborative network are solved using a kind of iterative linear algorithm IL, comprise the following steps that：

S1：Selection initial value, θ=θ^(m-1), m=1 must according to formula (*)With the value of a；

S2：Make θ^(m)=θ^(m-1)+ Δ θ, solves the upper dividing value of SPEB

${\tilde{P}}^U\left(p_k\right)=\left[\begin{array}{cc}{\widetilde{\mathrm{\μ}}}_1^{\left(m\right)}& n\end{array}\right]*{\left[\begin{array}{cc}{\mathrm{\μ}}_2^{(m-1)}& {\mathrm{\μ}}_1^{(m-1)}\\ a^2& a\end{array}\right]}^{-1}*\left[\begin{array}{c}n\\ 1\end{array}\right]$

WhereinIt is μ₁In a first order Taylor expansion value of zonule R, by object functionReplace withIt is right to representThe result after linear process is carried out, wherein making θ=[P^T,B^T,f^T]^TTo simplify symbol Expression, and add constraints on this basis | | θ-θ^(m-1)| |=| | Δ θ | |≤R, thus set up and complete a MILP mould Type problem ζ '₄；

S3：ζ ' is solved using YALMIP instruments₄, export θ^(m)=θ^(m-1)+ Δ θ, updates m=m+1；

S4：Whether≤R is restrained for | | the Δ θ | | that judges the condition of convergence, and S5 is jumped to if setting up, and otherwise returns to S2；

S5：Export the result of SPEB and θ.