Matlab to solve finite difference for Ode
2 ビュー (過去 30 日間)
古いコメントを表示
% Finite difference method
% Approximate the solution of y"=(-2/x)y'+(2/x^2)y+ sin(lnx)/x^2
% for 1<=x<=2 with y(1)=1 and y(2)=2.
p = @(x) (-2/x);
q = @(x) (2/x^2);
r = @(x) (sin(log(x))/x^2);
aa = 1; bb = 2; alpha = 1; beta = 2; n=29; % h = (bb-aa)/(n+1); h=0.1
fprintf(' x w \n');
h = (bb-aa)/(n+1);
a = zeros(1,n+1);
b = zeros(1,n+1);
c = zeros(1,n+1);
d = zeros(1,n+1);
l = zeros(1,n+1);
u = zeros(1,n+1);
z = zeros(1,n+1);
w = zeros(1,n+1);
x = aa+h;
a(1) = 2+h^2*q(x);
b(1) = -1+0.5*h*p(x);
d(1) = -h^2*r(x)+(1+0.5*h*p(x))*alpha;
m = n-1;
for i = 2 : m
x = aa+i*h;
a(i) = 2+h^2*q(x);
b(i) = -1+0.5*h*p(x);
c(i) = -1-0.5*h*p(x);
d(i) = -h^2*r(x);
end
x = bb-h;
a(n) = 2+h^2*q(x);
c(n) = -1-0.5*h*p(x);
d(n) = -h^2*r(x)+(1-0.5*h*p(x))*beta;
l(1) = a(1);
u(1) = b(1)/a(1);
z(1) = d(1)/l(1);
for i = 2 : m
l(i) = a(i)-c(i)*u(i-1);
u(i) = b(i)/l(i);
z(i) = (d(i)-c(i)*z(i-1))/l(i);
end
l(n) = a(n)-c(n)*u(n-1);
z(n) = (d(n)-c(n)*z(n-1))/l(n);
w(n) = z(n);
for j = 1 : m
i = n-j;
w(i) = z(i)-u(i)*w(i+1);
end
i = 0;
fprintf('%5.4f %11.8f\n', aa, alpha);
for i = 1 : n
x = aa+i*h;
fprintf('%5.4f %11.8f\n', x, w(i));
end
i = n+1;
fprintf('%5.4f %11.8f\n', bb, beta);am trying to run d code but it keeps writing Error: Invalid text character. The text ' ' contains an unsupported non-ASCII whitespace character. Please where am I doing wrong?
0 件のコメント
回答 (2 件)
Viranch Patel
2021 年 7 月 9 日
I tried running the code on my machine and it's working fine. You can paste this code to your machine and see if its working. I am attaching the output too. You might be using some wrong character for any particular character which is not supported.
% Finite difference method
% Approximate the solution of y"=(-2/x)y'+(2/x^2)y+ sin(lnx)/x^2
% for 1<=x<=2 with y(1)=1 and y(2)=2.
p = @(x) (-2/x);
q = @(x) (2/x^2);
r = @(x) (sin(log(x))/x^2);
aa = 1; bb = 2; alpha = 1; beta = 2; n=29; % h = (bb-aa)/(n+1); h=0.1
fprintf(' x w \n');
h = (bb-aa)/(n+1);
a = zeros(1,n+1);
b = zeros(1,n+1);
c = zeros(1,n+1);
d = zeros(1,n+1);
l = zeros(1,n+1);
u = zeros(1,n+1);
z = zeros(1,n+1);
w = zeros(1,n+1);
x = aa+h;
a(1) = 2+h^2*q(x);
b(1) = -1+0.5*h*p(x);
d(1) = -h^2*r(x)+(1+0.5*h*p(x))*alpha;
m = n-1;
for i = 2 : m
x = aa+i*h;
a(i) = 2+h^2*q(x);
b(i) = -1+0.5*h*p(x);
c(i) = -1-0.5*h*p(x);
d(i) = -h^2*r(x);
end
x = bb-h;
a(n) = 2+h^2*q(x);
c(n) = -1-0.5*h*p(x);
d(n) = -h^2*r(x)+(1-0.5*h*p(x))*beta;
l(1) = a(1);
u(1) = b(1)/a(1);
z(1) = d(1)/l(1);
for i = 2 : m
l(i) = a(i)-c(i)*u(i-1);
u(i) = b(i)/l(i);
z(i) = (d(i)-c(i)*z(i-1))/l(i);
end
l(n) = a(n)-c(n)*u(n-1);
z(n) = (d(n)-c(n)*z(n-1))/l(n);
w(n) = z(n);
for j = 1 : m
i = n-j;
w(i) = z(i)-u(i)*w(i+1);
end
i = 0;
fprintf('%5.4f %11.8f\n', aa, alpha);
for i = 1 : n
x = aa+i*h;
fprintf('%5.4f %11.8f\n', x, w(i));
end
i = n+1;
fprintf('%5.4f %11.8f\n', bb, beta);
Output:
x w
1.0000 1.00000000
1.0333 1.03067949
1.0667 1.06155254
1.1000 1.09262608
1.1333 1.12390423
1.1667 1.15538887
1.2000 1.18708018
1.2333 1.21897697
1.2667 1.25107701
1.3000 1.28337728
1.3333 1.31587416
1.3667 1.34856355
1.4000 1.38144105
1.4333 1.41450201
1.4667 1.44774163
1.5000 1.48115505
1.5333 1.51473732
1.5667 1.54848354
1.6000 1.58238883
1.6333 1.61644834
1.6667 1.65065735
1.7000 1.68501118
1.7333 1.71950528
1.7667 1.75413522
1.8000 1.78889666
1.8333 1.82378540
1.8667 1.85879736
1.9000 1.89392857
1.9333 1.92917520
1.9667 1.96453354
2.0000 2.00000000
0 件のコメント
Walter Roberson
2024 年 12 月 1 日
S = 'aa = 1; bb = 2; alpha = 1; beta = 2; n=29; % h = (bb-aa)/(n+1); h=0.1'
idx = find(S>128)
S(idx)
S(idx)+0
Your code contains a number of char(160), which is the Unicode Non-Breaking Whitespace character https://www.compart.com/en/unicode/U+00A0
The first group of NBSP is after the n=29;
S(40:50)