How to code following equation?
7 views (last 30 days)
Show older comments
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1108960/image.png)
My code for the above equation is given below.
But looks like it don't work well. Could anyone have a look at whether I implemented the above formula well?
clear all;clc;close all;
syms j l k p q s K
m = 3;
z=linspace(0.0001,5000);
w0=0.05;
Omega = 30;
k_w = 4.05*10^6;
gama=(1./(w0.^2))+((1i.*k_w)./(2.*z));
gama_star = subs(gama, 1i, -1i);
con = 1./(16.*pi.*pi.*z.*z);
total = 0;
for j = 0:m
J = nchoosek(m,j).*((1i).^(m-j));
sum = 0;
for l = 0:m
L = nchoosek(m,l).*((1i).^(m-l));
for k = 0:1/2:l/2
KK = (((-1).^k).*factorial(l))./(gamma(k+1).*factorial(l-(2.*k)));
for q= 0:(l-2*k)
Q = nchoosek((l-(2.*k)),q);
for p= 0:1/2:(m-l)/2
P = (((2.*1i)./(sqrt(gama))).^(m-(2.*p)-(2.*k)));
for s = 0:(m-l-2*p)
a = (hermiteH(m-j+s,((Omega.*1i)./sqrt(gama)))).*(hermiteH(j+q,((Omega.*1i)./sqrt(gama))));
b = (exp(Omega).*((Omega).^(m-s-q-(2.*p)-(2.*k)))) - (exp(-Omega).*((-Omega).^(m-s-q-(2.*p)-(2.*k))));
S = nchoosek((m-l-(2.*p)),s).*a.*b.*(1./((2.*1i.*sqrt(gama_star)).^(m+s+q)));
end
end
end
end
end
sum = sum + S.*P.*Q.*KK.*L;
end
total = sum.*con;
F = vpa(real(total))
4 Comments
Answers (1)
Image Analyst
on 11 Sep 2022
Many times with complicated equations it's difficult to code and often it's because the parentheses were misplaced. I recommend you break it up into smaller terms, like term1, term2, term3, sum1, sum2, sum3, etc., and then combine them.
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1122080/image.jpeg)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1122085/image.jpeg)
5 Comments
Image Analyst
on 11 Sep 2022
That's because you have not defined z:
Unrecognized function or variable 'z'.
Error in test8 (line 2)
c=1./(16.*pi.^2.*z.^2)
Torsten
on 11 Sep 2022
m = 3;
z=linspace(0.0001,5000);
w0=0.05;
Omega = 30;
k_w = 4.05*10^6;
Gamma=(1./(w0.^2))+((1i.*k_w)./(2.*z));
Gamma_star = (1./(w0.^2))+((-1i.*k_w)./(2.*z));
F = 0;
c=1./(16.*pi.^2.*z.^2);
for l = 0:m
L = ((1i).^(m-l)).*nchoosek(m,l);
sum_j = 0;
for j = 0:m
J= ((1i).^(m-j)).*nchoosek(m,j);
sum1 = 0.0;
for k = 0:1/2:l/2
faktor1_K = (((-1).^k).*factorial(l))./(gamma(k+1).*factorial(l-(2.*k)));
for q = 0:l-2*k
faktor1_Q = nchoosek((l-(2.*k)),q);
for p = 0:1/2:(m-l)/2
faktor1_P = (((2.*1i)./(sqrt(Gamma))).^(m-(2.*p)-(2.*k)));
sum_inner = 0;
for s = 0:m-l-2*p
E1 = ((exp(Omega).*(Omega.^(m-s-q-2.*p-2.*k))) - (exp(-Omega).*(-Omega.^(m-s-q-2.*p-2.*k)))).*hermiteH(j+q,((1i.*Omega)./(sqrt(Gamma)))).*hermiteH(m-j+s,((1i.*Omega)./(sqrt(Gamma))));
sum_inner = sum_inner + nchoosek((m-l-(2.*p)),s).*(1./(((2.*1i.*sqrt(Gamma_star))).^(m+q+s))).*E1;
end
end
end
sum1 = sum1 + faktor1_K.*faktor1_P.*faktor1_Q.*sum_inner;
end
sum_j= sum_j + J.*(sum1);
end
F = F + (L.*sum_j);
end
F0 = c.*F
See Also
Categories
Find more on Calculus in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!