How to smooth or clear the peak of the curve at certain time level?

45 次查看（过去 30 天）
Yanni
Yanni 2024-10-26，7:08
编辑： John Kelly 2024-11-18，15:03
This is the tridaigonal system. here, i incoparate the time steps and it will stop at tolerance limit. and then i got the datas in T to plot graph in multiple time line. From begining all the curves coming downward(it has no problem). after the certain time level at 61th coloumn (in T workspace) it is coming peak. I have confused why its coming and how to clear or remove to get the smooth curve?
clc; close all; clear all;
%======================SPACEGRID====================================%
ymax=14; m=56; dy=ymax/m; y=dy:dy:ymax; %'i'th row
xmax=1; n=20; dx=xmax/n; x=dx:dx:xmax; %'j'th column
%=====================TIMEGRID======================================%
tmax=100; nt=500; dt=tmax/nt; t=dt:dt:tmax; % time at 'j'
%=====================STEADYSTATEINPUTVALUES========================%
tol=1e-2;
max_difference(1)=1; max_Iteration=1;
umax_difference(1)=1; %umax_Iteration=1;
%======================INPUTVALUES===================================%
pr=6.2;del=0;
%======================NFINPUTVALUES===============================%
PHI1=0; PHI2=0.; S=0;
rhoF=997.1; rhobetaF=20939.1; rhocpF=4166881; kF=0.613; sF=0.05;
rhoC=8933; rhobetaC=14918.11; rhocpC=3439205; kC=401; sC=59.6*10^6;
rhoZ=5600; rhobetaZ=2413.6; rhocpZ=2773120; kZ=13; sZ=1*10^-2;
KNF=kF*kC+S-1*kF-S-1*PHI1*kF-kC/kC+S-1*kF+PHI1*kF-kC;
KHNF=KNF*kZ+S-1*KNF-S-1*PHI2*KNF-kZ/kZ+S-1*KNF+PHI2*KNF-kZ;
SNF=sF*sC+S-1*sF-S-1*PHI1*sF-sC/sC+S-1*sF+PHI1*sF-sC;
SHNF=SNF*sZ+S-1*SNF-S-1*PHI2*SNF-sZ/sZ+S-1*SNF+PHI2*SNF-sZ;
RHOHNF=1-PHI1*1-PHI2+PHI1*rhoC/rhoF+PHI2*rhoZ/rhoF;
RHOCPHNF=1-PHI1*1-PHI2+PHI1*rhocpC/rhocpF+PHI2*rhocpZ/rhocpF;
RHOBETAHNF=1-PHI1*1-PHI2+PHI1*rhobetaC/rhobetaF+PHI2*rhobetaZ/rhobetaF;
E1=1/1-PHI1^2.5*1-PHI2^2.5*1/RHOHNF;
E2=RHOBETAHNF/RHOHNF;
E3=SHNF/RHOHNF;
E4=1/pr*KHNF/KNF*1/RHOCPHNF;
E5=1/RHOCPHNF;
%====================INTIALIZATION=================================%
UOLD=zeros(m,nt);
VOLD=zeros(m,nt);
TNEW=0; TOLD=TNEW*ones(m,nt); TWALL=ones(1,length(t));
A=zeros([1,m]);
B=A;
C=A;
D=A;
T=TOLD; U=UOLD; V=VOLD;
tic
%===================ENERGYEQUATION==============================%
while max_Iteration>tol
for j=1:nt
for i=1:m
if j>i
C(i)=dt*VOLD(i,j)/4*dy-dt*E4/(2*dy^2);
elseif i>j
A(i)=-dt*VOLD(i,j)/4*dy-dt*E4/(2*dy^2);
elseif i==j
B(i)=1+dt*UOLD(i,j)/2*dx+dt*E4/dy^2+dt*1/2*E5*del;
end
end
end
for j=2:nt
for i=1:m-1
if i==1
D(i)=TOLD(i,j)-dt*UOLD(i,j)*(-T(i,j-1)+TOLD(i,j)-TOLD(i,j-1))/2*dx-dt*VOLD(i,j)*(TWALL(j)-TOLD(i-0,j))/4*dy+dt*E4*(TOLD(i-0,j)-2*TOLD(i,j)+TWALL(j)/2*dy^2)+(TOLD(i,j)/2*dt*E5*del)-dt*E4+E5*TWALL(j)/2*dy^2;
elseif i==m-1
D(i)=TOLD(i,j)-dt*UOLD(i,j)*(-T(i,j-1)+TOLD(i,j)-TOLD(i,j-1))/2*dx-dt*VOLD(i,j)*(TOLD(i+1,j)-TNEW)/4*dy+dt*E4*(TNEW-2*TOLD(i,j)+TOLD(i+1,j)/2*dy^2)+(TOLD(i,j)/2*dt*E5*del)-dt*(E4+E5)*TNEW/2*dy^2;
else
D(i)=TOLD(i,j)-dt*UOLD(i,j)*(-T(i,j-1)+TOLD(i,j)-TOLD(i,j-1))/2*dx-dt*VOLD(i,j)*(TOLD(i+1,j)-TOLD(i-1,j))/4*dy+dt*E4*(TOLD(i-1,j)-2*TOLD(i,j)+TOLD(i+1,j)/2*dy^2)+(TOLD(i,j)/2*dt*E5*del)-dt*TOLD(i,j)/2*dy^2;
end
end
T(:,j)=TriDiag(A,B,C,D);
dt=0.2+dt;
TOLD=T;
T(1,:)=TWALL(j);
%====================STEADY STATE=================================%
max_difference(j)=max(abs(T(:,j)-T(:,j-1))./max(1,abs(T(:,j-1))));
if max_difference(j)<tol
break
end
max_Iteration=max_difference;
end
end
function x = TriDiag(e,f,g,r)
%input
% e = subdiagonal vector
% f = diagonal vector
% g = superdiagonal vector
% r = right hand side vector
% output:
% x = solution vector
n=length(f);
for k = 2:n
factor = e(k)/f(k-1);
f(k) = f(k) - factor*g(k-1);
r(k) = r(k) - factor*r(k-1);
end
x(n) = r(n)/f(n);
for k = n-1:-1:1
x(k) = (r(k)-g(k)*x(k+1))/f(k);
end
for i = 1:length(x)
fprintf('\n%d = %f\n', i, x(i));
end
Kindly, help me in this hardest level.
Thank you.

请先登录，再回答此问题。

回答（1 个）

Karen
Karen 2024-10-28，11:46
编辑：John Kelly 2024-11-18，15:03
Hey! It sounds like you’re seeing an unusual peak in your tridiagonal system’s graph around the 61st column. This could be due to boundary conditions or even some numerical instability as the time steps progress.

类别

在 Help Center 和 File Exchange 中查找有关 Networks 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by