Need a better method to plot this function
显示 更早的评论
Essentially I am trying to find a better way of solving this:
clear all; close all; home
N=10;
M=100;
L=400;
ki=zeros(L-1,floor(2/3*M+1));
kr=zeros(L-1,M-floor(2/3*M+1)+1);
options=optimset('MaxIter', 1e8, 'TolFun', 1e-6);
for i=1:floor(2/3*M+1)
p=pi/M*(i-1);
g_p=2*cos(p/2);
for j=1:L-1
g=@(x) sin(x*N)+g_p*sin(x*(N+1));
[ki(j,i), feval1, exitflag1(j,i)]=fsolve(g,pi/L*j,options);
if (exitflag1(j,i)~=1)
ki(j,i)=0;
end
end
end
for i=(floor(2/3*M+1)+1):(M+1)
m=i-floor(2/3*M+1);
p=pi/M*(i-1);
g_p=2*cos(p/2);
for j=1:L-1
g=@(x) sin(x*N)+g_p*sin(x*(N+1));
[kr(j,m), feval2, exitflag2(j,m)]=fsolve(g,pi/L*j,options);
if (exitflag2(j,m)~=1)
kr(j,m)=0;
end
end
end
which is nothing more than a set of nonlinear equations. The problem is that it takes too long for those given parameters. I know the system has the following properties: if i<=2/3*M+1 then we have N different real solutions; if i>2/3*M+1 then we have N-1 different real solutions. In both cases, the solutions have to lie between 0 and pi, without the boundary as possible solutions. The next piece of code just gets rid of the invalid solutions and sorts them in crescent order.
for j=1:floor(2/3*M+1)
ki(:,j)=sort(ki(:,j));
end
for j=(floor(2/3*M+1)+1):(M+1)
m=j-floor(2/3*M+1);
kr(:,m)=sort(kr(:,m));
end
Ki=zeros(L,floor(2/3*M+1));
Kr=zeros(L-1,M-floor(2/3*M+1)+1);
for j=1:floor(2/3*M+1)
m=0;
for i=1:L-1
if ((abs(Ki(:,j)-ki(i,j))>1e-5)&(ki(i,j)>0)&(ki(i,j)<3.1416))
m=m+1;
Ki(m,j)=ki(i,j);
end
end
end
for j=(floor(2/3*M+1)+1):(M+1)
n=j-floor(2/3*M+1);
m=0;
for i=1:L-1
if ((abs(Kr(:,n)-kr(i,n))>1e-5)&(kr(i,n)>0)&(kr(i,n)<3.1416))
m=m+1;
Kr(m,n)=kr(i,n);
end
end
end
And effectively, we get N different values for each column of Ki and N-1 different values for each column of Kr. I am seeking a faster and neater way of solving this problem making sure that for a given N (the actual problem requires N=30), it gets all possible solutions. Thanks.
Note: More information about the system I am dealing with can be found here http://iopscience.iop.org/article/10.1088/1468-6996/11/5/054504/meta, appendix B, look for the dispersion relation.
回答（1 个）
Alan Weiss
2015-11-4
1 个投票
Without understanding what you wrote in any detail, it seems to me that you are looking at a large number of real roots of a real function. If that is true, perhaps you should try chebfun.
Alan Weiss
MATLAB mathematical toolbox documentation
类别
在 帮助中心 和 File Exchange 中查找有关 Programming 的更多信息
2015-11-4
2015-11-4
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
