WoGu Wiki

y=[0.67 0.67 0.67 0.67 0.67 0.02 0.02 0.02 0.02];
 
temp = 727;
R = 8.3144;
Q = 140000;
temp_k = temp + 273;
d0= 0.000041;
D = d0 * exp( -Q /( R * temp_k)) * 1 * 10^10;
 
dt = 1;
dx = 1;
 
tmp = (D*dt)/dx^2;
 
n = 100;
[s_h s_w] = size(y);
 
y2 = zeros(1, s_w);
 
for j = 1:n
    for i = 1:s_w
 
        im = i-1;
        ip = i+1;
        if im < 1
            im = 1;
        end
        if ip > s_w
            ip = s_w;
        end
        y2(i) = (1 - 2*tmp)*y(i) + tmp*(y(im) + y(ip));
 
    end
 
    y = y2;
end
 
plot(y2, 'DisplayName', 'y2', 'yDataSource', 'y2');

temp=727;
Q = 140000;
R = 8.3144;
temp_k = temp + 273;
d0 = 0.000041;
D=d0*exp(-Q /(R * temp_k))*1E10;
dt = 0.01;
dx = 0.1;
pN = 0.1;
%Kryterium stabilności Neumanna:
if (D*dt)/dx^2<=0.5
%(D*dt)/dx^2<=0.5
 
%Obszar rozwiązania: 100 mm
%pierwszy i ostatni element powielamy
    K(6:100) = 0.02;
    K(1:5) = 0.67;
    ksi = 6;
    for k = 0:dt:10
 
        K = [K(1), K, K(length(K))];
        K1 = K;
 
        for i = 2:1:ksi+1%length(K)-1
           K1(i) = D*(dt/dx^2)*(K(i+1)+K(i-1)) + K(i)*(1 - 2*D*(dt/dx^2));
        end
 
        K = K1;
        K = K(2:1:length(K)-1);
 
        plot(K)
        hold on
 
        y=-0.004074074074074075*temp_k+3.7155555555555555;
        if K(ksi) >= y && ksi < length(K)-1
            ksi = ksi + 1  
        end
        temp_k = temp_k + pN * dt;
        D=d0*exp(-Q /(R * temp_k))*1E10;
 
    end
   % disp(K);
end