0

我正在尝试使用流函数绘制流线：

psi = Y.*(F.*M.*cosh((1 + phi.*cos(2.*pi.*X)).*M)+ (1 + F.*M.^2*beta).*sinh((1 + phi.*cos(2.*pi.*X)).*M))-(F + (1 + phi.*cos(2.*pi.*X))).*sinh(M.*Y)/(1 + phi.*cos(2.*pi.*X)).*M.*cosh((1 + phi.*cos(2.*pi.*X)).*M)+ (-1 + (1 + phi.*cos(2.*pi.*X)).*M.^2.*beta).*sinh((1 + phi.*cos(2.*pi.*X)).*M )+ (alpha.*((F + (1 + phi.*cos(2.*pi.*X))).*Y.*cosh(2.*(1 + phi.*cos(2.*pi.*X)).*M) + Y.*(-1 + 2.*(1 + phi.*cos(2.*pi.*X)).^2.*M.^2- 2.*M.*Y.*cosh(M.*Y).*((1 + phi.*cos(2.*pi.*X)).*M.*cosh((1 + phi.*cos(2.*pi.*X)).*M)+ (-1 + (1 + phi.*cos(2.*pi.*X)).*M.^2.*beta).*sinh((1 + phi.*cos(2.*pi.*X)).*M))- 2.*(1 + phi.*cos(2.*pi.*X)).*M.*sinh(2.*(1 + phi.*cos(2.*pi.*X)).*M) + 2.*((1 + phi.*cos(2.*pi.*X)).*M.*(Y + (1 + phi.*cos(2.*pi.*X)).^2.*M.^2.*beta).*cosh((1 + phi.*cos(2.*pi.*X)).*M) + (-Y + (1 + phi.*cos(2.*pi.*X)).*M.*2.*((1 + phi.*cos(2.*pi.*X)).^2 - (1 + phi.*cos(2.*pi.*X)).*beta + Y.*beta)).*sinh((1 + phi.*cos(2.*pi.*X).*M)).*sinh(M.*Y)))./8.*((1 + phi.*cos(2.*pi.*X)).*M.*cosh((1 + phi.*cos(2.*pi.*X)).*M) + (-1 + (1 + phi.*cos(2.*pi.*X)).*M.^2.* beta).*sinh((1 + phi.*cos(2.*pi.*X)).*M)).^2));

我写了一个代码：

M = 0;
phi = 0.4;
beta = 0.0;
theta = 0.65;
F = theta - 1;
alpha = 0.0;

[X,Y]= meshgrid(linspace(-0.5,0.5),linspace(0,1.5));
Z = Y.*(F.*M.*cosh((1 + phi.*cos(2.*pi.*X)).*M)+ (1 + F.*M.^2*beta).*sinh((1 + phi.*cos(2.*pi.*X)).*M))-(F + (1 + phi.*cos(2.*pi.*X))).*sinh(M.*Y)/(1 + phi.*cos(2.*pi.*X)).*M.*cosh((1 + phi.*cos(2.*pi.*X)).*M)+ (-1 + (1 + phi.*cos(2.*pi.*X)).*M.^2.*beta).*sinh((1 + phi.*cos(2.*pi.*X)).*M )+ (alpha.*((F + (1 + phi.*cos(2.*pi.*X))).*Y.*cosh(2.*(1 + phi.*cos(2.*pi.*X)).*M) + Y.*(-1 + 2.*(1 + phi.*cos(2.*pi.*X)).^2.*M.^2- 2.*M.*Y.*cosh(M.*Y).*((1 + phi.*cos(2.*pi.*X)).*M.*cosh((1 + phi.*cos(2.*pi.*X)).*M)+ (-1 + (1 + phi.*cos(2.*pi.*X)).*M.^2.*beta).*sinh((1 + phi.*cos(2.*pi.*X)).*M))- 2.*(1 + phi.*cos(2.*pi.*X)).*M.*sinh(2.*(1 + phi.*cos(2.*pi.*X)).*M) + 2.*((1 + phi.*cos(2.*pi.*X)).*M.*(Y + (1 + phi.*cos(2.*pi.*X)).^2.*M.^2.*beta).*cosh((1 + phi.*cos(2.*pi.*X)).*M) + (-Y + (1 + phi.*cos(2.*pi.*X)).*M.*2.*((1 + phi.*cos(2.*pi.*X)).^2 - (1 + phi.*cos(2.*pi.*X)).*beta + Y.*beta)).*sinh((1 + phi.*cos(2.*pi.*X).*M)).*sinh(M.*Y)))./8.*((1 + phi.*cos(2.*pi.*X)).*M.*cosh((1 + phi.*cos(2.*pi.*X)).*M) + (-1 + (1 + phi.*cos(2.*pi.*X)).*M.^2.* beta).*sinh((1 + phi.*cos(2.*pi.*X)).*M)).^2));
contour(X,Y,Z,100)

[C,h] = contour(X,Y,Z);
set(h,'ShowText','on','TextStep',get(h,'LevelStep')*2)
colormap cool

我正进入（状态

Warning: Matrix is singular to working precision.
Warning: Contour not rendered for non-finite ZData
> In contour>parseargs at 192
  In contour at 69
Warning: Contour not rendered for non-finite ZData
> In contour>parseargs at 192
  In contour at 69

问题是什么？谢谢。

4

1 回答 1

1

你正在划分一个奇异矩阵：

(F + (1 + phi.*cos(2.*pi.*X))).*sinh(M.*Y)/ (1 + phi.*cos(2.*pi.*X))

我认为您在斜线之前缺少一个点：

M = 1;
phi = 0.4;
beta = 0.0;
theta = 0.65;
F = theta - 1;
alpha = 0.0;
[X,Y]= meshgrid(linspace(-0.5,0.5),linspace(0,1.5));
Z = Y.*(F.*M.*cosh((1 + phi.*cos(2.*pi.*X)).*M)+ (1 + F.*M.^2*beta).*sinh((1 + phi.*cos(2.*pi.*X)).*M))-(F + (1 + phi.*cos(2.*pi.*X))).*sinh(M.*Y)./(1 + phi.*cos(2.*pi.*X)).*M.*cosh((1 + phi.*cos(2.*pi.*X)).*M)+ (-1 + (1 + phi.*cos(2.*pi.*X)).*M.^2.*beta).*sinh((1 + phi.*cos(2.*pi.*X)).*M )+ (alpha.*((F + (1 + phi.*cos(2.*pi.*X))).*Y.*cosh(2.*(1 + phi.*cos(2.*pi.*X)).*M) + Y.*(-1 + 2.*(1 + phi.*cos(2.*pi.*X)).^2.*M.^2- 2.*M.*Y.*cosh(M.*Y).*((1 + phi.*cos(2.*pi.*X)).*M.*cosh((1 + phi.*cos(2.*pi.*X)).*M)+ (-1 + (1 + phi.*cos(2.*pi.*X)).*M.^2.*beta).*sinh((1 + phi.*cos(2.*pi.*X)).*M))- 2.*(1 + phi.*cos(2.*pi.*X)).*M.*sinh(2.*(1 + phi.*cos(2.*pi.*X)).*M) + 2.*((1 + phi.*cos(2.*pi.*X)).*M.*(Y + (1 + phi.*cos(2.*pi.*X)).^2.*M.^2.*beta).*cosh((1 + phi.*cos(2.*pi.*X)).*M) + (-Y + (1 + phi.*cos(2.*pi.*X)).*M.*2.*((1 + phi.*cos(2.*pi.*X)).^2 - (1 + phi.*cos(2.*pi.*X)).*beta + Y.*beta)).*sinh((1 + phi.*cos(2.*pi.*X).*M)).*sinh(M.*Y)))./8.*((1 + phi.*cos(2.*pi.*X)).*M.*cosh((1 + phi.*cos(2.*pi.*X)).*M) + (-1 + (1 + phi.*cos(2.*pi.*X)).*M.^2.* beta).*sinh((1 + phi.*cos(2.*pi.*X)).*M)).^2));
contour(X,Y,Z,100)
[C,h] = contour(X,Y,Z);
set(h,'ShowText','on','TextStep',get(h,'LevelStep')*2)
colormap cool

另外，请务必将 M 更改为不同于 0 的值

在此处输入图像描述

于 2012-04-07T01:00:04.537 回答