Thomas algorithm tri diagonal matrix

16 次查看(过去 30 天)
Asha
Asha 2023-12-21
In this code i want a tridiagonal matrix with different diagonal elements. but after some point i got stuck with the elments.i have attached the form of the matrix in image. please help me with this code.
h = 0.01;
x0 = 0;x1 = 1;
b = 1;
t0 = 0; t1 = 0.1;
k = 0.01;
x = x0: h: x1;
sizex = x;
t = t0: k: t1;
M = length(x);
N = length(t);
function YYp = Tridiag( p0,h,k,M)
alf=0.5;
B1=@(x,t)2;
B1x=@(x,t)0;
B2=@(x,t)4;
B2x=@(x,t)0;
R1=-3*B1x/h+6*B1/h^2;
R2=-3*B2x/h+6*B2/h^2;
S1=-12*B1/h^2;
S2=-12*B2/h^2;
T1=6*B1/h^2+3/h*B1x;
T2=6*B2/h^2+3/h*B2x;
A = zeros(M,M);
%Mvalue = M;
R(1)=0;
P=zeros(1,N);
Q=zeros(1,N-1);
R=zeros(1,N);
Y=zeros(1,N-1);
L1 = 0;
for i = 1:N
X(i) = (i-1) * h;
if X(i) <= alf
L1 = i;
end
end
for j=1:M
T(j)=(j-1)*k;
end
for j=2:M
for i = 2:L1-1
P(i, i) = S1;
end
for i = 2:L1-2
Q(i, i+1) = R1;
end
for i = 2:L1-2
R(i+1, i) = T1;
end
for i = L1+2:N
P(i, i) = S2;
end
for i = L1+2:N-1
Q(i, i+1) = R2;
end
for i = L1+2:N-1
R(i+1, i) = T2;
end
i=L1
P(i, i) = pi*S1;
Q(i, i+1) =6* R1;
R(i+1, i) = pi*T1;
i=L1+1
P(i, i) = pi*S2;
Q(i, i+1) =6* R2;
R(i+1, i) = pi*T2;
end
A = diag([4*R1 P(2:L1-1) R2-T2, 0]) + diag([0 R(2:L1-2)], 1) + diag([Q(2:L1-2) 0], -1);
ft = p0;
alpha0 = 0;
alpha1 = 0;
yp0 = alpha0;
yp1 = alpha1;
dm = h^2*ft;
dm(1) = dm(1) - yp0*(P(1));
dm(M) = dm(M) - yp1*(P(M));
Amatrix = A;
dm = dm';
sizeA = size(A);
sizedm = size(dm);
cc = A\dm;
ccm1 = 0 - 4*cc(1) - cc(2);
ccmp1 = 0 - 4*cc(M) - cc(M-1);
ccc = [ccm1 cc' ccmp1];
YYp = ccc(1:M) + 4*ccc(2:M+1) + ccc(3:M+2);
end

回答(1 个)

John D'Errico
John D'Errico 2023-12-21
编辑:John D'Errico 2023-12-21
DON'T WRITE YOUR OWN CODE FOR A PROBLEM LIKE THIS.
Instead, just create the tridiagonal matrix, as a sparse matrix. Then use backslash to do the solve. All of that is literally trivial, requiring not much more than two calls, one to spdiags, and a second to backslash. Instead, you spent probably hours writing messy, heavily looped, unreadable, uncommented, and un-debuggable code that you got wrong. and now cannot fix.
Use MATLAB as it is designed to be used.
Create an nx3 array of the elements on the three diagonals. Call spdiags. Use backslash, or perhaps linsolve or decomposition as you wish.

类别

Help CenterFile Exchange 中查找有关 Linear Algebra 的更多信息

产品


版本

R2020a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by