I don't get a reasonable output for my code?

AI Thread Summary
The user is attempting to plot the probability of error using a complex equation in Matlab but is encountering incorrect y-axis values, which should range between 0 and 1.2 but instead fall between 0.492 and 0.5. The code utilizes the "trapz" function for numerical integration, but the user is advised to validate intermediate results and debug the code by breaking it into smaller parts. There is a specific concern regarding the definition of the variable 'z', which is set to a range that may not align with the intended calculations. The discussion emphasizes the importance of checking variable dimensions and ensuring that all components of the code are functioning correctly.
deema_master
Messages
13
Reaction score
0
I'm trying to plot the probability of error for the following equation using Matlab software, i want to use the command "trapz" for the numerical integration, the problem is that i get a fine shape for the plot, but the values in the y-axis are wrong, the whole curve should be between 0 and 1.2 but it is between 0.492 and 0.5! Can anyone just tell me what is wrong in my code, or just give me a hint? I really need help. Here is my formula that i need to plot written using Maketex:$$P(e)=0.5-\frac{.5}{\sqrt{\pi}}\sum_{\alpha=1-k_1/2}^{\infty}\sum_{\beta=1-k_2/2}^{\infty}\sum_{\eta=0}^{\infty}\sum_{N=0}^{\infty}\sum_{M=0}^{\infty}\sum_{Q=0}^{\infty}\sum_{i=0}^{\eta}\sum_{j=0}^{N}\sum_{A=0}^{N-j}v^{\eta+N+1/2} C
\int_0^{\infty}\exp(-z*v(1+1.5/v))z^{\eta+N-.5} \frac{1}{(z+1)^{A+Q+Nrkr/2} (z+.5)^{M+i+j+Nsks/2}} dz$$

where:
$$C=0.25*\exp \left(-\frac{\lambda_1}{2}\right)*\left(\frac{\lambda_1^2}{4}\right)^{\alpha/2}*\left(\frac{\lambda_2^2}{4}\right)^{\beta/2}*\exp \left(-\frac{\lambda_2}{2}\right)*\left(\frac{\lambda_1}{4 *em *v}\right)^{\eta}*\left(\frac{\lambda_2}{4 *em* v}\right)^N*\exp (-\frac{\lambda_r*Nr}{2})*0.25^{Ns*ks/4-0.5}*\exp (-\frac{\lambda_s*Ns}{2})*0.25^{Nr*kr/4-0.5}*(Ns*\lambda_s/4)^M*(Nr*\lambda_r/4)^Q*\frac{1}{\eta! M! N! Q! \Gamma (M+Ks*Ns/2) \Gamma (Q+Nr*kr/2) \Gamma (\beta+N+1) \Gamma (\eta+\alpha+1)}*{em}^A {\eta\choose{i}}{N\choose{j}}{{N-j}\choose{A}}(em+1)^{N-j-A} \Gamma (A+Q+Nr*kr/2) \Gamma (i+j+M+Ns*ks/2)*\left(2^{A+Q+Nr*kr/2)}\right.$$And here is my Matlab code:
Code:
close all; clear;clc;
Nr=2;Ns=2;
lmda1=.3; lmda2=.3;
lmdas=.1; lmdar=.1;

z= 0.0001:1:40;
k1=2;k2=2;
kr=2.*Nr;ks=2.*Ns;
ax=0;
avg=0.0001:1:40;
em=1;
ch=2;

for alp=1-k1.*.5:ch
    for beta=1-k2.*.5:ch
        for eta=0:ch
            for N=0:ch
                for M=0:ch
                    for Q=0:ch
                        for id=0:eta
                            for jd=0:N
                                for A=0:N-jd
                                    %
                                  
                                  
                                  
                                  
                                    up=.25.*exp(-lmda1./2).*(lmda1./2).^(alp).*(lmda2.^2./(4)).^(beta./2).*exp(-lmda2./2).*(lmda1./(4.*em.*avg)).^eta.*(lmda2./(4.*em.*avg)).^N.*exp(-lmdas.*Ns.*.5).*.25.^(ks.*.25-.5).*exp(-lmdar.*Nr.*.5).*.25.^(kr.*.25-.5).*(Ns.*lmdas.*.25).^M.*(Nr.*lmdar.*.25).^Q;
                                  
                                    cy=up.*(1./(factorial(eta).*factorial(N).*factorial(M).*factorial(Q).*gamma(eta+alp+1).*gamma(N+beta+1).*gamma(M+ks.*.5).*gamma(Q+kr.*.5)));
                                    cj=cy.*(factorial(eta)./(factorial(id).*factorial(eta-id))).*(factorial(N)./(factorial(jd).*factorial(N-jd))).*gamma(M+id+jd+ks.*.5);
                                    f1=(cj.*(factorial(N-jd)./(factorial(A).*factorial(N-jd-A))).*em.^A.*(((em+1).^(N-jd-A))).*gamma(kr.*.5+Q+A));
                                  
                                    f2=f1.*(2.^(kr.*.5+Q+A)).*avg.^(eta+N);
                                    ax=ax+f2;
                                  
                                  
                                end
                            end
                        end
                    end
                end
            end
        end
    end
end
q2=2;n2=2;N2=1;eta2=1;

fun2 = exp(-z.*avg.*(1+1.5./avg)).*z.^(eta2+N2-1./2).*(1./((1+z).^(q2).*(1./2+z).^(n2)));

out= trapz(z,fun2);

b=.5.*(1-ax.*(1./sqrt(pi)).*out.*avg.^(1./2));
plot(avg,b);grid;
 
Physics news on Phys.org
I fixed your latex as you needed $$ to bracket the expressions
 
jedishrfu said:
I fixed your latex as you needed $$ to bracket the expressions
great! Thank you.
 
This is probably something we can't help you with directly. My suggestion is to pick apart your program and validate that intermediate results are coming out correctly.

Programmers will often break up code into smaller pieces to make it easier to debug and would add print statements to show intermediate results.
 
What is going on with z (and avg)?

You have:
z=0.0001:1:40

Is this what you intend?
z = 0.0001 1.0001 2.0001 3.0001 ...etc... 39.0001 (length of z is 40)
 
Ye
lewando said:
What is going on with z (and avg)?

You have:
z=0.0001:1:40

Is this what you intend?
z = 0.0001 1.0001 2.0001 3.0001 ...etc... 39.0001 (length of z is 40)
Yes this is what i want...i know z should go to a larger value...but i get a dimension error so i made their size the same.
 

Similar threads

Back
Top