我在一个数值分析课上,我正在做一个家庭作业问题。这是Timothy Sauer的数值分析,在第二个建议活动部分。我一直在和我的教授谈论这个代码,似乎错误和近似是错误的,但我们都不能找出原因。下面的代码是我正在使用的,这是在MatLab中。有谁对欧拉·伯努利梁有足够的了解,有谁能帮上忙?
function ebbeamerror %This is for part three
disp(['tabe of errors at x=L for each n'])
disp([' n ',' Aprox ',' Actual value',' Error'])
disp(['======================================================='])
format bank
for k = 1:11
n = 10*(2^k);
D = sparse(1:n,1:n,6*ones(1,n),n,n);
G = sparse(2:n,1:n-1,-4*ones(1,n-1),n,n);
F = sparse(3:n,1:n-2,ones(1,n-2),n,n);
S = G+D+F+G'+F';
S(1,1) = 16;
S(1,2) = -9;
S(1,3) = 8/3;
S(1,4) = -1/4;
S(2,1) = -4;
S(2,2) = 6;
S(2,3) = -4;
S(2,4) = 1;
S(n-1,n-3)=16/17;
S(n-1,n-2)=-60/17;
S(n-1,n-1)=72/17;
S(n-1,n)=-28/17;
S(n,n-3)=-12/17;
S(n,n-2)=96/17;
S(n,n-1)=-156/17;
S(n,n)=72/17;
E = 1.3e10;
w = 0.3;
d = 0.03;
I = w*d^3/12;
g = -9.81;
f = 480*d*g*w;
h = 2/10;
L = 2;
x = (h^4)*f/(E*I);
x1 = ones(n ,1);
b = x*x1;
size (S);
size(b);
pause
y = S\b;
x=2;
a = (f/(24*E*I))*(x^2)*(x^2-4*L*x+6*L^2);
disp([n y(n) a abs(y(n)-a)])
end
end
1条答案
按热度按时间332nm8kg1#
我不知道任何关于伯努利的梁,但我知道很多关于数值分析。这可能不是你的问题,但它非常好的阅读我写的东西,无论你有多久一直在编程。
当除以整数时,计算机的进位位可能会导致错误的答案并产生舍入错误。我看到你到处使用/17,其他地方使用24,12,3和4。你必须始终使这些值浮动,以便它们在代码中显示为17.0,24.0,12.0,3.0和4.0。如果公式中的所有内容都是整数,结果可以被转换为整数。
此外,如果数字很小,并且您试图在数值计算中创建公差和收敛数字,舍入误差可能会变得非常显著,导致错误的答案。因此,理想情况下,绝对使用浮点数或双精度数,包括乘法。
特别要注意的是,在你的公式中,E是巨大的,因此1/E是微小的。这对我来说,舍入误差是你的近似和误差的一个很可能的原因。