フィルターのクリア

Use second-order and fourth-order Runge-Kutta methods

2 ビュー (過去 30 日間)
Cassidy Holene
Cassidy Holene 2021 年 3 月 6 日
回答済み: darova 2021 年 3 月 9 日
The one-dimensional linear convection equation for a scalar T is,
𝜕𝑇/𝜕𝑡 + 𝑎*𝜕𝑇/𝜕𝑥 = 0 0 ≤ 𝑥 ≤ 1
With the initial condition and boundary conditions of, (𝑥, 0) = 𝑒 ^(−200(𝑥−0.25)^ 2) , 𝑇(0,𝑡) = 0, 𝑇(1,𝑡): outflow condition Take 𝑎 = 1, ∆𝑥 = 0.01
Plot the results for t = 0.25, .5, 0.75 for both of them.
  2 件のコメント
darova
darova 2021 年 3 月 6 日
DO you have any attempts? Do you have difference scheme?
Cassidy Holene
Cassidy Holene 2021 年 3 月 6 日
a = 1; %initial condition
f = %given function
yprime = f(x,y) %direvative of f
xRange = [0,1]; %where the solution is sought on 0<=x<=1
yInitial = exp(-200*(xRange-0.25).^2); %column vector of initial values for y at x=0
numSteps = 100; %number of equally-sized steps to take from x1 to x2
t = 0.25; %first value of t changes to 0.5 and 0.75
x=zeros(numSteps+1,1); %initalzie x
x(1) = xRange(1); %set the first x value
h = ( xRange(2) - xRange(1) ) / numSteps; %time step
y(1,:) = (yInitial)';
%first try
for k = 1 : numSteps
x(k) = ??? ; %first x at k
y(k) = ??? ; %first y at k
x(k+1) = x(k) + h; %x at next k value
y(k+1,:) = y(k,:) + h * (feval( ??? ))'; %y at next k value
end
%second try
for i=1:(length(x)-1)
k1 = h*f(x(i),y(i)); %1st order
k2 = h*f(x(i)+0.5*k1,y(i)+0.5*h); %2nd order
k3 = h*f((x(i)+0.5*k2),(y(i)+0.5*h)); %3rd order
k4 = h*f((x(i)+k3),(y(i)+h)); %4th order
y(i+1) = y(i) + (1/6)*(k1+2*k2+2*k3+k4)*h; % main equation
end
figure, plot(x, y) % To see the solution results

サインインしてコメントする。

回答 (1 件)

darova
darova 2021 年 3 月 9 日
I'd try simple euler scheme first. The scheme i used
𝜕𝑇/𝜕𝑡 + 𝑎*𝜕𝑇/𝜕𝑥 = 0
clc,clear
a = 1;
x = 0:0.01:1;
t = 0:0.01:0.75;
U = zeros(length(t),length(x));
U(1,:) = exp(-200*(x-0.25).^2);
U(:,1) = 0;
dx = x(2)-x(1);
dt = t(2)-t(1);
for i = 1:size(U,1)-1
for j = 1:size(U,2)-1
U(i+1,j) = U(i,j) + a*dt/dx*(U(i,j+1)-U(i,j));
end
end
surf(U)

カテゴリ

Help Center および File ExchangeVector Fields についてさらに検索

タグ

Community Treasure Hunt

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

Start Hunting!

Translated by