嵌入式数值计算

2024-06-11

嵌入式数值计算（精选四篇）

嵌入式数值计算篇1

随着矿山资源开发技术的日益发展,我国大部分的金属与非金属矿山开采深度不断增大,深部采矿已成为我国矿业发展不可避免的趋势。由于竖井与矿体的位置关系,开拓井巷的长度不断增加,工作量日益增大,因为掘进过程中凿岩和爆破产生的粉尘和炮烟等有毒、有害物质不能及时有效地排出而使掘进巷道的工作环境也越来越恶劣[1]。当前,掘进巷道通风工作主要依赖于局部通风,一个良好的通风环境在为矿山带来良好的经济效益的同时,也为职工的健康提供安全保证[2,3]。所以对掘进巷道通风除尘的研究对掘进工作有着重大的实际意义。

掘进巷道内的局部压入式通风,因为受巷道有限空间的限制,可以将其看作有限空间受限贴附射流[4,5,6]。而Fluent是一款能高效模拟复杂流动的流体模拟软件,它含有各种计算模型,并且有精度高和计算速度快等优点[7]。所以本文采用Fluent软件对掘进巷道压入式局部通风进行模拟,并对其风流流场分布规律进行分析,为矿山在巷道掘进通风除尘的实际工作提供数据参考。

2 数值模型

2.1 基本控制方程

在遵循质量守恒、动量守恒和能量守恒等基本物理守恒定律的前提下,假设掘进巷道内无热源,且不考虑掘进巷道周围围岩的热辐射,从而建立了该数值模型的控制方程组[8,9,10],主要方程式如下。

连续性方程:

式中:ρ——流体密度,kg/m3;

t—时间,s;

——速度矢量。

动量守恒方程:

式中:p——微元体上的压力,N;

τxx、τxy、τxz——粘性应力τ的分量;

Fx、Fy、Fz——微元体上的体积力,N。假设体积力只受重力作用,且z轴竖直向上,则Fx=0,Fy=0,Fz=ρg。

2.2 湍流模型

k-ε模型是当前模拟三维复杂湍流运用最多,也是最成功的模型之一,但经过许多专家学者进一步大量的模拟对比分析,发现对高应变率和流线弯曲程度较大的流动进行模拟时,RNG k-ε模型比标准的k-ε模型更合适,即能模拟射流撞击、旋流以及二次流等复杂的流动问题,并取得良好效果。掘进工作面风流流动时多数为湍流状态,且流线弯曲也较大,因此,采用RNG k-ε模型来模拟掘进工作面三维复杂湍流的运动。

RNG k-ε模型为:

2.3 计算区域及初始条件

以三星拱巷道为实例,巷道实际尺寸为:净断面面积S=4.21m2(宽为x=2.1m,墙高为h1=1.5m,拱高为h2=0.7m,其中风筒出口中心距巷道底面高为h3=1.1m)的计算区域(掘进巷道),风筒布置在巷道左侧壁面中部,风筒直径为d0=0.4m,设立风筒出口与掘进巷道工作面间距离为z=10m,风筒出口的平均风速为u0=6m/s。采用非结构化网格来离散,Pressure Based三维隐式求解,RNG k-ε模型,壁面无滑移,SIMPLEC压力速度耦合方式。

3 数值模拟结果与分析

3.1 掘进巷道压入式通风风流场分析

通过Fluent数值模拟计算,掘进巷道局部压入式通风风流场分布如图1和图2所示。

从图1中给出的不同X-Z水平截面速度矢量图得出:①图1a、图1b、图1c上的情形大致相同,风流从风筒中激射而出,向前运动过程中,由于风流与巷道壁面的粘附作用,在靠近壁面区域的风流形成了附壁射流区。②当射流到巷道尽头后,风流依然有较大速度,因前进方向受工作面限制,风流反向形成回流,并且回流区的最大速度值出现在与射流区相对的壁面。③另外从图中可明显看出在巷道中心区域有较明显的涡流存在,并且X-Z水平截面位于掘进巷道底板壁面时,涡流向射流区偏移,此时风流流动主要回流风速远大于射流的风速。

从图2中给出的不同截面速度矢量图可以看出:①与图1中射流区有所不同,在风筒出风口附近的Y-Z垂直截面,射流有自由射流的特性,因为风筒前方顶部与靠近巷道中心区域受巷道壁面影响较小,风流射流区自风筒出口射出之后吸卷周围气流开始逐渐扩大,如图2c和图2d中射流产生明显的扩展效应,与图1a的特征较相似。②另外,在巷道中心处,掘进巷道风流流场出现了明显的分界线,巷道靠风筒一侧主要以射流为主,相对壁面一侧以回流为主,如图2h所示该断面上几乎整个断面都为回流。③涡流主要存在于巷道的中下部区域,与射流区和回流区相比较不明显。④另外可以发现在射流区正下方区域风流流速非常低,尤其在风筒出口中心轴线下方。

3.2 掘进巷道内断面等速度云图分析

通过对建立的物理模型进行数值模拟计算,得出掘进工作面距射流出口不同距离的X-Y截面的等速度云图,如图3所示。

从图3中给出掘进工作面等速度云图得出:①从Z=0.0m和Z=1.0m断面云图中清晰显示出射流区对周围气流有明显的吸卷效应。②射流经过一段距离的发展,因为受巷道壁面的影响,射流流体的截面速度等值线会由圆形开始向半圆形转变,当距风筒出口距离约Z=6m时,射流截面几乎变成近似半圆形,同时射流开始析出空气,射流断面开始收缩,并且在这个过程中速度不断减小。③涡流主要集中在掘进巷道的中心部位,而回流主要集中在射流对面的壁面附近中下部区域内,而且距离射流区越远风速越大。④另外从图中也可以发现在风筒出口的射流区正下方有一小块区域和巷道顶部有部分区域风速很低。

4 结论

利用Fluent软件,并在适宜的简化基础上,对掘进工作面压入式通风风流流场进行数值模拟,定性地指出了掘进工作面压入式通风风流运动规律以及射流(风流)风筒出口中心轴线速度变化规律等流场分布,模拟结果得出以下结论。

(1)掘进巷道采用压入式通风时,风流流场存在4个明显的分区:附壁射流区、冲击射流区、回流区及涡流区。涡流主要分布在巷道中心部位,回流区与射流区形成对角分布,回流主要集中在风筒对面巷道壁面的中下部附近,而且距离射流区域越远风速越大。

(2)射流由风筒出口射出后,不是一直向外扩展,而是首先卷吸周围气流不断扩大,由于受到掘进巷道空间的限制及回流的影响,射流范围扩展到一定程度后,将向外析出空气,射流体断面开始不断收缩直至遇到掘进断面产生回流。

(3)巷道内风流分布很不均匀,在射流区正下方有一小块区域和巷道顶部部分区域内风速非常低,可能会出现通风盲区,不利于通风除尘,应加强通风管理。

参考文献

[1]张国枢.通风安全学[M].徐州:中国矿业大学出版社,2000.

[2]汪德淇.矿井通风防尘[M].北京:冶金工业出版社,1994.

[3]邬长福,汤民波,古鹏等.独头巷道局部通风数值模拟研究[J].有色金属科学与工程,2012,3(3):71-73.

[4]王海桥,施式亮,刘荣华等.压入式受限贴附射流流场特征及参数计算[J].黑龙江科技学院学报,2001,12(3):4-7.

[5]王淑芳,王剑波,张丽,王汝琳.局部通风机调速控制系统的研究[J].煤炭学报,2006,12(6):813-818.

[6]赵成龙,陈喜山.基于Fluent的独头巷道通风降温的数值模拟研究[J].现代矿业,2011,(2):33-36.

[7]温正,任毅如.FLUENT流体计算应用教程[M].北京:清华大学出版社,2009.

[8]周雪漪.计算水力学[M].北京:清华大学出版社,1995.

[9]陶文铨.数值传热学(第二版)[M].西安:西安交通大学出版社,2001.

[10]郭鸿志.传输过程数值模拟[M].北京:冶金工业出版社,1998.

[11]张兆顺,崔桂香,许春晓.湍流理论与模拟[M].北京:清华大学出版社,2005.163-165.

数值分析计算实习题篇2

姓名：

学号：

班级：

第二章

1、程序代码

Clear;clc;

x1=[0.2

0.4

0.6

0.8

1.0];

y1=[0.98

0.92

0.81

0.64

0.38];

n=length(y1);

c=y1(:);

for

j=2:n

%求差商

for

i=n:-1:j

c(i)=(c(i)-c(i-1))/(x1(i)-x1(i-j+1));

end

syms

df(1)=1;d(1)=y1(1);

for

i=2:n

%求牛顿差值多项式

df(i)=df(i-1)*(x-x1(i-1));

d(i)=c(i-1)*df(i);

end

P4=vpa(sum(d),5)

%P4即为4次牛顿插值多项式,并保留小数点后5位数

pp=csape(x1,y1,＇variational＇);%调用三次样条函数

q=pp.coefs;

q1=q(1,:)*[(x-.2)^3;(x-.2)^2;(x-.2);1];

q1=vpa(collect(q1),5)

q2=q(1,:)*[(x-.4)^3;(x-.4)^2;(x-.4);1];

q2=vpa(collect(q2),5)

q3=q(1,:)*[(x-.6)^3;(x-.6)^2;(x-.6);1];

q3=vpa(collect(q3),5)

q4=q(1,:)*[(x-.8)^3;(x-.8)^2;(x-.8);1];

q4=vpa(collect(q4),5)%求解并化简多项式

2、运行结果

0.98*x

0.3*(x

0.2)*(x

0.4)

0.625*(x

0.2)*(x

0.4)*(x

0.6)

0.20833*(x

0.2)*(x

0.4)*(x

0.8)*(x

0.6)

0.784

1.3393*x^3

0.80357*x^2

0.40714*x

1.04

1.3393*x^3

1.6071*x^2

0.88929*x

1.1643

1.3393*x^3

2.4107*x^2

1.6929*x

1.4171

1.3393*x^3

3.2143*x^2

2.8179*x

1.86293、问题结果

4次牛顿差值多项式=

0.98*x

0.3*(x

0.2)*(x

0.4)

0.625*(x

0.2)*(x

0.4)*(x

0.6)

0.20833*(x

0.2)*(x

0.4)*(x

0.8)*(x

0.6)

0.784

三次样条差值多项式

第三章

1、程序代码

Clear;clc;

x=[0

0.1

0.2

0.3

0.5

0.8

1];

y=[1

0.41

0.5

0.61

0.91

2.02

2.46];

p1=polyfit(x,y,3)%三次多项式拟合p2=polyfit(x,y,4)%四次多项式拟合y1=polyval(p1,x);

y2=polyval(p2,x);%多项式求值

plot(x,y,＇c--＇,x,y1,＇r:＇,x,y2,＇y-.＇)

p3=polyfit(x,y,2)%观察图像，类似抛物线，故用二次多项式拟合。

y3=polyval(p3,x);

plot(x,y,＇c--＇,x,y1,＇r:＇,x,y2,＇y-.＇,x,y3,＇k--＇)%画出四种拟合曲线

2、运行结果

-6.6221

12.8147

-4.6591

0.9266

2.8853

-12.3348

16.2747

-5.2987

0.9427

3.1316

-1.2400

0.73563、问题结果

三次多项式拟合P1=

四次多项式拟合P2=

二次多项式拟合P3=

第四章

1、程序代码

1)建立函数文件f.m:

function

y=fun(x);

y=sqrt(x)*log(x);

2)编写程序：

a.利用复化梯形公式及复化辛普森公式求解：

Clear;clc;

h=0.001;%h为步长,可分别令h=1,0.1,0.01,0.001

n=1/h;t=0;s1=0;s2=0;

for

i=1:n-1

t=t+f(i*h);

end

T=h/2*(0+2*t+f(1));T=vpa(T,7)

%梯形公式

for

i=0:n-1

s1=s1+f(h/2+i*h);

end

for

i=1:n-1

s2=s2+f(i*h);

end

S=h/6*(0+4*s1+2*s2+f(1));S=vpa(S,7)

%辛普森公式

a’复化梯形公式和复化辛普生公式程序代码也可以是：

Clear;clc;

x=0:0.001:1;

%h为步长,可分别令h=1,0.1,0.01,0.001

y=sqrt(x).*log(x+eps);

T=trapz(x,y);

T=vpa(T,7)

（只是h=1的运行结果不一样，T=1.110223*10^(-16)，而其余情况结果都相同）

Clear;clc;

f=inline(＇sqrt(x).*log(x)＇,x);

S=quadl(f,0,1);

S=vpa(S,7)

b.利用龙贝格公式求解：

Clear;clc;

m=14;%m+1即为二分次数

h=2;

for

i=1:m

h=h/2;n=1/h;t=0;

for

j=1:n-1

t=t+f(j*h);

end

T(i)=h/2*(0+2*t+f(1));%梯形公式

end

for

i=1:m-1

for

j=m:i+1

T(j)=4^i/(4^i-1)*T(j)-1/(4^i-1)*T(j-1);

%通过不断的迭代求得T（j)，即T表的对角线上的元素。

end

vpa(T(m),7)

2、运行结果

-0.4443875

-0.4444345

ans

-0.44444143、问题结果

a.利用复合梯形公式及复合辛普森公式求解：

步长h

0.1

0.01

0.001

梯形求积T=

[1.110223*10^(-16)]

-0.4170628

-0.4431179

-0.4443875

辛普森求积S=

-0.3267527

-0.4386308

-0.4441945

-0.4444345

b.利用龙贝格公式求解：

通过15次二分，得到结果：-0.4444414

第五章

1、程序代码

（1）LU分解解线性方程组：

Clear;clc;

A=[10

2.099999

2];

b=[8;5.900001;5;1];

[m,n]=size(A);

L=eye(n);

U=zeros(n);

flag=＇ok＇;

for

i=1:n

U(1,i)=A(1,i);

end

for

r=2:n

L(r,1)=A(r,1)/U(1,1);

end

for

i=2:n

for

j=i:n

z=0;

for

r=1:i-1

z=z+L(i,r)*U(r,j);

end

U(i,j)=A(i,j)-z;

end

abs(U(i,i))

flag=＇failure＇

return;

end

for

k=i+1:n

m=0;

for

q=1:i-1

m=m+L(k,q)*U(q,i);

end

L(k,i)=(A(k,i)-m)/U(i,i);

end

y=L＼b;x=U＼y

detA=det(L*U)

（2）列主元消去法：

function

gauss(A,b);

A=[10

1;-3

2.099999

2;5

-1;2

2];

b=[8;5.900001;5;1];

[n,n]

size(A);

zeros(n,1);

Aug

[A,b];

%增广矩阵

for

1:n-1

[piv,r]

max(abs(Aug(k:n,k)));

%找列主元所在子矩阵的行r

列主元所在大矩阵的行

r>k

temp=Aug(k,:);

Aug(k,:)=Aug(r,:);

Aug(r,:)=temp;

end

Aug(k,k)==0,error(‘对角元出现0’),end

把增广矩阵消元成为上三角

for

k+1:n

Aug(p,:)=Aug(p,:)-Aug(k,:)*Aug(p,k)/Aug(k,k);

end

解上三角方程组

Aug(:,1:n);

Aug(:,n+1);

x(n)

b(n)/A(n,n);

for

n-1:-1:1

x(k)=b(k);

for

p=n:-1:k+1

x(k)

x(k)-A(k,p)*x(p);

end

x(k)=x(k)/A(k,k);

end

detA=det(A)

2、运行结果

1）

LU分解解线性方程组

1.0e+006

0.0000

-0.0000

0.0000

-2.5000

0.0000

-2.4000

0.0000

1.0e+007

0.0000

-0.0000

0.0000

-0.0000

0.0000

1.5000

0.5750

0.0000

-0.0000

-1.0000

1.0000

detA

-762.0001

2）列主元消去法

detA

762.0001

ans

0.0000

-1.0000

1.0000

1.00003、问题结果

1）

LU分解解线性方程组

x=（0.0000，-1.0000，1.0000，1.0000）T

detA=-762.001

2)列主元消去法

x=（0.0000，-1.0000，1.0000，1.0000）T

detA=762.001

第六章

1、程序代码

（1）Jacobi迭代

Clear;clc;

%也可取n=8,10

hilb(n);

ones(n,1);

0.00001;

for

1:n

H(i,i)==0

＇对角元为零，不能求解＇

return

end

zeros(n,1);

kend

10000;

while

k<=kend

r>e

for

1:n

for

1:i

H(i,j)

x0(j);

end

for

1:n

H(i,j)

x0(j);

end

x(i)

b(i)

H(i,i)

H(i,i);

end

norm(x

x0,inf);

end

k>kend

＇迭代不收敛，失败＇

else

＇求解成功＇

end

（2）SOR迭代

1）程序代码

function

SOR(n,w);

hilb(n);

H*ones(n,1);

0.00001;

for

1:n

H(i,i)==0

‘对角线为零，不能求解’

return

end

zeros(n,1);

kend

10000;

while

k<=kend

r>e

for

1:n

for

1:i

H(i,j)

x(j);

end

for

1:n

H(i,j)

x0(j);

end

x(i)

x0(i)

H(i,i)

(b(i)

s);

end

norm(x

x0,inf);

end

k>kend

＇迭代不收敛，失败＇

else

＇求解成功＇

end

2）

从命令窗口中分别输入：

SOR（6,1）

SOR（8,1）

SOR（10,1）

SOR（6,1.5）

SOR（8,1.5）

SOR（10，1.5）

2、运行结果

Jacobi迭代：

ans

迭代不收敛，失败

SOR迭代：

第七章

1、程序代码

(1)不动点迭代法

1）建立函数文件：g.m

function

f=g(x)

f(1)=20/(x^2+2*x+10);

2)建立函数文件：bdd.m

function

[y,n]

bdd(x,eps)

nargin==1

eps=1.0e-8;

elseif

nargin<1

error

return

end

g(x);

while

(norm(x1-x)>=eps)&&(n<=10000)

x1;

g(x);

end

3)从命令窗口输入：bdd(0)

（2）牛顿迭代

clear;clc;

format

long;

m=8;

%m为迭代次数，可分别令m=2,4,6,8,10

x=sym(＇x＇);

f=sym(＇x^3+2*x^2+10*x-20＇);

df=diff(f,x);

FX=x-f/df;

%牛顿迭代公式

Fx=inline(FX);

disp(＇x=＇);

x1=0.5;

disp(x1);

Eps=1E-8;

k=0;

while

x0=x1;

k=k+1;

x1=feval(Fx,x1);

%将x1代入牛顿迭代公式替代x1

disp(x1);

%在屏幕上显示x1

k==m

break;

end

k,x12、运行结果

（1）不动点迭代法

bdd(0)

ans

1.3688

（2）

牛顿迭代

1.***

1.37***21

1.***

3、问题结果

(1)不动点迭代法

x=1.3688

n=25

收敛太慢

（2）牛顿迭代

初值取0

迭代次数k=8时，k

x(k)

1.4666666

1.3715120

1.3688102

1.3688081

嵌入式数值计算篇3

关键词：数值计算方法；创新意识；计算平台

作者简介：张俊丽（1980-），女，山东菏泽人，内蒙古民族大学数学学院，讲师。（内蒙古?通辽?028000）

中图分类号：G642.0?????文献标识码：A?????文章编号：1007-0079（2012）28-0087-01

随着科技的飞速发展和计算机技术的广泛应用，数值计算方法已成为重要的桥梁和工具深入到航天航空、地质勘探、汽车制造、桥梁设计、天气预报等各个领域，成为每一位科研人员和工程技术人员所必备的知识。为了满足社会需求，数值计算方法现已成为高等院校理工类学生的一门专业必修课程，其目的是让学生掌握设计数值算法的基本方法，培养学生分析问题与解决问题的能力，为以后用计算机解决科学计算问题打下坚实的基础。

一、“数值计算方法”课程的特点与教学现状

数值计算方法，简称计算方法，又叫数值分析，是一门研究数学问题的近似解并利用计算机进行数值实现的学科，是数学分析、高等数学、高等代数、概率统计等数学基础课的后续课程，它既有数学理论上的抽象性与严谨性，又有实验性与应用性的数值特征。计算方法课程的内容包括插值和拟合、数值微分和数值积分、求解线性方程组的数值方法（直接法和迭代法）、非线性方程数值解、矩阵特征值计算及常微分方程初值问题数值解法等；[2]它的计算对象是数学中的微积分、线性代数、常微分方程，只是它不像别的数学课程那样只是研究纯粹的数学理论，而是把数学理论与计算相结合，重点探讨数学问题的数值解法及应用；它的课程要求是在掌握算法原理的前提下设计算法编程实现。

“数值计算方法”是一门介绍科学计算的核心理论和基本方法的数学课程，它对培养学生的科学计算能力和解决实际问题的能力具有不可替代的作用。从20世纪80年代起，“数值计算方法”相继成为各高等院校数学及其他理工科（如物理、计算机等）专业本科生的一门专业基础课。但内蒙古民族大学（以下简称“我校”）的数值计算方法课程只在应用数学、信息与计算科学两个专业开设必修课，一般开设在第三或第四学期，理论课48学时，上机实验16学时，在别的学院（如物理、计算机等）没有开设该课程。该课程普遍存在的教学现状是：理论课内容多，学时少，各部分内容不连贯，公式繁多，枯燥乏味，使得学生产生厌学情绪；上机课时间紧，且一般集中上机，与理论课内容脱节，失去了上机实验操作的意义；很多时候这门课程的学习都结束了，学生还不清楚这门课程与原来的课程有什么联系，学习这门课有什么用，更无从谈起培养学生的创新能力；而且“数值计算方法”课程教学过程中还存在着教学内容陈旧、教学方式落后及考试形式单一等问题。针对该课程目前的教学现状，如何对该课程教学进行教学改革，是值得深入思考的问题。

二、关于“数值计算方法”课程改革的若干建议

根据前文分析可知，目前“数值计算方法”课程教学中存在着一些不容忽视的问题。那么如何进行教学改革，培养学生的实际应用能力，体现该课程在工程科学中的价值和意义，是值得数学界思考的问题。根据近年来我校师生在该课程教学中出现的问题，本文对“数值计算方法”课程教学改革提出以下几点建议：

1.优化教学内容，选择合适教材

“数值计算方法”课程讲授时既要强调它的理论结构与使用价值，又要注重提升它与计算机使用密切结合的实用性特点，所以该门课程对教材的要求很高。然而现行教材有的理论偏深，不适合普通本科生使用；有的内容陈旧，与实际联系缺乏；有的实用性强，但与实践结合的算例较少；[3]再加上该课程内容抽象，知识连贯性不强，定理和公式较多，推导过程烦琐，从而导致学生对该课程的学习没有兴趣，只是为了应付考试机械性地记忆公式。按照教育部关于“数值计算方法”课程在教学过程中应把握“重概念、重方法、重应用、重能力”的培养要求，对该课程的教学内容应灵活把握，知识点讲解应详略得当，不同专业的学生对该课程的要求不同，讲解的侧重点也应有所不同，最好选用的教材也不同。对数学类的学生来说，理论与实践应并重，而对于别的理工科的学生来说，不在于理论的论证与推导，而应侧重算法原理与实际应用。当选定教材后，在实际教学过程中还需要对教学内容灵活整合，对于一些复杂且后继课程将会深入学习的内容（例如微分方程的数值解法等），[4]可以略讲甚至不讲。不同地区的高校对该课程的教学要求也略有不同，例如我校处少数民族地区，学生的基础知识相对较差，在该课程授课时更应减少烦琐公式的推导，重在加强学生对知识点的掌握与实际应用能力的培养。鉴于该课程对以后学习和工作的重要性，我校建议除了数学与信息类的学生以外，别的理工科（如物理，计算机、信息工程等）的学生也应开设数值计算方法课程的选修课。我院本专业教师在包玉兰教授的带领下，根据我校学生的状况及多年积累的教学经验，编写了比较适合少数民族地区学生特点的数值计算方法教材，现已经出版在我校试用。该教材内容较浅，并配备一定量的习题和上机实验题，要求理论学时50～60学时（包含习题课），上机实验16～20学时，并且标注了一些选讲的内容，不同专业的学生可以针对性地学习，[5]基本上满足了我校学生对该课程教材的要求。

2.转变教学模式，活跃课堂气氛

“数值计算方法”是一门对数学问题进行数值求解的课程，主要培养学生算法思想与科学计算能力。因此，在教学过程中，面对烦琐的公式推导、累积的误差、多次迭代与数据处理等问题，教师必须改变传统的黑板加粉笔的教学方式。如果将计算机多媒体教学恰当地引入数值计算方法课堂，利用多媒体技术生动、形象、鲜明的特点，既可以保留传统教学中教师与学生面对面交流的优势，又会使某些抽象、枯燥、难以理解的概念、理论及冗长公式推导变得直观、形象。例如，对于Runge振荡现象，传统教学很难说清楚怎么回事，如果借用多媒体，通过选取不同的等距插值节点，将相应的插值图形动态地描述出来，学生马上就能理解振荡的原因。[4]另外，还可以介绍一些与“计算方法”课程相关的学术热点问题以及数学竞赛、数学模型中的典型算例，例如对于不同的数学模型，有不同类型的数值算法，即使同一模型也有多种数值解法，这些算法都有各自的应用背景，在教学中应根据背景、目的和设计出发点的不同引导学生积极思考各种算法。[6]

3.提升计算平台，改变考试方式

“数值计算方法”是一门连接传统数学理论和实际应用的课程，该课程教学关键的一步就是对算法进行上机实现，培养学生实际应用能力。因此，该课程的算法思想应该随着计算机的发展不断更新。但目前该课程很多版本的教材偏重理论，实际算例少，用传统的计算语言编写算法思想及算法图解，不能及时建立与先进算法衔接的平台。现阶段计算方法教材里的程序多是采用C语言编写的，学生上机时难度大，积极性不高，再加上我校学生对计算机语言掌握得不够，上机实验课的时间不能得到有效利用。随着MATLAB等计算平台的出现，传统C语言编写的程序渐已退出，因此不必对数值计算方法教材里的算法细节过多考虑，更多考虑对计算结果的影响因素（如初值的选取、参数的选择、迭代次数等）。[7]实验课程的上机时间不易集中安排在理论课结束以后，应该在每一章节理论课结束后结合具体算例带领学生上机操作，这样既可以巩固理论知识，又可培养学生实践能力，激发学习兴趣。

目前“数值计算方法”课程仍然采用传统的一卷定音的考试形式，这样不能调动学生上机实验的积极性，更不能全面考查学生对该门课程内容掌握的真实情况。根据近年来的教学体验，笔者认为应选用适合本课程特点的科学的考试方式，具体实施如下：若作为专业必修课，考核内容包括平时考核、上机实验与期末考试，比例分别为10%、30%和60%；若为选修课，考试方式可以灵活掌握，也许一个上机实例便会成为一份很好的考核学生能力的试卷，不一定拘泥于闭卷考试的形式，这样更有利于学生的学习积极性和综合能力的提高。

三、结束语

随着素质教育呼声的不断高涨，在高等学校进行教学改革已是时代发展的必然，在“数值计算方法”课程教学改革过程中，应该以培养适应时代发展的新型人才为目标，鼓励教师不断完善知识体系、创新教学方式，充分调动学生的主观能动性，有计划有步骤地提高该课程的整体教学水平。

参考文献：

[1]任铭,李振平,余亚辉.关于数值计算方法课程教学改革的思考与探讨[J].中国校外教育,2009,(S4).

[2]李庆扬,王能超,易大义.数值分析[M].北京:清华大学出版社,2003.

[3]栾孟杰.数值计算方法课程的教学改革研究[J].林区教学,2011,(4).

[4]李小林.关于数值计算方法课程教学改革的探讨[J].重庆文理学院学报(自然科学版),2010,(2).

[5]包玉兰,吉日木图,肖丽霞.关于计算方法课程教学改革的探讨[J].内蒙古民族大学学报(自然科学版),2007,(6).

[6]杜廷松.关于《数值分析》课程教学改革研究的综述和思考[J].大学数学,2007,(2).

[7]李秀梅.计算方法课程教学改革初探[J].中国电力教育,2009,(6).

嵌入式数值计算篇4

关键词：压入式通风,掘进工作面,数值模拟

引言

矿井掘进工作面处于地下暴露的岩体中,工作面的粉尘浓度大和有害气体浓度高。这些灾害严重影响采掘工人的作业,若长期在比较差的环境下工作,必然会影响到人体健康,并容易引发事故和灾害[1、2]。因此,对掘进工作面进行通风数值模拟对粉尘和有害气体扩散和降低工作面温度有重要意义,为改善矿井掘进工作面通风提供一些指导。

目前,局扇通风是独头巷道通风最常用的一种通风方式。局扇通风又分为压入式通风、抽出式通风和混合式通风。其中压入式通风是局部通风机及其附属装置安装在离掘进巷道口10m以为的进风侧,将新鲜风流经风筒输送到掘进工作面,污风沿掘进巷道排出。

1 数值模型的建立

1.1 数学模型

1.1.1 基本假设

在确定数学模型之前,作出如下假设:

(1)通风气流可视为不可压缩流体,并忽略由流体粘性力做功所引起的耗散热。同时壁面绝热,等温通风;

(2)流动的紊流粘性具有各向同性,紊流粘性系数μt作为标量处理;

(3)流动为稳态紊流,满足Boussinesq假设。紊流模型采用k-ε双方程模型。

1.1.2 数学模型的建立

基于以上假设条件,可以得到模型的控制方程:

(1)质量守恒方程(Mass Conservation Equation)为

(2)动量守恒方程(Monentum Conservation Equation)为

式中:xj为坐标系中j方向上的空间位置;uj为坐标系中j方向上的速度;p为静压;Tij为应力张量;gi和Fi分别为i方向上的重力和外部体积力。

(3)RNGk-ε方程为

式中,k、ε、xi、ρ、μ和μt分别表示紊流动能(m2/s2)、紊流动能耗损率(m2/s3)、坐标、空气密度、层流动力粘性系数(Pa·s)和紊流粘性系数(Pa·s)。Gk为平均速度梯度引起的紊动能产生项。方程中各经验常数的取值为:C1ε=1.44,C2ε=1.92,σk=1.00,σε=1.30。

1.2 物理模型

本次数值模拟选择局扇压入式通风,为便于分析,选一截面下部为宽为2.2m,高为1.7m的矩形上部为直径2.2m的半圆拱型的水平独头巷道作为数值模拟对象。风筒出口分别位于巷道右侧,风筒中心距离底板1.5m,巷道底部和巷道正上方距顶板。风筒直径为0.5 m。模型计算区域为风筒出口到独头巷道迎头(掘进工作面)。根据经验公式在此模型中利用压入式通风时,风筒出口距工作面的距离为10 m。如图1所示。采用压入式通风时,风流的有效射程满足式:

式中S—巷道断面,m2

1.3 模型网格划分和边界条件设置

Fluent是目前国际上比较流行的商业CFD(Computational Fluid Dynamics)计算流体动力学软件包之一[3]。它是用于模拟具有复杂外形流体流动和热传导的计算机程序。CFD软件包含了3个主要的功能部分:前处理器、求解器和后处理器。用GAMBIT软件建立三维掘进工作面的模型并对模型进行网格划分计算。三个模型分别划分成88218、89188和99895个网格单元。然后使用FLUENT软件进行方程组的计算。

结合数学模型和FLUENT的模拟方法[4,5,6],确定数值模拟的各参数及边界条件。

FLUENT软件选用隐式分离三维稳定流求解器,湍流模型(Viscous Model)选用k-epsilon模型(k-ε两方程模型);求解器参数中压力速度耦合方式设置为SIMPLEC,压力离散方式设置为标准格式,离散格式设置为一阶迎风格式,收敛标准为0.001。

设定风筒出口处为模型的入口边界,入口速度为17.8 m/s,类型为Velocity-inlet。

设定出口边界为巷道自由断面处,没有相对压力,类型为Outflow。

设定壁面边界条件:所有壁面采用无滑动边界条件,类型为No Slip。

设定水力直径和湍流强度:水力直径为2.49m,湍流强度为3.72%。其中水力直径按计算公式:

式中,A为风流断面面积,/m2;

S为流体与固体接触周长,/m。

湍流强度计算:

ReH按水力直径计算的Reynolds数。

设定材料属性为Air(空气),其中密度为1.225 kg/m3、黏性系数1.81×10-5Pa·s。

2 数值模拟结果及分析

根据模型边界条件,用FLUENT软件进行模拟,然后用Tecplot软件进行后处理。得到不同情况下独头巷道工作面风速稳定后的风速分布。如图,图分别为抽出式风筒不同位置Y=1.5m的高度的巷道内速度矢量图。

由图中可以看出,风流从风筒出口流出时是以逐渐扩大的自由风流状态射向工作面沿自由风流的轴线方向风速逐渐降低到一定距离之后再反向流至巷道口。当风筒处于顶板和巷道底部时,距离地面1.5m高度的巷道工作面附近空气形成紊流和循环涡流区,工作面的风速不利于瓦斯和粉尘扩散;当风筒中心处于Y=1.5m的时候,距离地面1.5m高度的巷道工作面附近新鲜风清洗工作面后,稳定流出巷道出口。

3 结语

(1)基于Fluent软件的掘进巷道风流状态模拟,能较准确的模拟风流的状态,结果可信、直观、可视性强。

(2)模拟结果显示,由于受独头巷道局限空间的限制和风流的连续性,射流体到达巷道掘进工作面,由于受壁面的影响,射流开始冲击并附壁回转。

(3)在风筒中心处于离底板1.5m的位置,距离地面1.5m高度的巷道工作面附近新鲜风清洗工作面后,稳定流出巷道出口。给工作面工人提供新鲜风流。

(4)影响独头巷道通风效果的因素有很多,对于其他影响因素可以结合有害气体扩散和掘进工作面温度场做进一步研究。

参考文献

[1]刘钊春,柴军瑞,贾晓梅等.压入式通风掘进面有害气体浓度扩散数值模拟[J].岩土力学,2009,30(2).

[2]刘荣华,王海桥,施式亮等.压入式通风掘进工作面粉尘分布规律研究[J].煤炭学报,2002(3)

[3]王瑞金,张凯,王刚.Fluent技术基础与应用实例[M].北京:清华大学出版社,2007.

[4]韩占忠,王敬,兰小平.FLUENT-流体工程仿真计算实例与应用[M].北京:北京理工大学出版社,2004

[5]王晓珍,蒋仲安,王善文等.煤巷掘进过程中粉尘浓度分布规律的数值模拟[J].煤炭学报,2007,(4)

本文来自 360文秘网(www.360wenmi.com)，转载请保留网址和出处