Fixed point iteration

Pages: 12
I am implementing an analytical model. In the model, I need to use fixed point iteration to find a variable "gamma". I write the following code, but the loop doesn't work. It only repeats two times even though the condition is not satisfied in the while loop.

void NormalizedOfferedLoad(double a){ //a is just use for identifying the study case,
//i.e. if a =1, uniform distribution, a=2, with SPR
double lamda=0;
double v1;
double v2; //store the result of the nCr function
double Lpay= 4000; //packet payload in bits
double Rdata = 11 * pow(10, 6); //Tranmission rate
double σ = 20 * pow(10, -6);
double Pb=0,Ps=0,P1=0,P2=0,PE=0;
double T,Tv=0, Ts, Tc, Tdata, Tsifs, Tack,Tdifs;
double R,b=0; //γ is the collision probability experencied by a tagged node
double z,x,c; // The probability of choosing the channel to transmit the packet
double S=0; //Average Service Time
double L =4576; //packet size in bits
double gamma=0.5; //Let 1 be the initial guess of γ
double gamma2=0.5;
double sum=0;
double tau=0;
double ATPS=0; //ATPS is the attempt rate per slot, i.e. τ ′;
double normalizedOL=0.1;
double ρ =1;
int count=0;
Tack = 304/(11*pow(10, 6));
Tsifs = 10 * pow(10, -6);
Tdifs = 50 * pow(10, -6);
Tdata =4.16* pow(10, -4);
//Tack = 3.04 * pow(10, -4);

while(fabs(gamma2)-fabs(gamma)<=0.001){
lamda =(normalizedOL * Rdata)/(N*Lpay);
gamma =gamma2;

count++;

if(a==1){
PE = 1- (M* pow((1/M), 2));}
else{
PE= 1- pow(z, 2) + pow(x, 2) + pow(c, 2);
}
v1 = CPascal(N, 1);
v2 = CPascal(N, 2);
for(int i=0;i<=K;i++){ //R(γ)
R += pow(gamma,i);
}

for(int j=0;j<=K;j++){
if (j==0)
b= W/2;
else if(j>=1 && j<=m-1)
b=pow(2, j)* (W/2)* pow(gamma,j);
else if(j>=m &&j<=K){
b = pow(2, m) * (W/2) * pow(gamma,j);}
sum =sum+b;


}
ATPS = R/sum;

Ts = Tdata + Tsifs + Tack + Tdifs;
Tc = Tdata+Tack+Tsifs+Tdifs;
Tv = ((1-Pb)* σ)+ (Pb*Ps*(Ts+σ))+(Pb*(1-Ps)*(Tc+σ));
S = sum *Tv;
ρ = lamda*S;
if(ρ>1)
{
ρ=1;
}
tau = ρ* ATPS;

P1= v1 * tau * pow((1-tau), N-1);
P2= v2 * pow(tau, 2) * pow((1-tau), (N-2)) * PE;
Pb = 1-(pow((1-tau),N));
Ps = (P1 + P2)/ Pb;


gamma2 = 1- (pow(1-tau, N-1) - (N-1) * tau *pow(1-tau, N-2)*PE);

}
Last edited on
There is something badly wrong with your "continue-to-iterate" condition
while(fabs(gamma2)-fabs(gamma)<=0.001)
Half-way down your loop you set gamma = gamma2 and I don't think (difficult to tell, because your code indentation is non-existent) that you don't change either one again in the loop. Thus your loop will run once.

From what single-point iteration does I presume that your condition ought to be
while(fabs(gamma2-gamma) > 0.001)
since you are presumably testing the difference between successive iterates.
However, this won't loop at all unless you set different values of gamma and gamma2 at the start. You will also have to update one of them WITHIN the loop.


Other things:
- Please use code tags: it makes it much easier for people to help you;
- Please supply complete, compileable code: it makes it much easier for people to help you;
- Please don't write horrible expressions like
Tack = 304/(11*pow(10, 6));
Just use scientific notation:
Tack = 304/11.0e6;
- Presumably N is a global variable; I can't fine it defined anywhere in what you have supplied.


So,
- check your "continue-to-iterate" condition and fix it;
- update gamma or gamma2 correctly WITHIN the loop;
- repost FULL code, with INDENTATION and CODE TAGS.
#include <stdio.h>
#include <iostream>
#include <cmath>
#include <fstream>
using namespace std;
double m=5; //determine the maximum backoff window, CWmax= pow(2,m)*CWmin
//since CWmax=1024, m equals 5
double W =32; //minimum contention window
double K=7; //K is the retransmission limit
double N=5; //number of node
double M=3; //number of non-zero power avaliable


int CPascal(int n , int r) { //The nCr function
int v = 1;
for (int i = 1; i <= r; i++)
v = v * (n + 1 - i) / i;
return v;
}




void NormalizedOfferedLoad(double a){ //a is just use for identifying the study case,
//i.e. if a =1, uniform distribution, a=2, with SPR
double lamda=0;
double v1;
double v2; //store the result of the nCr function
double Lpay= 4000; //packet payload in bits
double Rdata = 11 * pow(10, 6); //Tranmission rate
double σ = 20 * pow(10, -6);
double Pb=0,Ps=0,P1=0,P2=0,PE=0;
double T,Tv=0, Ts, Tc, Tdata, Tsifs, Tack,Tdifs;
double R,b=0; //γ is the collision probability experencied by a tagged node
double z,x,c; // The probability of choosing the channel to transmit the packet
double S=0; //Average Service Time
double L =4576; //packet size in bits
double gamma=0.2; //Let 1 be the initial guess of γ
double gamma2=0.2;
double sum=0;
double tau=0;
double ATPS=0; //ATPS is the attempt rate per slot, i.e. τ ′;
double normalizedOL=0.1;
double p ;
int count=0;
Tack = 304/(11*pow(10, 6));
Tsifs = 10 * pow(10, -6);
Tdifs = 50 * pow(10, -6);
Tdata =4.16* pow(10, -4);
//Tack = 3.04 * pow(10, -4);

while(fabs(gamma2-gamma)<0.0001){

lamda =(normalizedOL * Rdata)/(N*Lpay);
gamma =gamma2;

count++;

if(a==1){
PE = 1- (M* pow((1/M), 2));}
else{
PE= 1- pow(z, 2) + pow(x, 2) + pow(c, 2);
}
v1 = CPascal(N, 1);
v2 = CPascal(N, 2);
for(int i=0;i<=K;i++){ //R(γ)
R += pow(gamma,i);
}

for(int j=0;j<=K;j++){
if (j==0)
b= W/2;
else if(j>=1 && j<=m-1)
b=pow(2, j)* (W/2)* pow(gamma,j);
else if(j>=m &&j<=K){
b = pow(2, m) * (W/2) * pow(gamma,j);}
sum =sum+b;


}
ATPS = R/sum;

Ts = Tdata + Tsifs + Tack + Tdifs;
Tc = Tdata+Tack+Tsifs+Tdifs;
Tv = ((1-Pb)* σ)+ (Pb*Ps*(Ts+σ))+(Pb*(1-Ps)*(Tc+σ));
S = sum *Tv;
p = lamda*S;
if(p>1)
{
p=1;
}

tau = p* ATPS;


P1= v1 * tau * pow((1-tau), N-1);
P2= v2 * pow(tau, 2) * pow((1-tau), (N-2)) * PE;
Pb = 1-(pow((1-tau),N));
Ps = (P1 + P2)/ Pb;

gamma2 = 1- (pow(1-tau, N-1) - (N-1) * tau *pow(1-tau, N-2)*PE);




}

/*if((γ2-γ)/γ2<=0.0001){
goto End;
}

else{
γ =γ2;
goto Loop;}*/
cout<<"Number of iteration: "<<count<<endl;
cout<<"gamma2: "<<gamma2<<endl;
cout<<"gamma: "<<gamma<<endl;
cout<<"Tv: "<<Tv<<endl;
cout<<"P1: "<<P1<<endl;
cout<<"P2: "<<P2<<endl;
cout<<"Pb: "<<Pb<<endl;
cout<<"Ps: "<<Ps<<endl;
cout<<"S: "<<S<<endl;
cout<<"p: "<<p<<endl;
cout<<"R: "<<R<<endl;
cout<<"Sum: "<<sum<<endl;
cout<<"Attempt rate per slot: "<<ATPS<<endl;
cout<<"τ: "<<tau<<endl;


T = ((L*P1) + (2*L*P2))/Tv;
cout<<"System throughput: "<<T<<" bits per seond"<<endl;


//λ =(normalizedOL * Rdata)/(N*Lpay);
cout<<"Normalized Offered Load: "<<normalizedOL<<endl;
cout<<"λ: "<<lamda <<endl;


}




int main(int argc, const char * argv[]) {



NormalizedOfferedLoad(1);


return 0;
}
Here is the full code. Please check it.
Please go back and edit it. Put it in code tags (the <> in the format menu to the right) so that other people can read it and can see the loop structure from the indentation.

It is still using an incorrect condition for iterating (why would you want to keep iterating once the gamma values are close to each other?) and you are still failing to update gamma2 within the loop.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
#include <stdio.h>
#include <iostream>
#include <cmath>
#include <fstream>
using namespace std;
double m=5; //determine the maximum backoff window, CWmax= pow(2,m)*CWmin
//since CWmax=1024, m equals 5
double W =32; //minimum contention window
double K=7; //K is the retransmission limit
double N=5; //number of node
double M=3; //number of non-zero power avaliable


int CPascal(int n , int r) { //The nCr function
int v = 1;
for (int i = 1; i <= r; i++)
v = v * (n + 1 - i) / i;
return v;
}




void NormalizedOfferedLoad(double a){ //a is just use for identifying the study case,
//i.e. if a =1, uniform distribution, a=2, with SPR
double lamda=0;
double v1;
double v2; //store the result of the nCr function
double Lpay= 4000; //packet payload in bits
double Rdata = 11 * pow(10, 6); //Tranmission rate
double σ = 20 * pow(10, -6);
double Pb=0,Ps=0,P1=0,P2=0,PE=0;
double T,Tv=0, Ts, Tc, Tdata, Tsifs, Tack,Tdifs;
double R,b=0; //γ is the collision probability experencied by a tagged node
double z,x,c; // The probability of choosing the channel to transmit the packet
double S=0; //Average Service Time
double L =4576; //packet size in bits
double gamma=0.2; //Let 1 be the initial guess of γ
double gamma2=0.2;
double sum=0;
double tau=0;
double ATPS=0; //ATPS is the attempt rate per slot, i.e. τ ′;
double normalizedOL=0.1;
double p ;
int count=0;
Tack = 304/(11*pow(10, 6));
Tsifs = 10 * pow(10, -6);
Tdifs = 50 * pow(10, -6);
Tdata =4.16* pow(10, -4);
//Tack = 3.04 * pow(10, -4);

while(fabs(gamma2-gamma)<0.0001){

lamda =(normalizedOL * Rdata)/(N*Lpay);
gamma =gamma2;

count++;

if(a==1){
PE = 1- (M* pow((1/M), 2));}
else{
PE= 1- pow(z, 2) + pow(x, 2) + pow(c, 2);
}
v1 = CPascal(N, 1);
v2 = CPascal(N, 2);
for(int i=0;i<=K;i++){ //R(γ)
R += pow(gamma,i);
}

for(int j=0;j<=K;j++){
if (j==0)
b= W/2;
else if(j>=1 && j<=m-1)
b=pow(2, j)* (W/2)* pow(gamma,j);
else if(j>=m &&j<=K){
b = pow(2, m) * (W/2) * pow(gamma,j);}
sum =sum+b;


}
ATPS = R/sum;

Ts = Tdata + Tsifs + Tack + Tdifs;
Tc = Tdata+Tack+Tsifs+Tdifs;
Tv = ((1-Pb)* σ)+ (Pb*Ps*(Ts+σ))+(Pb*(1-Ps)*(Tc+σ));
S = sum *Tv;
p = lamda*S;
if(p>1)
{
p=1;
}

tau = p* ATPS;


P1= v1 * tau * pow((1-tau), N-1);
P2= v2 * pow(tau, 2) * pow((1-tau), (N-2)) * PE;
Pb = 1-(pow((1-tau),N));
Ps = (P1 + P2)/ Pb;

gamma2 = 1- (pow(1-tau, N-1) - (N-1) * tau *pow(1-tau, N-2)*PE);




}

cout<<"Number of iteration: "<<count<<endl;
cout<<"gamma2: "<<gamma2<<endl;
cout<<"gamma: "<<gamma<<endl;
cout<<"Tv: "<<Tv<<endl;
cout<<"P1: "<<P1<<endl;
cout<<"P2: "<<P2<<endl;
cout<<"Pb: "<<Pb<<endl;
cout<<"Ps: "<<Ps<<endl;
cout<<"S: "<<S<<endl;
cout<<"p: "<<p<<endl;
cout<<"R: "<<R<<endl;
cout<<"Sum: "<<sum<<endl;
cout<<"Attempt rate per slot: "<<ATPS<<endl;
cout<<"τ: "<<tau<<endl;


T = ((L*P1) + (2*L*P2))/Tv;
cout<<"System throughput: "<<T<<" bits per seond"<<endl;


//λ =(normalizedOL * Rdata)/(N*Lpay);
cout<<"Normalized Offered Load: "<<normalizedOL<<endl;
cout<<"λ: "<<lamda <<endl;


}




int main(int argc, const char * argv[]) {



NormalizedOfferedLoad(1);


return 0;
}
Since I'm trying to do a fixed point iteration. I just set the tolerance that if the new gamma is close to the previous gamma which is not exceeding the tolerance the new gamma is the final gamma. How should I update the gamma2 so that I won't end it after the first time?
Your code. With indentation. And gamma and gamma2 sufficiently different at the start that the loop will run.

It produces something - for you to check. However, your CPascal() function will fail due to integer division. You still have that to fix.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
#include <stdio.h>
#include <iostream>
#include <cmath>
#include <fstream>
using namespace std;
double m=5; //determine the maximum backoff window, CWmax= pow(2,m)*CWmin
//since CWmax=1024, m equals 5
double W =32; //minimum contention window
double K=7; //K is the retransmission limit
double N=5; //number of node
double M=3; //number of non-zero power avaliable


int CPascal(int n , int r) 
{
   int v = 1;
   for (int i = 1; i <= r; i++)
   v = v * (n + 1 - i) / i;
   return v;
}




void NormalizedOfferedLoad(double a)
{
   //i.e. if a =1, uniform distribution, a=2, with SPR
   double lamda=0;
   double v1;
   double v2; //store the result of the nCr function
   double Lpay= 4000; //packet payload in bits
   double Rdata = 11 * pow(10, 6); //Tranmission rate
   double s = 20 * pow(10, -6);
   double Pb=0,Ps=0,P1=0,P2=0,PE=0;
   double T,Tv=0, Ts, Tc, Tdata, Tsifs, Tack,Tdifs;
   double R,b=0; //? is the collision probability experencied by a tagged node
   double z,x,c; // The probability of choosing the channel to transmit the packet
   double S=0; //Average Service Time
   double L =4576; //packet size in bits
   double gamma=0.2; //Let 1 be the initial guess of ?
   double gamma2=0.2 + 1;
   double sum=0;
   double tau=0;
   double ATPS=0; //ATPS is the attempt rate per slot, i.e. t ';
   double normalizedOL=0.1;
   double p ;
   int count=0;

   Tack = 304/11e6;
   Tsifs = 10e-6;
   Tdifs = 50e-6;
   Tdata =4.16e-4;
   
   while ( fabs( gamma2-gamma ) > 0.0001 )
   {
      lamda =( normalizedOL * Rdata ) / ( N * Lpay );
      gamma =gamma2;
   
      count++;
   
      if ( a==1 ) PE = 1- ( M * pow( ( 1 / M ), 2 ) );
      else        PE= 1 - z * z + x * x + c * c;
      v1 = CPascal(N, 1);
      v2 = CPascal(N, 2);
      for ( int i = 0; i <= K; i++ ) R += pow( gamma, i );
   
      for ( int j = 0; j <= K; j++ )
      {
         if      ( j == 0               )  b = W / 2;
         else if ( j >= 1 && j <= m - 1 )  b = pow( 2, j ) * ( W / 2 ) * pow( gamma, j );
         else if ( j >= m && j <= K     )  b = pow( 2, m ) * ( W / 2 ) * pow( gamma, j );
         sum += b;
      }
      ATPS = R/sum;
   
      Ts = Tdata + Tsifs + Tack + Tdifs;
      Tc = Tdata+Tack+Tsifs+Tdifs;
      Tv = ((1-Pb)* s)+ (Pb*Ps*(Ts+s))+(Pb*(1-Ps)*(Tc+s));
      S = sum *Tv;
      p = lamda*S;
      if ( p > 1 ) p = 1;
      
      tau = p * ATPS;
      
      P1= v1 * tau * pow( ( 1 - tau ), N - 1 );
      P2= v2 * pow( tau, 2 ) * pow( ( 1 - tau ), ( N - 2 ) ) * PE;
      Pb = 1 - ( pow( ( 1 - tau ), N ) );
      Ps = ( P1 + P2 )/ Pb;
      
      gamma2 = 1 - ( pow( 1 - tau, N - 1) - ( N - 1 ) * tau * pow( 1 - tau, N - 2 ) * PE );
   }
   
   cout << "Number of iteration: " << count << endl;
   cout << "gamma2: " << gamma2 << endl;
   cout << "gamma: " << gamma << endl;
   cout << "Tv: " << Tv << endl;
   cout << "P1: " << P1 << endl;
   cout << "P2: " << P2 << endl;
   cout << "Pb: " << Pb << endl;
   cout << "Ps: " << Ps << endl;
   cout << "S:  " << S  << endl;
   cout << "p:  " << p  << endl;
   cout << "R:  " << R  << endl;
   cout << "Sum:" << sum<<endl;
   cout << "Attempt rate per slot: " << ATPS << endl;
   cout << "t: " << tau << endl;
   
   
   T = ( ( L * P1 ) + ( 2 * L * P2 ) ) / Tv;
   cout << "System throughput: " << T << " bits per seond" << endl;
   
   
   //? =(normalizedOL * Rdata)/(N*Lpay);
   cout << "Normalized Offered Load: " << normalizedOL << endl;
   cout << "?: " << lamda << endl;
}
   
   
   
   
int main(int argc, const char * argv[])
{
   NormalizedOfferedLoad(1);
   return 0;
}

Last edited on
thank you so much!
But now I have another question, since the variable sum accumulate every time. But I want it to accumulate. I add the code "sum=0;" at the beginning of the loop. However, the program goes wrong and cannot show the result correctly.
I don't know what your set of equations are, so I'm offering suggestions blind. I assume sum is created from the loop in lines 67-73. You need to set sum=0 just before that, not at the start of your program.

On another issue, make sure that you go back and fix your CPascal function. Test it separately - it will fail due to integer division.
since base on the program you modify for me, if gamma final approximately equals to 0.21, sum should not be so large.
ρ = λS
τ = min(1, ρ) τ ′( τ ′= ATPS in the code)
R(γ) = 1 + γ + · · · + γK
bi=W/2 for i=0,
bi= 2ib0 for 1<i<m-1
bi=2mb0 for m<i<K
W =b0 +γb1 +···+γKbK
τ ′= R(γ)/W
P[Ei ̸=Ej]=1−M* (1/M)2
Pb =1−(1−τ)N
Ps=(P1+P2)/Pb
P1= nC1 τ(1−τ)N−1
P2 = nC2 τ2(1−τ)N-2*P[Ei ̸=Ej]
Tv =(1−Pb)σ+PbPs(Ts +σ)+Pb(1−Ps)(Tc +σ)
Ts =Tdata +TSIFS +TACK +TDIFS
Tc = Tdata + Ttimeout + TDIFS
S=WTv
γ=1−(1−τ)N−1−(N−1)τ(1−τ)N−2P[Ei ̸= Ej]
T= (LP1+2LP2)/Tv
here is the equations , hope it won't bother you too much to take a look.
Have you made the changes that I recommended?

- fix the CPascal function - at present it will fail due to integer division.
- set sum=0 just before line 67 above, NOT at the start of the program - otherwise you will just accumulate size every pass through the iteration.

Then post this code.

Also, cut down the unnecessary use of the pow() function.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
//
//  main.cpp
//  FYP4
//
//  Created by Lau Siu Chung on 9/11/2016.
//  Copyright © 2016年 Lau Siu Chung. All rights reserved.
//

#include <stdio.h>
#include <iostream>
#include <cmath>
#include <fstream>
using namespace std;
double m=5; //determine the maximum backoff window, CWmax= pow(2,m)*CWmin
//since CWmax=1024, m equals 5
double W =32; //minimum contention window
double K=7; //K is the retransmission limit
double N=5; //number of node
double M=3; //number of non-zero power avaliable


/*int CPascal(int n , int r)
{
    int v = 1;
    for (int i = 1; i <= r; i++)
        v = v * (n + 1 - i) / i;
    return v;
}*/




void NormalizedOfferedLoad(double a)
{
    //i.e. if a =1, uniform distribution, a=2, with SPR
    double lamda=0;
    double Lpay= 4000; //packet payload in bits
    double Rdata = 11e6; //Tranmission rate
    double s = 20e-6;
    double Pb=0,Ps=0,P1=0,P2=0,PE=0;
    double T,Tv=0, Ts, Tc, Tdata, Tsifs, Tack,Tdifs;
    double R,b=0; //? is the collision probability experencied by a tagged node
    double z,x,c; // The probability of choosing the channel to transmit the packet
    double S=0; //Average Service Time
    double L =4576; //packet size in bits
    double gamma=0.2; //Let 1 be the initial guess of ?
    double gamma2=0.2 +0.1;
    //double sum=0;
    double tau=0;
    double ATPS=0; //ATPS is the attempt rate per slot, i.e. t ';
    double normalizedOL=0.5;
    double p ;
    int count=0;
    
    Tack = 304/11e6;
    Tsifs = 10e-6;
    Tdifs = 50e-6;
    Tdata =4.16e-4;
    
    while ( fabs( gamma2-gamma ) > 0.0001 )
    {
        
        lamda =( normalizedOL * Rdata ) / ( N * Lpay );
        gamma =gamma2;
        
        count++;
        
        if ( a==1 )
            PE = 1- ( M * pow( ( 1 / M ), 2 ) );
        else
            PE= 1 - z * z + x * x + c * c;
        //v1 = CPascal(N, 1);
        //v2 = CPascal(N, 2);
        for ( int i = 0; i <= K; i++ )
            R += pow( gamma, i );
        
        double sum =0;
        for ( int j = 0; j <= K; j++ )
        {
            if  ( j == 0)
                b = W / 2;
            else if ( j >= 1 && j <= m - 1 )
                b = pow( 2, j ) * ( W / 2 ) * pow( gamma, j );
            else if ( j >= m && j <= K     )
                b = pow( 2, m ) * ( W / 2 ) * pow( gamma, j );
            sum += b;
        }
        
        ATPS = R/sum;
        Ts = Tdata + Tsifs + Tack + Tdifs;
        Tc = Tdata+Tack+Tsifs+Tdifs;
        Tv = ((1-Pb)* s)+ (Pb*Ps*(Ts+s))+(Pb*(1-Ps)*(Tc+s));
        S = sum *Tv;
        p = lamda*S;
        if ( p > 1 )
            p = 1;
        
        tau = p * ATPS;
        
        P1= 5 * tau * pow( ( 1 - tau ), N - 1 );
        P2= 10 * pow( tau, 2 ) * pow( ( 1 - tau ), ( N - 2 ) ) * PE;
        Pb = 1 - ( pow( ( 1 - tau ), N ) );
        Ps = ( P1 + P2 )/ Pb;
        
        gamma2 = 1 - ( pow( 1 - tau, N - 1) - ( N - 1 ) * tau * pow( 1 - tau, N - 2 ) * PE );
        
    }
    
    
    cout << "Number of iteration: " << count << endl;
    cout << "gamma2: " << gamma2 << endl;
    cout << "gamma: " << gamma << endl;
    cout << "Tv: " << Tv << endl;
    cout << "P1: " << P1 << endl;
    cout << "P2: " << P2 << endl;
    cout << "Pb: " << Pb << endl;
    cout << "Ps: " << Ps << endl;
    cout << "S:  " << S  << endl;
    cout << "p:  " << p  << endl;
    cout << "R:  " << R  << endl;
    //cout << "Sum:" << sum<<endl;
    cout << "Attempt rate per slot: " << ATPS << endl;
    cout << "t: " << tau << endl;
    
    
    T = ( ( L * P1 ) + ( 2 * L * P2 ) ) / Tv;
    cout << "System throughput: " << T << " bits per seond" << endl;
    
    
    //? =(normalizedOL * Rdata)/(N*Lpay);
    cout << "Normalized Offered Load: " << normalizedOL << endl;
    cout << "Lamda: " << lamda << endl;
}




int main(int argc, const char * argv[])
{
    NormalizedOfferedLoad(1);
    return 0;
}
I delete the CPasal function as the value won't change in the program. I just replace it by directly multiply 5 and 10.
I add the code:double sum =0; in line 77, but the result is wrong, it show:
Number of iteration: 20
gamma2: nan
gamma: -9.12792e+54
Tv: -4.98015e+64
P1: nan
P2: nan
Pb: nan
Ps: nan
S:  nan
p:  nan
R:  nan
Attempt rate per slot: nan
t: nan
System throughput: nan bits per seond
Normalized Offered Load: 0.5
Lamda: 275
Program ended with exit code: 0
I also doubt if the order of the equations matter?
Try this (and READ MY NOTES AT THE BOTTOM)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
//
//  main.cpp
//  FYP4
//
//  Created by Lau Siu Chung on 9/11/2016.
//  Copyright © 2016? Lau Siu Chung. All rights reserved.
//

#include <stdio.h>
#include <iostream>
#include <cmath>
#include <fstream>
using namespace std;
double m = 5;   // determine the maximum backoff window, CWmax= pow(2,m)*CWmin
                // since CWmax=1024, m equals 5
double W = 32;  //minimum contention window
double K = 7;   //K is the retransmission limit
double N = 5;   //number of node
double M = 3;   //number of non-zero power avaliable


/*int CPascal(int n , int r)          // *** FIX THIS **** wrong due to integer division
{
    int v = 1;
    for (int i = 1; i <= r; i++)
        v = v * (n + 1 - i) / i;
    return v;
}*/




void NormalizedOfferedLoad(double a)
{
    //i.e. if a = 1, uniform distribution, a = 2, with SPR
    double lamda =0;
    double Lpay = 4000;           // packet payload in bits
    double Rdata = 11.0e6;        // transmission rate
    double sigma = 20.0e-6;
    double Pb = 0, Ps = 0, P1 = 0, P2 = 0, PE = 0;
    double T, Tv = 0, Ts, Tc, Tdata, Tsifs, Tack, Tdifs;
    double R; 
    double z, x, c;               // the probability of choosing the channel to transmit the packet
    double S;                     // Average Service Time
    double L = 4576;              // packet size in bits
    double ATPS = 0;              // ATPS is the attempt rate per slot, i.e. t ';
    double normalizedOL = 0.5;
    double rho;
    double tau;                                                                                                                      
    double sum;
    int count = 0;
    double gamma = 0.2;           // initial guess
    double gamma2 = 0.2 + 0.1;    // anything different from gamma

    Tack  = 304 / 11.0e6;                                          // did you mean 304 / Rdata?
    Tsifs = 10.0e-6;
    Tdifs = 50.0e-6;
    Tdata = 4.16e-4;
    
                                                                   // Moved the following out of the loop, as they are constant
    lamda = ( normalizedOL * Rdata ) / ( N * Lpay );               // 
    Ts = Tdata + Tsifs + Tack + Tdifs;                             //           These seem to be the same?
    Tc = Tdata + Tack + Tsifs + Tdifs;                             //           *** DOESN'T CORRESPOND TO EQUATIONS
    if ( a == 1 ) PE = 1 - 1 / M;                                  //           Simplified equation
    else          PE = 1 - z * z + x * x + c * c;                  //           CHECK THIS - NOT TESTED HERE


    while ( fabs( gamma2 - gamma ) > 0.0001 )
    {
        gamma = gamma2;
        count++;
        
        //v1 = CPascal(N, 1);                                      // *** FIX THIS ***
        //v2 = CPascal(N, 2);                                      // *** FIX THIS ***

        R = 0;                                                     // *** CORRECTION *** added initialiser
        for ( int i = 0; i <= K; i++ ) R += pow( gamma, i );
        
        sum = 1;
        for ( int j = 1; j <= m - 1; j++ ) sum += pow( 2.0 * gamma, j );             // simplified and CORRECTED sums *****
        for ( int j = m; j <= K    ; j++ ) sum += pow( 2.0, m ) * pow( gamma, j );    
        sum = sum * W / 2;

        ATPS = R / sum;
        Tv = ( 1 - Pb ) * sigma + Pb * Ps * ( Ts + sigma ) + Pb * ( 1 - Ps ) * ( Tc + sigma );      // simplified

        S = sum * Tv;
        rho = lamda * S;   if ( rho > 1 ) rho = 1;
        tau = rho * ATPS;
        
        P1 = 5 * tau * pow( ( 1 - tau ), N - 1 );                       // *** FIX nCr
        P2 = 10 * pow( tau, 2 ) * pow( ( 1 - tau ), ( N - 2 ) ) * PE;   // *** FIX nCr
        Pb = 1 - ( pow( ( 1 - tau ), N ) );
        Ps = ( P1 + P2 ) / Pb;
        
        gamma2 = 1 - pow( 1 - tau, N - 1) - ( N - 1 ) * tau * pow( 1 - tau, N - 2 ) * PE;    // *** CORRECTED TO MATCH EQUATIONS ***
        
    }
    
    
    cout << "Number of iteration: " << count << endl;
    cout << "gamma2: " << gamma2 << endl;
    cout << "gamma: " << gamma << endl;
    cout << "Tv: " << Tv << endl;
    cout << "P1: " << P1 << endl;
    cout << "P2: " << P2 << endl;
    cout << "Pb: " << Pb << endl;
    cout << "Ps: " << Ps << endl;
    cout << "S:  " << S  << endl;
    cout << "rho:" << rho<< endl;
    cout << "R:  " << R  << endl;
    cout << "Sum:" << sum<< endl;
    cout << "Attempt rate per slot: " << ATPS << endl;
    cout << "t: " << tau << endl;
        
    T = ( ( L * P1 ) + ( 2 * L * P2 ) ) / Tv;
    cout << "System throughput: " << T << " bits per seond" << endl;
    
    
    //? =(normalizedOL * Rdata)/(N*Lpay);
    cout << "Normalized Offered Load: " << normalizedOL << endl;
    cout << "Lamda: " << lamda << endl;
}


int main(int argc, const char * argv[])
{
    NormalizedOfferedLoad(1);
    return 0;
}



- Your code for Tc doesn't match the equations you gave
- I have corrected a number of things (LOOK AT THE CODE) - most important was probably initialising the sum for R within the loop
- DO NOT USE GREEK SYMBOLS IN YOUR CODE (write out in full: e.g. rho, sigma)
- if you want some code for the nCr combination then you can cut and paste the nCr() function from below:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
#include<iostream>
#include<cmath>
using namespace std;

int nCr( int n, int r )
{
   double value = 1.0;
   for ( int i = 0; i < r; i++ ) value = value * ( n - i ) / ( r - i );
   return floor( value + 0.5 );
}


int main()
{
   int n, r;
   cout << "Input n and r: ";   cin >> n >> r;
   cout << "nCr = " << nCr( n, r );
}

You can then correct your code.

Last edited on
closed account (48T7M4Gy)
j <= K where j is int and K is double can lead to spurious results due to rounding/precision problems with the inequality relation.

This (OP) code should be written in such a way that it can be unit-tested by appropriate debugging and/or variable trace methods etc.

Having just one function as is the case here makes the testing task well nigh impossible.

Along with all the admirable advice from lastchance the program should be written in a well thought out and structured way instead of a vast continuous list of calculations. Test data should be documented and used for each stage.
I thought I was stretched just to make the program run! Signal processing isn't my field!

I think you can safely (and appropriately) make m, W, K and N into ints rather than doubles, but check for any integer divisions. You should eventually make M an int, but you will have to prevent integer division in line 64.

As kemort has suggested, develop your code slowly (and write out variable values, or use a debugger, to check intermediate results), consider breaking down into subfunctions (especially to do those summations) and test it scrupulously.

Individual functions should be tested separately.
Pages: 12