APPLICATIONS OF MATLAB ON ENGINEERING PROBLEMS

12/08/2006

林書如：學生期中考成績分佈情形

因為脫離物理這門學科太久了
所以老師網頁上的題目都讓我覺得很陌生....
因此，我是以網頁內跟統計較相關的程式"度數分配程式"做了一些修正以及增加一些小指令，以敘述統計方式呈現學生期中考成績分佈情形。

12/01/2006

http://r95622024.blogspot.com/2006/11/blog-post_2921.htm

前言:
常常在攝影的人都知道
用單眼數位拍出來的照片通常需要 “後製”
也就是對照片作修改

因為常常用到像photoimpact或photoshop等的編輯軟體
突然想到可不可以利用
Matlab的強大功能作類似的處理

經過摸索之後找到可以把照片載入matlab的方法
並進行幾項簡單的修改

介紹:

在這裡我選擇了幾項簡單的功能並製作了一個選單
基本的功能包括
1. 觀看原圖 View photo
2. 顏色濾鏡功能 RGB filter
3. 反白 Counter image
4. 彩色照片轉黑白照 Color ro B/W

功能說明：
1. View photo
觀看指定的照片並有放大縮小的功能

2. RGB filter
就像photoshop的濾鏡功能，分成R、G、B 濾鏡
點選之後會輸入三個數值，r、g、b 分別是0~1的數值
1 代表顏色原始呈現 0 代表完全被遮住。

3. Counter image
也就是反白功能，能將明暗對調

4. Color ro B/W
把彩色照片轉換成黑白照片

使用方法：
把想處理的照片 xxx.xxx (如 w.jpg) 放在matlab的 Current Directory。
輸入指令 : Photo
打入 xxx. xx (記得打副檔名)
就進入選單模式

範例：
選這張照片進行處理 “ w.jpg “

5. (1) View photo
會顯示原照片

 http://photos1.blogger.com/x/blogger2/2907/859880749988154/1600/963015/w.jpg

6. (2) RGB filter

R=1
G=0
B=0

http://photos1.blogger.com/x/blogger2/2907/859880749988154/1600/36517/r.jpg

R=0
G=1
B=0

http://photos1.blogger.com/x/blogger2/2907/859880749988154/1600/388570/g.jpg
R=0.5
G=0.5
B=0.5
http://photos1.blogger.com/x/blogger2/2907/859880749988154/1600/500246/05.jpg

R=0
G=1
B=1
http://photos1.blogger.com/x/blogger2/2907/859880749988154/1600/859043/bg.jpg

7. Counter image

http://photos1.blogger.com/x/blogger2/2907/859880749988154/1600/872523/count.jpg

8. Color ro B/W
http://photos1.blogger.com/x/blogger2/2907/859880749988154/1600/673361/gra.jpg

結語：
目前只是最基本的修改功能
未來計畫要加入的包括進接的修改功能及目前富士最新的臉部辨識功能

程式碼

file=input('select your photo (xxx.xxx)\n:','s');
w=imread(file);
k=menu(' Matphoto v.1 by Ching','View photo','RGB filter','Counter image','Color ro B/W');
switch k
case 1
imview(w);
case 2
r=input('Input your R value(0~1)\n:');
g=input('Input your G value(0~1)\n:');
b=input('Input your B value(0~1)\n:');
w2(:,:,1)=w(:,:,1)*r;
w2(:,:,2)=w(:,:,2)*g;
w2(:,:,3)=w(:,:,3)*b;
imshow(w2);
case 3
w2=256-w;
imshow(w2);
case 4
w2=(w(:,:,1)+w(:,:,2)+w(:,:,3))/3;
imshow(w2);
end

熱傳模組建立（續）

繼上次的穩態一般邊界
這次加入穩態所有邊界情形

11/28/2006

小程式月曆產生器


把程式加入介面


  GUIDE程式說明
    首先在command window下打guide就會出現下圖

  接下來選Blank GUI(Default)，就會開啟新的window，如下圖

接下來你會看到左邊有很多控制物件，
你只要以滑鼠點選其中一項
並把它拖到你想要的位置就可以.
    月曆這個程式我有用到[push button]二個(run 和run)、
九個[text](如萬年曆、月、年和藍底的框都是)，
而其它的空白框則都是[edit](可以key in數字的)。
如下圖
(小格子是中文字，因為我7版的不能用，所以我就用6.5版的
，在此要注意，應該會有些指令在6.5不能用)。


而要是想更改其顏色，字體等…則在其物件點擊二下則會出現下圖(這是以紅色的run為例子)，
或是mark其物件，選View，再選Property Inspector也可以
  
當你選擇好你要用的物件，即可存檔。此時你會存成二個檔
一個叫 xxx.fig的檔，另一個是xxx.m的檔。接下來則是程式編輯的部份了,
點選View，再選M-file Editor即可開始編輯程式。
當你一開啟M-file你就會看到密密麻麻的字，但在這個程式會用到的只有幾個而已。
      大部份我們都是在function 物件名稱_callback(hObject, eventdata, handles)
下編輯其程式，如下圖：
這個程式的[edit物件名稱]空白處的程式碼為
function  edit物件名稱_Callback(hObject, eventdata, handles)
% hObject    handle to edit8 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
newstrval=get(hObject,'string');  % 截取輸入的字串
newval=str2double(newstrval);  % 變成數字
handles.edit物件名稱=newval;  % 存儲起來
guidata(hObject,handles);
% Hints: get(hObject,'String') returns contents of edit8 as text
%        str2double(get(hObject,'String')) returns contents of edit8 as a double
以下為pushbutton的圖
  這個程式的[pushbutton物件名稱]的程式碼為  
Function  pushbutton物件名稱_Callback(hObject, eventdata, handles)  
% hObject    handle to pushbutton3 (see GCBO)  
% eventdata  reserved - to be defined in a future version of MATLAB  
% handles    structure with handles and user data (see GUIDATA)  
Y=handles.edit6; % 截取年的值
M=handles.edit7; % 截取月的值
K=handles.edit8; % 截取日的值  
% 以下為判斷你輸入的日子是星期幾  
if rem(Y,100)==0 & rem(Y,400)==0 | rem(Y,4)==0
    d1=[31 29 31 30 31 30 31 31 30 31 30 31];
else
    d1=[31 28 31 30 31 30 31 31 30 31 30 31];
end
s=Y-1;
c=d1(1:(M-1));
c1=sum(c)+K;
ss=s+fix(s/4)-fix(s/100)+fix(s/400)+c1;
dd=rem(ss,7)+1;
di=['日','一','二','三','四','五','六']; 
ggo=sprintf('你查詢的%d年%d月%d日是星期%s',Y,M,K,di(dd)); 
set(handles.text12,'string',ggo);  % 把它show出來在藍色框框中如下圖


另一個[pushbutton物件名稱]的程式碼為
% --- Executes on button press in pushbutton4.
function  pushbutton4_Callback(hObject, eventdata, handles)
% hObject    handle to pushbutton4 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
Y=handles.edit9; % 讀取年的值
M=handles.edit10; % 讀取用的值
ans=[ ];
if rem(Y,100)==0 & rem(Y,400)==0 | rem(Y,4)==0
    d1=[31 29 31 30 31 30 31 31 30 31 30 31];
else
    d1=[31 28 31 30 31 30 31 31 30 31 30 31];
end
s=Y-1;
c=d1(1:(M-1));
% 以下為建月曆的運算式
for k=1:d1(M)
    c1=sum(c)+k;
    ss=s+fix(s/4)-fix(s/100)+fix(s/400)+c1;
    ans(k)=rem(ss,7);        
end
ans(ans==0)=7;
pp=ans(1);
ans1=[ ];
k1=0;
for n=1:6;
    for m1=pp:7
        ans1(n,m1)=k1;
        k1=k1+1;    
    end
    pp=1;
end
ans1(ans1>d1(M))=0; % 結果
 % 以下就是把所出來的值給代入小格子
(這個東東，還可以修改)  
set(handles.e1,'string','日');
set(handles.e2,'string','一');
set(handles.e3,'string','二');
set(handles.e4,'string','三');
set(handles.e5,'string','四');
set(handles.e6,'string','五');
set(handles.e7,'string','六');  
set(handles.e8,'string',num2str(ans1(1,1)));
set(handles.e9,'string',num2str(ans1(1,2)));
set(handles.e10,'string',num2str(ans1(1,3)));
set(handles.e11,'string',num2str(ans1(1,4)));
set(handles.e12,'string',num2str(ans1(1,5)));
set(handles.e13,'string',num2str(ans1(1,6)));
set(handles.e14,'string',num2str(ans1(1,7)));  
set(handles.e15,'string',num2str(ans1(2,1)));
set(handles.e16,'string',num2str(ans1(2,2)));
set(handles.e17,'string',num2str(ans1(2,3)));
set(handles.e18,'string',num2str(ans1(2,4)));
set(handles.e19,'string',num2str(ans1(2,5)));
set(handles.e20,'string',num2str(ans1(2,6)));
set(handles.e21,'string',num2str(ans1(2,7)));  
set(handles.e22,'string',num2str(ans1(3,1)));
set(handles.e23,'string',num2str(ans1(3,2)));
set(handles.e24,'string',num2str(ans1(3,3)));
set(handles.e25,'string',num2str(ans1(3,4)));
set(handles.e26,'string',num2str(ans1(3,5)));
set(handles.e27,'string',num2str(ans1(3,6)));
set(handles.e28,'string',num2str(ans1(3,7))); 
set(handles.e29,'string',num2str(ans1(4,1)));
set(handles.e30,'string',num2str(ans1(4,2)));
set(handles.e31,'string',num2str(ans1(4,3)));
set(handles.e32,'string',num2str(ans1(4,4)));
set(handles.e33,'string',num2str(ans1(4,5)));
set(handles.e34,'string',num2str(ans1(4,6)));
set(handles.e35,'string',num2str(ans1(4,7)));  
set(handles.e36,'string',num2str(ans1(5,1)));
set(handles.e37,'string',num2str(ans1(5,2)));
set(handles.e38,'string',num2str(ans1(5,3)));
set(handles.e39,'string',num2str(ans1(5,4)));
set(handles.e40,'string',num2str(ans1(5,5)));
set(handles.e41,'string',num2str(ans1(5,6)));
set(handles.e42,'string',num2str(ans1(5,7)));  
set(handles.e43,'string',num2str(ans1(6,1)));
set(handles.e44,'string',num2str(ans1(6,2)));
set(handles.e45,'string',num2str(ans1(6,3)));
set(handles.e46,'string',num2str(ans1(6,4)));
set(handles.e47,'string',num2str(ans1(6,5)));
set(handles.e48,'string',num2str(ans1(6,6)));
set(handles.e49,'string',num2str(ans1(6,7)));    
  下面為顯示結果圖

等程式都key in完，就可以按run，即可執行程式
月曆產生器檔案下載

11/24/2006

熱傳模組建立

this a test

11/13/2006

簡單的濾波器製作

　　對於擷取肌電圖（EMG）時，擷取到的資料對於我們而言
並不是都是有意義的，有時候常常會收取到一些雜訊，例如直
流電所產生的60Hz的雜訊，或是其實我們收取EMG的訊號以取
樣頻率1000為例，大多介於20-500之間。
　　因此濾波對於我們肌電訊號的分析格外重要。這裡撰寫了
一個非常簡單的濾波功能。有興趣的話可以點我看看。或是按
標題也可以連結。

11/11/2006

請問您有關MATLAB的問題楊淑媛

老師您好

抱歉打擾您了

想請問您有關"比較資料作累加動作:的問題

data2是k*3的矩陣

現在的比較方式是

當第三行的其中一個數比另一個數小

其對應的第二行資料卻比另一個數對應的第二行資料大時

則p+1

相反的

當第三行的其中一個數比另一個數大

其對應的第二行資料也比另一個數對應的第二行資料大時

則m+1

附上的是我寫的程式

基本的debug是沒錯

但它跑很久

(雖然說k=1000)

因此我怕這是一個無窮迴圈

請老師幫我看一下

謝謝老師

抱歉打擾您

麻煩您了

學生淑媛

11/03/2006

黃世榮:經緯儀的定位

前言:早前出差去使用經緯儀，希望能用兩種角度去作定位
問題:如何在只有兩組經緯儀的角度下作空間點群的定位，在建立公式以及程式後，發覺程式當時剛學技巧有需要改進且似乎這樣的計算是有錯的，希望大家幫忙想
程式流程:
clear
format long
matrix=xlsread('diow.xls','sheet1','A2:L3'); %% read the data
[m,n]=size(matrix); %% get row and column
L=input('Input the L you get==> ')
%%request the distance between two measure point

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% to build the alpha's array
for i=1:m
A(i)=matrix(i,1)
%% to build the alpha's array
A2(i)=matrix(i,2)
%% to build the alpha's minute array
A3(i)=matrix(i,3)
%% to build the alpha's second array

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% to build the theta's array
S(i)=matrix(i,4)
%% to build the theta's array
S2(i)=matrix(i,5)
%% to build the theta's minute array
S3(i)=matrix(i,6)
%% to build the theta's second array

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% to build the Beta's array
B(i)=matrix(i,7)
%% to build the Beta's array
B2(i)=matrix(i,8)
%% to build the Beta's minute array
B3(i)=matrix(i,9)
%% to build the Beta's second array

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% to build the PHI's array

F(i)=matrix(i,10)
%% to build the PHI's array
F2(i)=matrix(i,11)
%% to build the PHI's minute array
F3(i)=matrix(i,12)
%% to build the PHI's second array
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% to make their unit coherence
for i=1:m
A(i)=A(i)+(A2(i)/60)+(A3(i)/3600)
S(i)=S(i)+(S2(i)/60)+(S3(i)/3600)
B(i)=B(i)+(B2(i)/60)+(B3(i)/3600)
F(i)=F(i)+(F2(i)/60)+(F3(i)/3600)
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%find that limit it or not is no matter
for i=1:m
% if B(i)<90
C(i,1)=L*sin(B(i)/180*pi)/sin((B(i)-A(i))/180*pi)
% to build X's array (X,Y are not the real coordinate)
C(i,2)=L*sin(A(i)/180*pi)/sin((B(i)-A(i))/180*pi)
% to build Y's array
% end
% if B(i)>90
% C(i,1)=L*sin(B(i)/180*pi)/sin((B(i)+A(i))/180*pi)
% to build X's array
% C(i,2)=L*sin(A(i)/180*pi)/sin((B(i)+A(i))/180*pi)
% to build Y's array
% end
end

%% to build coordinates' matrix
for i=1:m
H(i,1)=C(i,1)*sin(A(i)/180*pi)
H(i,2)=C(i,1)*cos(A(i)/180*pi)
H(i,3)=C(i,1)*tan(S(i)/180*pi)
H(i,4)=C(i,2)*sin(B(i)/180*pi)
H(i,5)=C(i,2)*cos(B(i)/180*pi)+L
H(i,6)=C(i,2)*tan(F(i)/180*pi)
end

%output the data to excel
xlswrite('test',H,'sheet1','A2');

for i=1:m
plot3(H(i,1),H(i,2),H(i,3),'o');
plot3(H(i,4),H(i,5),H(i,6),'o');
hold on
end

吳子青: 試作搜尋引擎( I )

前言：
本來構想是希望作一個可以將老師部落格上面
關於Matlab程式碼斷行問題解決的程式

利用copy指令將部落格的程式碼複製下來
會發現變成一整行老師在上課的時候也不太方便
構想:
因為學到cell() 本來希望將整個文件按照字元存成
cell, 再去對矩陣作處理, 把斷行點找出來之後
重新輸出成正確的文件
但是...
找不到可以斷行的依據

所以轉換方向作一個網頁或是文件的搜尋器 (google..??)
目前進度緩慢，所以分階段進行

問題
如何將字元讀出並進行搜尋動作？

----------------------------------------------------------
答案:
題旨與分析
將文件讀入並將文件裡面的每個字元存入記憶體宣告成矩陣
(或者文件太大以動態搜尋方式不先將整個文件全部存入記憶體）
對矩陣作搜尋動作, 找出關鍵字對應的位置
將搜尋到的字元或位置進行後置處理
(列出一整行,或者標示不同顏色？）
作出使用者介面

程式流程
因為尚未完成所以先簡述如下

/tra.m/
[fid,message]=fopen('test2.txt','r');
mydata=fread(fid);
if fid==-1
disp(message);
end
fclose all;

%如此就將文件文字全部以二進位儲存方式
%儲存成矩陣mydata

/dis.m/
disp(char(mydata'))

/run.m/
[m,n]=find(mydata==keyword);
%找出矩陣中的keyword的二進位元,在mydata中作搜尋
找出對應的位置
A={mydata(m,n-2) mydata(m,n-1) mydata(m,n) mydata(m,n+1) mydata(m,n+2) };
disp(char(A));
%將關鍵字元的二進位碼轉換, 並顯示前兩字後兩字

執行結果
test2.txt 內容 ghuodjgpobinhindpohinjestpohdbvuisneptogindbuodb
關鍵字p
結果 :
jgpob , ndpoh , stpoh , nepto

討論
這目前只是初步的概念 ,只是針對字串中的字母搜尋關鍵字母並作處理
未來要克服的包括關鍵字串的設定, 文件檔案讀取更general, 及加入
if 等規則等等等

吳堉光：EMG之Onset時間點分析

物理治療研究中，肌電圖（electromyography）是非常重要的研究工具之一，以下是我利用所教的while loop找出肌電訊號開始的時間點
EMG時間點

11/02/2006

林書如：您的體重過重嗎

前言：
BMI指數可以反映出體重是否處於理想狀態。過去已有許多文獻指出BMI指數過高往往是造成糖尿病、高血壓，甚至是心血管疾病的危險因子。因此，適度地監控自己的BMI指數是否處於正常值內是很重要的。

問題：
試寫出一個函數，當輸入體重和身高的數值後，能計算出其身體質量指數BMI值。並且依據其BMI值判斷其是否有體重過重或過輕的問題。(18.5<=BMI<24時顯示"體重標準"，BMI>24時顯示"體重過重"，BMI<18.5時顯示"體重過輕")

指令：

一開始我寫了一個比較簡單敘述檔案(Script file)

function [A]=BMI(w,h)
% BMI(Body Mass Index) calculation
A=w/h^2; % A means BMI
if 18.5<=A & A<24
fprintf('BMI=%0.2f 體重標準\n',A)
elseif A<18.5
fprintf('BMI=%0.2f 體重過輕\n',A)
else
fprintf('BMI=%0.2f 體重過重\n',A)
end

模擬輸出結果：

>> BMI(60,1.5)
BMI=26.67 體重過重

ans =

26.6667

結果顯示，一位體重60公斤，身高150公分的人，其BMI指數為26.67，有體重過重的傾向

==============================================

但是，為了能夠加上交談的功能，讓使用者更能清楚地了解應該輸入的內容
例如，輸入之體重應以公斤為單位，輸入之身高應以公尺為單位
因此我又改寫成另一個函數檔(function file)，如下

% BMI2.m
% BMI(Body Mass Index) calculation
a=input('Please keyin your body weight? (kg)');
if isnumeric(a)
b=input('Please keyin your body height? (m)');
B=a/b^2; % B means BMI
if 18.5<=B & B<24
fprintf('BMI=%0.2f 體重標準\n',B)
elseif B<18.5
fprintf('BMI=%0.2f 體重過輕\n',B)
else
fprintf('BMI=%0.2f 體重過重\n',B)
end
end

模擬輸出結果：

>> BMI2
Please keyin your body weight? (kg)60
Please keyin your body height? (m)1.5
BMI=26.67 體重過重

結果同上，一位體重60公斤，身高150公分的人，其BMI指數為26.67，有體重過重的傾向

討論：

我是以簡單的if-elseif指令來撰寫這兩個小程式
程式中如果還有可以補強的地方(例如發現程式的bias，或是其實還有更簡單的寫法)，還麻煩大家提供意見 :p

10/25/2006

吳子青: 民意調查真的公正嗎

前言：
最近又到了選舉期了在電視媒體或報章雜誌常常可以看到
“ 某某參選人支持率多少”…
民意支持率各家都宣稱客觀公正，但其實是可以操作的
在這裡舉出一個以前碰過的問題讓大家思考一下　

問題：

若台北市某候選人的真實 “民意支持率”為15%，找出使調查樣本的
民意支持率超過0.18 的出現機率不大於 10% 的最小樣本數 N。

----------------------------------------------------------------------
答案：

題旨與分析
這問題的思考流程是這樣的
將全部台北市民當成一個族群，而調查的對象代表抽樣的樣本
“某候選人的真實民意支持率為15%” 這句話的意思就是代表族群的支持率
而調查樣本的民意支持率就是樣本的支持率。

將這問題視為一個二項式分佈的問題，也就是支持和不支持兩種情況
支持的機率 P=0.15 不支持的機率 Q=0.85
舉例說明
如果某一次的抽樣從台北市民隨機抽10個人 (Sample Size m=10)
10個人中支持候選的人有2個人，不支持的有8個人
這個事件出現的機率為

 http://photos1.blogger.com/blogger2/2907/859880749988154/1600/formula1.jpg

所以如果問題是 ”找出使調查樣本的民意支持率超過0.18”
如果抽樣人數維持10人，支持人數大於1.8的都符合條件
也就是支持人數 2 到10人分別算出的機率加總。

題目所求是 “出現機率不大於 10% 的最小樣本數 N。”所以從n=1
開始算到符合出現機率不大於10％為止。

程式流程
Step 1 : 宣告function p1(m,n) 為抽樣m個人，有n個人支持的機率

Step 2 : 計算A為抽樣m個人，支持率超過0.18的機率總和

Step 3 : 從m=1 開始增加用For 跑迴圈，決定停止條件為A < >程式內容


/ p1.m / (* 這行代表m-files p1.n 開始 *)

function output=p1(m,n)

if n>0.5*m
a=m-n;
output=prod((m-a+1):m)/prod(1:a)*(0.15)^n*(0.85)^(m-n);

else
output=prod((m-n+1):m)/prod(1:n)*(0.15)^n*(0.85)^(m-n);
end

/ pp1.m / (* 這行代表m-files pp1.n 開始 *)
clear;
max=400;
ppp1=zeros(max,2);

for m=1:max
k=ceil(m*0.18);
A=0;

if A<=0.1 fprintf('最小樣本數%d\t出現機率%d',m,A);
break
end
fprintf('最小樣本數%d\t出現機率%d\n',m,A);
clc;
end
plot (ppp1(1:max,2), 'DisplayName', 'ppp1(1:max,2)', 'YDataSource', 'ppp1(1:max,2)');
figure(gcf)
xlabel('樣本數目');
ylabel('出現機率');
/ Command windows / (* 這行代表command windows輸入 *)
pp1

執行結果
最小樣本數217 出現機率9.559103e-002
這是最後跑出來的結果
N=217為民意支持率超過0.18 的出現機率不大於 10% 的最小樣本數 N

這是程式結果

http://photos1.blogger.com/blogger2/2907/859880749988154/1600/hy2.jpg

這是圖

http://photos1.blogger.com/blogger2/2907/859880749988154/1600/hydro.jpg

討論
(1)因為老師上週的上課內容剛交到For 和If的用法，所以在這裡練習一下

(2) 程式中碰到最大的問題在於 prod () 連乘函數括號內的矩陣不能太大
從1＊2＊3……..＊171 ＝prod(1:171) 這樣結果會爆掉，可是在算二項式分佈
的機率時就一定會用到。

所以在這裡用了一個小技巧，在p1裡當 n > 0.5 m 時令n=m-n
因為在排列組合裡 C10取2=C10取8，10!/8!2! =10!/2!8!
在這裡可以避免連乘太多項，可以有效避免算到M=171時，
結果一定會出現NaN。但是當m 算到300多時
仍然會出現同樣的問題。或許會有更有效率的寫法。

10/19/2006

PROBLEM #4.3

氣體在特定容器下，其壓力、體積與溫度間之關係如下：

P=nRT/V

這是波義耳的氣體公式，但萬得爾(van der Vaals)公式則加入分子間之引力因素。其關係如下：

P=nRT/(V-nb)-an²/V²

其中，a,b分別為氣體之特定常數。理想氣體在n=1摩爾數下，溫度在0C(273.2K)下，氣壓為一大氣壓時，具有體積為22.4公升，其Ｒ值為0.0826。對於氯氣而言，上述之a,b值分別為a=6.49，b=0.0562。試就上兩公式計算其壓力差。

PROBLEM #4.2

若採用定存的方式，每個月存放一定的數額（例如每個月五百元），銀行年息5%。

1. 問何時可以達到10,000,000元？

2. 若改為每週存入原每月存放額度之四分之一，則何時可以達到10,000,000元？

PROBLEM #4.1

有一個很好玩的數字遊戯是：假設第一天到銀行存一分錢，第二天存二分錢，第三天三分錢，如此累積，問何時可以累積到10,000,000元？（設銀行的年息為5%）

撰寫一程式，試算上述之過程，並將其每達1000元倍數之日數同時印出，直到達到所期望之金額。
若每天所存改為五分錢、十分錢、十五分錢，。。。，試比較達到期望金額之日數。
實際生活上能否執行此項累積方式，試說明你的想法。

10/11/2006

PROBLEM #3.4：九九乘法表

請使用最簡單之指令，製作一個９X９之九九乘法表，其中之項目格式為 3 x 4 = 12　3 x5 = 15 ...。

10/10/2006

PROBLEM #3.3：平均氣溫

五月份的台北天氣已經開始轉熱，其整月的日平均溫度如下：

T=[25 26 24 23, 18 17 20 30 29, 25 19 19 21 20, 23 24 23 25 23,…
27 28 24 26 27, 23 25 24 21 19, 20]C

1. 試撰寫一程式，以上述之溫度資料為輸入值，求其月平均、最高及最低溫度。

2. 試撰寫另一程式，以交談式輸入之方式，按月份輸入每日之溫度，（必須分清楚每月應有之輸入日數）然後呼叫第一項之程式進行分析。

註：請同學針對本題目直接回應(COMMENT)，請勿採用EMAIL方式提交。

PROBLEM #3.2：抛物線運動

在抛物運動中，設初值分別為速度ν及仰角θ，試其投射之最高高度、最遠距離。已知ν= [ 10 15 20 25 ] m/s，其對應θ=[30 40 50 60]度。
1. 試撰寫一函數，並繪出其投射之軌跡。
2. 同上題，若發射點高於平面100m時，求其對應之狀況。
3. 同上題，求物體自發射後，其著地之時間。

註：請同學針對本題目直接回應(COMMENT)，請勿採用EMAIL方式提交。

PROBLEM #3.1：三角函數表

試做一個三角函數表，第一欄為度數，其區間由0-90度，每間隔為一度，第二欄為正弦值(SINE)，第三欄為正切值(TANGENT)，兩者之有效數字為小數點後四位數。
針對第一項主題撰寫一函數檔案，可以隨時改變度數之範圍及小數點之位數。

註：請同學針對本題目直接回應(COMMENT)，請勿採用EMAIL方式提交。

10/07/2006

problem#2林冠宏

第一題
題意:
(A)把B矩陣單位由華氏溫度轉為攝氏溫度
(B)把溫度超過(>38℃)定義為"1",而沒有超過(<38℃)定義為"0"
分析:
(A)先定義矩陣B,把B用公式轉換成攝氏溫度
(B)把換算出來的攝氏溫度,再和38℃比較,再定義超過(>38℃)定義為"1",而沒有超過
(<38℃)定義為"0"
程式流程:
(A)>>B=[98 100 102 104 98.4 98.2 98.5 101 102 99.5]; %Their temperatures
measured
>>C=(B-32)*5/9 %Converse the unit of B into Centigrade

C =

36.6667 37.7778 38.8889 40.0000 36.8889 36.7778 36.9444
38.3333 38.8889 37.5000
(B)>> D=C>38 %design Dto show that "1" for the one who has a fever (>38℃
) and "0" for the normal

D =

0 0 1 1 0 0 0 1 1 0
結果:
(A)C =36.6667 37.7778 38.8889 40.0000 36.8889 36.7778 36.9444
38.3333 38.8889 37.5000
(B)D =0 0 1 1 0 0 0 1 1 0
討論:
若是在每個矩陣中的數字可以定義其所代表的名字,若是其中有人超過,就可以立即知道
其名
第二題
題意:
比較A、B矩陣的相關:(A)D=(A>B)(B)D=(A>5)(C)D=A+B
分析:
比較其矩陣的各項和題意的相互關係
程式流程:
(A)>> A=[-10 8 6 4 -5 20];
>>B=[2 8 5 10 -6 3];
>> D=(A>B) %當A大於B時,定義當A>B為真時,表為"1",為假時表為"0"

D =

0 0 1 0 1 1

(B)>> D=(A>5) %當A大於5時,定義當A>B為真時,表為"1",為假時表為"0

D =

0 1 1 0 0 1

(C)>> D=A+B %把A、B矩陣相加

D =

-8 16 11 14 -11 23
結果:
(A)D =0 0 1 0 1 1
(B) D =0 1 1 0 0 1
(C) D =-8 16 11 14 -11 23
討論:
學會如何使用比大小
第三題
題意:
做矩陣加減乘除法和三角函數運算
分析:
先確定矩陣的大小是否符合加減乘除所須之維度,在做運算\
程式流程:
>> x=[10 20 30];
>>y=[1 4 6];
>> A=3*x+y

A =

31 64 96

>> B=5*y./x

B =

0.5000 1.0000 1.0000

>> C=(4*x).*(2*y)

C =

80 640 1440

>> D=sin(x).*cos(y)

D =

-0.2939 -0.5967 -0.9487 %三角函數以逕度量進行運算

>> E=(5*x).*sin(2*y)

E =

45.4649 98.9358 -80.4859

討論:是矩陣中的單元對應的相加減乘除,而不是矩陣的乘除,所以要注意其"."位子
第四題:
題意:在M-file中做運算,紀錄其工作列,在運算視窗中直接叫出
分析:按File後New後M-file在其中做串聯和並聯電阻的運算,在紀錄在Work裡,在Matlab
中把檔案叫出
程式流程:
1- R=[10 30 200 400];
2- series=sum(R)%串聯電阻
3- parallel=1/(sum(1./R))%並聯電阻
(紀錄在Untitled.m)
>> untitled

series =

640

parallel =

7.1006
討論:在M-file中做運算,若是非常大的運算,是非常好用的方法
第五題:
題意:把一個A矩陣(1-24)變成3*8,6*4,2*12矩陣大小
分析:使用reshape(X,Y,Z)
結果:
>> A=1:24

A =

Columns 1 through 13

1 2 3 4 5 6 7 8 9 10 11 12 13

Columns 14 through 24

14 15 16 17 18 19 20 21 22 23 24

>> B=reshape(A,3,8)

B =

1 4 7 10 13 16 19 22
2 5 8 11 14 17 20 23
3 6 9 12 15 18 21 24

>> C=reshape(A,6,4)

C =

1 7 13 19
2 8 14 20
3 9 15 21
4 10 16 22
5 11 17 23
6 12 18 24

>> D=reshape(A,2,12)

D =

1 3 5 7 9 11 13 15 17 19 21 23
2 4 6 8 10 12 14 16 18 20 22 24

_________________________________________________________________
MSN 交友和 Match.com 為您線上尋找您週圍的單身貴族！

http://match.tw.msn.com/match/mt.cfm?pg=join&tcid=281200

訂閱：文章 (Atom)