前言
在天津大学流体力学实验室中,研究人员常采用粒子图像测速技术(PIV,Particle Image Velocimetry)对湍流边界层进行细致的流动测量。当对精度要求不断提高时,还需要确定流体的粘度。
粘度分为两种:动力粘度μ和运动粘度ν。其中动力粘度定义为流体切应力与切应变之间的比例系数,量纲为
,单位为,而 运动粘度
,量纲为
,单位为
。
本文对动力粘度 μ 进行窄温度范围内的二次函数拟合。
拟合方法
在网上能够查询到水动力粘度与温度之间的关系表格,如图。但是不能够得到连续温度的粘度,常用的方法有插值和最小二乘曲线拟合,本文介绍最小二乘拟合二次函数的方法。
图:水动力粘度温度对照表,来源:https://wenwen.sogou.com/question/q827107145.htm
实验室内温度一般不低于14℃,不高于27℃,因此只需要这一区间内温度的平滑曲线即可。取14~27℃水的动力粘度,对二次函数进行最小二乘拟合(matlab提供最小二乘拟合函数lsqcurvefit)。拟合结果为:
拟合效果与误差如图。误差均在1‰以下。注意:此公式只适用于14~27℃的粘度计算,不同原始数据拟合结果不同。
图:二次函数拟合效果
图:每个温度的误差
拟合代码
%% 14~27℃水动力粘度的二次函数最小二乘拟合
%% 数据给定
T=[14,15,16,17,18,19,20,...
21,22,23,24,25,26,27]; % 温度列表
D_v=[1.1709e-3,1.1404e-3,1.111e-3,1.0828e-3,1.0559e-3,1.0299e-3,1.005e-3,...
0.981e-3,0.9579e-3,0.9358e-3,0.9142e-3,0.8937e-3,0.8737e-3,0.8545e-3]; % 动力粘度列表
%% 拟合
fun=@(W,T)W(1).*T.^2+W(2).*T+W(3); % 待拟合函数
W=[4.7e-07,-4.3e-05,1.7e-03]; % 给定初值
W=lsqcurvefit(fun,W,T,D_v); % 拟合结果为 W=[4.692308e-07,-4.346549e-05,1.686702e-03];
T0=min(T-0.5):0.1:max(T+0.5);
D_v0=fun(W,T0);
%% 拟合效果与误差
figure
plot(T0,D_v0)
hold on
plot(T,D_v,'o')
xlabel('T (℃)')
ylabel('μ (Pa·s)')
legend('data','parabolic fit')
figure
bar(T,abs((D_v-fun(W,T))./D_v*10000))
xlabel('T (℃)')
ylabel('err_μ (‱)')
%% 显示公式
fprintf('Dynamic Viscosity =\n %e * Temperature^2\n%e * Temperature\n+%e',W(1),W(2),W(3))
syms T
mu=vpa(W(1)*T^2+W(2)*T+W(3),7)