How to code following equation?
Show older comments
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
Torsten
on 11 Sep 2022
But looks like it don't work well.
Did you check using pencil and paper or how do you know ?
Athira T Das
on 11 Sep 2022
Torsten
on 11 Sep 2022
"not good enough" means "wrong" ?
Athira T Das
on 11 Sep 2022
Answers (1)
Image Analyst
on 11 Sep 2022
0 votes
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.
5 Comments
Athira T Das
on 11 Sep 2022
Image Analyst
on 11 Sep 2022
Start with the innermost look and compute the sum for that, like
sumOverS = 0
for s = 1 : 0:(m-l-2*p)
sumOverS = whatever.....
end
Then after the s loop but before the end of the p loop, use the sumOverS that you just got from the s loop.
Athira T Das
on 11 Sep 2022
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)
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
Categories
Find more on Matrix Indexing in Help Center and File Exchange
on 26 Aug 2022
on 11 Sep 2022
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
