實驗資料的數值積分問題
用不同溫度下固體鉛的定壓熱容資料,求298k時的絕對熵。資料見下表:
由熱力學方法可知,絕對熵與定壓熱容有如下關係:
採用matlab來計算,這裡將最小二乘b樣條擬合法和三次樣條插值函式計算作比較,命令列程式如下:
function entropy_int
clear all
clc% 讀入資料
t = [5. 10. 15.
20. 25. 30.
50. 70. 100.
150. 200. 250.
298.];
cp = [0.305 2.8 7.
0 10.8 14.1 16.
5 21.4 23.3 24.
5 25.4 25.8 26.
2 26.5];
cpt = cp./t;
t = [0 t]; cpt = [0 cpt];
% 用最小二乘b樣條擬合法和三次樣條插值函式計算
knots = 3; k = 3三次b樣條
sp = spapi(knots,k,t,cpt);
cs = csapi(t,cpt); % 生成三次樣條插值函式cs
% 繪製濃度擬合曲線
ti = linspace(t(1),t(end),299);
cpt_b = fnval(sp,ti);
cpt_c = fnval(cs,ti);
plot(t,cpt,'ro',ti,cpt_b,'b-',ti,cpt_c,'bl-.')
xlabel('t')
ylabel('cp/t')
legend('實驗值','b樣條擬合','立方樣條擬合')
% 進行數值積分
pp = fnint(cs);
s = fnval(pp,t)
s = 0 0.0433 0.8733 2.7820 5.3352 8.1104 10.9071
20.7236 28.2599 36.8085 46.9207 54.2919 60.0880 64.7233
在解決本問題時,需首先要考慮積分下限,即時的情形,因此時的為0/0。為此,可考慮先求出接近0k時的數值。可有多種方法此數值,如通過圖的樣條插值曲線讀出,或採用回歸擬合以及插值計算的方法。
這裡採用上述方法得到的數值是,。由於本題的溫度資料間隔不等,直接採用樣條插值積分是一較好的方法。運用imsl數學庫中的樣條插值csint和依據樣條插值係數的積分csitg,可以求得絕對熵,下面給出呼叫的主程式和計算結果。
c to obtain the absolute entropy of solid lead at 298 k from
c the experimental data of the heat capacity at constant pressure cp of
c solid lead at various temperatures up to and including 298 k
c by the cubic spline interpolant and integration from imsl library.
c the absolute entropy in crc handbook of chemistry and physics is 64.8 j k-1 mol-1.
integer ndata, i, nintv, nout
parameter (ndata=14)
real a, b, break(ndata), cscoef(4,ndata), csitg, error,
& exact, f, fdata(ndata), fi, float, value, x, xdata(ndata)
intrinsic float
external csint, csitg, umach
data (xdata(i),i=1,ndata)/0.05,5.,10.,15.,20.,25.,30.,50.,
70.,100.,150.,200.,250.,298./
data (fdata(i),i=1,ndata)/0.003,0.305,2.8,7.0,10.8,14.1,16.5,
21.4,23.3,24.5,25.4,25.8,26.2,26.5/
c compute cubic spline interpolant
do 10 i=1,ndata
fdata(i)=fdata(i)/xdata(i)
10 continue
c cubic spline interpolant
call csint (ndata, xdata, fdata, break, cscoef)
c cubic spline integration
a = 0.0
b = 298.
nintv = ndata - 1
value = csitg(a,b,nintv,break,cscoef)
exact = 64.8
error = exact - value
c get output unit number
call umach (2, nout)
c print the result
write (nout,99999) a, b, value, exact, error
99999 format (' on the closed interval (', f4.1, ',', f5.1,
we h**e :', /, 1x, 'computed integral = ', f10.5, /,
& 1x, 'exact integral = ', f7.2, /, 1x, 'error '
f10.5, /, /)
endoutput
on the closed interval (0.0,298.0) we h**e :
computed integral = 64.82892
exact integral = 64.8
error = -0.02892
查閱crc化學和物理手冊,可得到固體鉛298k時的絕對熵為64.8,與計算結果之間的差別很小。
本題可用多種方法來解決,如simpson數值積分法,或採用的形式,或先用資料點擬合方程得回歸引數,然後對回歸方程進行積分(回歸方程可用,,也可在此基礎上探索更好的回歸方程形式)。
MATLAB及其應用實驗指導書
一 實驗目的 學習matlab軟體的安裝過程,熟悉matlab軟體介面的組成及基本使用方法。理解陣列 array 的分類,及標量 scalar 向量 vector 和矩陣 matrix 的區別,熟悉陣列與矩陣的構造方法,掌握陣列與矩陣的基本運算法則。二 實驗要求 1 掌握matlab軟體的啟動與退出...
MATLAB在級數的應用
matlab語言 課程 姓名 楊小花 學號 12010245315 專業 通訊工程 班級 通訊工程 1 班 指導老師 湯全武 學院 物理電氣資訊學院 完成日期 2011.12.10 matlab在級數中的應用 楊小花 12010245315 2010級通訊班 摘要 matlab除了對數值計算外,還有...
實驗一Matlab軟體的使用
1 實驗目的和要求 掌握matlab軟體的使用方法 常用離散時間訊號的產生 顯示和運算。2 實驗內容 matlab軟體常用命令和工具箱的使用,學會簡單的程式設計 程式設計實現常用離散時間訊號 離散時間訊號的疊加 移位 線性卷積等基本運算。3 實驗原理 利用軟體生成數字訊號處理系統中所涉及的訊號及訊號...