Understanding the error in matlab code for the given double summation

Hello all, I had written the MATLAB code for the given equation that has double summation:

where is modified Bessel function of second kind,

Code:

% Simulation parameters
a_1=0.8;
a_2=0.2;
R_1=1.8;
R_2=2.4;
u_1=(2^(R_1))-1;
u_2=(2^(R_2))-1;
Omega_0=1; Omega_1=1; Omega_2=1;
m_0=2; m_1=2; m_2=2;
rho_1 = 20;
phi=(u_1/((a_1*rho_1)-(a_2*rho_1*u_1)));

p = (m_1)-1; % Outer summation limit
q = (m_2)-1; % Inner summation limit
prod = 0; % initialize to zero

for ii = 0:p % ii is representing i in equation
    
    for jj = 0:q % jj is representing j in equation
        a=2/((factorial(ii))*(factorial(jj)));
        A_1=(((m_1*phi)/Omega_1)^ii)*(((m_2*phi)/Omega_2)^jj);
        A_2=(1/gamma(m_0))*((m_0/Omega_0)^m_0);
        A_3=(((m_1/Omega_1)+(m_2/Omega_2))*((Omega_0)*(phi/m_0)))^((m_0-ii-jj)/2);
        A_4=besselk((m_0-ii-jj),(2*(sqrt(((m_1/Omega_1)+(m_2/Omega_2))*((m_0*phi)/Omega_0)))));
        prod_1=a*A_1*A_2*A_3*A_4;           
    end
    prod = prod+prod_1;
end
op_1 = 1-prod;

My query is that why the value of op_1 is coming as 0.9711 while it should be a very small value such as 0.065.

Any help in this regard will be highly appreciated.

Answer 1

[Torsten]

The summation

prod = prod+prod_1;

must be inside the loop over jj, not outside as you coded it.