2016-11-26 9 views
1

私は現在、分散ロードビームを与えられたMatlabの問題に取り組んでおり、数値的に(複合台形ルールを使用して)積分してせん断力と曲げモーメントを求めています。そこから、私は最大の瞬時値とこれが起こる位置を見つける。それから私は再び2倍に分けてせん断力と分散荷重を求めます。Matlabで数値的に統合されています。なぜ部門が増えると正確な出力が低下するのですか?

私は、数値積分と微分が正しく計算されていることを確認し、確認するために分析的に統合された関数を使用しています。私の問題は、計算を実行するのに13ポイントしか使用しない場合(x = 0:12;)、私の値は解析値の結果に近い値になるということです。 (この数字は、梁が12フィートであり、師団がそれぞれ1フィートであることに由来する)。私の数値積分の精度を上げるために分割数を増やすと、値は分析値からさらに遠ざかります。私はこれを引き起こした可能性があることを理解したいと思います。これまでのところ成功していません。

マイコード:2番目と3番目の行上記のコードで

clear; clc; 
% n=100; 
% x = linspace(0,12,n+1);   %dx (100 divisions) 
x=0:12;        %dx (12 divisions) 
w = 12.5.*x;       %distributed load 

%//Integrated Shear 
V(1)=300; 
for i = 2: length(w) 
    weight = [.5 ones(1,i-2) .5]; 
    V(i)=300-sum(weight.*w(1:i)); 
end 
figure(1); clf; 
plot(x,V,'--'); 
legend('Numerical Shear'); 
xlabel('Position (ft)'); ylabel('Shear Force (lb)'); 

%//Integrated Moment 
Mactual = 300.*x - ((25/12).*(x.^3));   %analytically integrated M 
M(1)=0; 
for i = 2: length(V) 
    weight = [.5 ones(1,i-2) .5]; 
    M(i)=sum(weight.*V(1:i)); 
end 
figure(2); clf; 
plot(x,M,'--'); hold on; 
plot(x,Mactual); hold off; 
legend('Numerical Moment' , 'Analytical Moment'); 
xlabel('Position (ft)'); ylabel('Bending Moment (lbf)'); 

%//Max Moment and Position at Max Moment 
[maxMValue, indexAtMaxM] = max(M); 
xValueAtMaxMValue = x(indexAtMaxM(1)); 
dispM = [ 'Maximum bending moment: ' , num2str(maxMValue) , ' lbf']; 
disp(dispM); 
dispX = [ 'Position at maximum bending moment: ' , num2str(xValueAtMaxMValue) , ' ft' ]; 
disp(dispX) 

%//Derived Shear 
dM(1)=300; 
for i = 1:length(x)-1   
    dM(i+1) = (M(i+1)-M(i));  
end 
figure(3); clf; 
plot(x,dM,'--'); hold on; 
Vactual = 300 - 6.25.*(x.^2);   %analytically integrated shear 
plot(x,Vactual); hold off; 
legend('Numerical Shear' , 'Analytical Shear'); 
xlabel('Position (ft)'); ylabel('Shear Force (lb)'); 

%//Derived Load 
dV(1)=0; 
for i = 1:length(x)-1   
    dV(i+1) = -(V(i+1)-V(i));  
end 
figure(4); clf; 
plot(x,dV,'--'); hold on; 
plot (x,w); hold off; 
legend('Numerical Load' , 'Analytical Load'); 
xlabel('Position (ft)'); ylabel('Distributed Load (lb/ft)'); 

は12部門で出力を表示するようにコメントアウトされています。これらの2行のコメントを外して4行目をコメントアウトすると、出力に100桁の値が表示されます。

クイックメモ:分析値は、積分モーメント、派生せん断、および派生分散荷重を比較としてプロットされます。分析値(MactualとVactual)を正しい出力値と見なしてください。

私はこの問題を解決する方法につながるどんなインプットも非常に高く評価されます。

答えて

0

正しい整数値を取得するには、n12と異なる場合dx=12/nまたは同一のdx=x(2)-x(1)で関数値の合計を乗算する必要があります。


(11/28)また、台形積分は、このように各工程における合計の重量と計算の再構築を回避

V(1)=300 
M(1)=0 
for i=2:n 
    V(i) = V(i-1) - 0.5*dx*(w(i-1)+w(i)) 
    M(i) = M(i-1) + 0.5*dx*(V(i-1)+V(i)) 
end 

ようにループを簡素化することができます。これにより、アルゴリズムのこの部分がO(n^2)からO(n)に減少します。

+0

私はこれを試して、動作させようとしましたが、2番目の積分(モーメント)と派生した両方の値の出力を歪ませます。私は、瞬間の積分のために私のループに別の問題があると信じていますが、これはまだ解決には至りません。 –

+0

私は 'sum'を' dx * sum'に変更しました。もちろん、あなたはまた、実際の商として差の商を計算しなければなりません、 'dM(i + 1)=(M(i + 1)-M(i))/ dx'と' dV(i + (V(i + 1)-V(i))/ dx; – LutzL