CFD作业：采用代数法生成结构网格

CFD本科作业，采用代数法生成结构网格

骆苏

1565人浏览 · 2022-04-20 20:15:28

骆苏 · 2022-04-20 20:15:28 发布

将变换关系定义为

${emξ−1em−1=x−xminxmax−xminenη−1en−1=y−y1(x)y2(x)−y1(x)\left\lbrace \begin{array}{l} \frac{e^{m\xi } -1}{e^m -1}=\frac{x-x_{min} }{x_{max} -x_{min} }\\ \frac{e^{n\eta } -1}{e^n -1}=\frac{y-y_1 \left(x\right)}{y_2 \left(x\right)-y_1 \left(x\right)} \end{array}\right.$

化简有
${x=xmin+(xmax−xmin)emξ−1em−1y=y1(x)+(y2(x)−y1(x))enη−1en−1\left\lbrace \begin{array}{l} x=x_{min} +\left(x_{max} -x_{min} \right)\frac{e^{m\xi } -1}{e^m -1}\\ y=y_1 \left(x\right)+\left(y_2 \left(x\right)-y_1 \left(x\right)\right)\frac{e^{n\eta } -1}{e^n -1} \end{array}\right.$

对于本题 $y1(x)y_1 \left(x\right)$ 、 $y2(x)y_2 \left(x\right)$ 均为常数， $y1(x)=yminy_1 \left(x\right)=y_{min}$ ， $y2(x)=ymaxy_2 \left(x\right)=y_{max}$ ，本题中，取为

${xmin=0xmax=10ymin=0ymax=10\left\lbrace \begin{array}{l} x_{min} =0\\ x_{max} =10\\ y_{min} =0\\ y_{max} =10 \end{array}\right.$

因此只需取值 $m$ 与 $n$ 来保证 $(x1,j−x0,j)∈(0.05,0.1)\left(x_{1,j} -x_{0,j} \right)\in \left(0.05,0.1\right)$ ， $(yi,1−yi,0)∈(0.05,0.1)\left(y_{i,1} -y_{i,0} \right)\in \left(0.05,0.1\right)$ 即可。

已经满足题目要求，即给出答案为

${a=10b=10m=1n=1\left\lbrace \begin{array}{l} a=10\\ b=10\\ m=1\\ n=1 \end{array}\right.$

但我们进一步分析一下这些参数的影响

为了有比对效果，我们先将x、y方向的分区数都提高到20，可以看到，间隔相应减小。
稠密化

再只调整y方向的参数而保持x方向的参数不变。

调整n的取值，发现随着n的增大，y方向上的不均匀性更强了。

参数影响

我们再随便画个拉瓦尔喷管试试
${y1(x)=x(10−x)10y2(x)=10−x(10−x)10\left\lbrace \begin{array}{l} y_1 \left(x\right)=\frac{x(10-x)}{10}\\ y_2 \left(x\right)=10-\frac{x(10-x)}{10} \end{array}\right.$

拉瓦尔喷管
程序还算有效👌。

%% 图1
%% 程序初始化
close;clear,clc

%% 调整参数
m = 1;
n = 1;

%% 选取分区数
count_x = 10;
count_y = 10;

%% 描述计算区域
x_max = 10;
x_min = 0;
y_1 = @(x) 0;
y_2 = @(x) 10;

%% 计算划分网格
xi = (0:1:count_x)/count_x;
eta = (0:1:count_y)/count_y;

for i = 1:(count_x+1)
    for j = 1:(count_y+1)
        x(i,j) = x_min + (x_max-x_min)*(exp(m*xi(i))-1)/(exp(m)-1);
        y(i,j) = y_1(x(i,j)) + (y_2(x(i,j))-y_1(x(i,j)))*(exp(n*eta(j))-1)/(exp(n)-1);
    end
end

%% 绘制划分结果
figure(1);
for i = 1:(count_x)
    for j = 1:(count_y)
        line([x(i,j) x(i,j+1)],[y(i,j) y(i,j+1)],'color','k','linewidth',2);
        line([x(i,j) x(i+1,j)],[y(i,j) y(i+1,j)],'color','k','linewidth',2);
    end
end
for i = 1:(count_x)
    line([x(i,count_y+1) x(i+1,count_y+1)],[y(i,count_y+1) y(i+1,count_y+1)],'color','k','linewidth',2);
end
for j = 1:(count_y)
    line([x(count_x+1,j) x(count_x+1,j+1)],[y(count_x+1,j) y(count_x+1,j+1)],'color','k','linewidth',2);
end
axis equal;
title({['m = ',num2str(m,'%d'),'，n = ',num2str(n,'%d'),'，x、y方向各分',num2str(count_x,'%d'),'、',num2str(count_y,'%d'),'个区域'],...
    ['x_{1,0} - x_{0,0} = ', num2str((x(2,1)-x(1,1)),'%.3f'),'，y_{0,1} - y_{0,0} = ', num2str((y(1,2)-y(1,1)),'%.3f')]});
xlabel('x');
ylabel('y');
xlim([-2.00 12.00])
ylim([-2.00 12.00])

%% 后面的几个图
%% 程序初始化
clear;

%% 调整参数
m = 1;
n = 1;

%% 选取分区数
count_x = 20;
count_y = 20;

%% 描述计算区域
x_max = 10;
x_min = 0;
y_1 = @(x) 0;
y_2 = @(x) 10;

%% 计算划分网格
xi = (0:1:count_x)/count_x;
eta = (0:1:count_y)/count_y;

for i = 1:(count_x+1)
    for j = 1:(count_y+1)
        x(i,j) = x_min + (x_max-x_min)*(exp(m*xi(i))-1)/(exp(m)-1);
        y(i,j) = y_1(x(i,j)) + (y_2(x(i,j))-y_1(x(i,j)))*(exp(n*eta(j))-1)/(exp(n)-1);
    end
end

%% 绘制划分结果
figure(2);
for i = 1:(count_x)
    for j = 1:(count_y)
        line([x(i,j) x(i,j+1)],[y(i,j) y(i,j+1)],'color','k','linewidth',2);
        line([x(i,j) x(i+1,j)],[y(i,j) y(i+1,j)],'color','k','linewidth',2);
    end
end
for i = 1:(count_x)
    line([x(i,count_y+1) x(i+1,count_y+1)],[y(i,count_y+1) y(i+1,count_y+1)],'color','k','linewidth',2);
end
for j = 1:(count_y)
    line([x(count_x+1,j) x(count_x+1,j+1)],[y(count_x+1,j) y(count_x+1,j+1)],'color','k','linewidth',2);
end
axis equal;
title({['m = ',num2str(m,'%d'),'，n = ',num2str(n,'%d'),'，x、y方向各分',num2str(count_x,'%d'),'、',num2str(count_y,'%d'),'个区域'],...
    ['x_{1,0} - x_{0,0} = ', num2str((x(2,1)-x(1,1)),'%.3f'),'，y_{0,1} - y_{0,0} = ', num2str((y(1,2)-y(1,1)),'%.3f')]});
xlabel('x');
ylabel('y');
xlim([-2.00 12.00])
ylim([-2.00 12.00])

%% 程序初始化
clear;

%% 调整参数
m = 1;
n = 10;

%% 选取分区数
count_x = 20;
count_y = 20;

%% 描述计算区域
x_max = 10;
x_min = 0;
y_1 = @(x) 0;
y_2 = @(x) 10;

%% 计算划分网格
xi = (0:1:count_x)/count_x;
eta = (0:1:count_y)/count_y;

for i = 1:(count_x+1)
    for j = 1:(count_y+1)
        x(i,j) = x_min + (x_max-x_min)*(exp(m*xi(i))-1)/(exp(m)-1);
        y(i,j) = y_1(x(i,j)) + (y_2(x(i,j))-y_1(x(i,j)))*(exp(n*eta(j))-1)/(exp(n)-1);
    end
end

%% 绘制划分结果
figure(3);
for i = 1:(count_x)
    for j = 1:(count_y)
        line([x(i,j) x(i,j+1)],[y(i,j) y(i,j+1)],'color','k','linewidth',2);
        line([x(i,j) x(i+1,j)],[y(i,j) y(i+1,j)],'color','k','linewidth',2);
    end
end
for i = 1:(count_x)
    line([x(i,count_y+1) x(i+1,count_y+1)],[y(i,count_y+1) y(i+1,count_y+1)],'color','k','linewidth',2);
end
for j = 1:(count_y)
    line([x(count_x+1,j) x(count_x+1,j+1)],[y(count_x+1,j) y(count_x+1,j+1)],'color','k','linewidth',2);
end
axis equal;
title({['m = ',num2str(m,'%d'),'，n = ',num2str(n,'%d'),'，x、y方向各分',num2str(count_x,'%d'),'、',num2str(count_y,'%d'),'个区域'],...
    ['x_{1,0} - x_{0,0} = ', num2str((x(2,1)-x(1,1)),'%.3f'),'，y_{0,1} - y_{0,0} = ', num2str((y(1,2)-y(1,1)),'%.3f')]});
xlabel('x');
ylabel('y');
xlim([-2.00 12.00])
ylim([-2.00 12.00])

%% 程序初始化
clear;

%% 调整参数
m = 1;
n = 1;

%% 选取分区数
count_x = 20;
count_y = 20;

%% 描述计算区域
x_max = 10;
x_min = 0;
y_1 = @(x) x*(10-x)/10;
y_2 = @(x) 10 - x*(10-x)/10;

%% 计算划分网格
xi = (0:1:count_x)/count_x;
eta = (0:1:count_y)/count_y;

for i = 1:(count_x+1)
    for j = 1:(count_y+1)
        x(i,j) = x_min + (x_max-x_min)*(exp(m*xi(i))-1)/(exp(m)-1);
        y(i,j) = y_1(x(i,j)) + (y_2(x(i,j))-y_1(x(i,j)))*(exp(n*eta(j))-1)/(exp(n)-1);
    end
end

%% 绘制划分结果
figure(4);
for i = 1:(count_x)
    for j = 1:(count_y)
        line([x(i,j) x(i,j+1)],[y(i,j) y(i,j+1)],'color','k','linewidth',2);
        line([x(i,j) x(i+1,j)],[y(i,j) y(i+1,j)],'color','k','linewidth',2);
    end
end
for i = 1:(count_x)
    line([x(i,count_y+1) x(i+1,count_y+1)],[y(i,count_y+1) y(i+1,count_y+1)],'color','k','linewidth',2);
end
for j = 1:(count_y)
    line([x(count_x+1,j) x(count_x+1,j+1)],[y(count_x+1,j) y(count_x+1,j+1)],'color','k','linewidth',2);
end
axis equal;
title({['m = ',num2str(m,'%d'),'，n = ',num2str(n,'%d'),'，x、y方向各分',num2str(count_x,'%d'),'、',num2str(count_y,'%d'),'个区域'],...
    ['x_{1,0} - x_{0,0} = ', num2str((x(2,1)-x(1,1)),'%.3f'),'，y_{0,1} - y_{0,0} = ', num2str((y(1,2)-y(1,1)),'%.3f')]});
xlabel('x');
ylabel('y');
xlim([-2.00 12.00])
ylim([-2.00 12.00])