2013-08-24 4 views
1

용지의 결과를 복제하려고합니다 (Immune Response). 즉, 병원균, 항체 및 혈장 세포의 농도와 면역 반응이 장기에 미치는 영향을 살펴 봅니다. 장기의 상태에 따라, 형질 세포의 농도에 영향을 미치는 변수가 변할 것입니다.Matlab의 ODE 변수 업데이트

저는 이것을 Matlab에서 모델링했습니다. 내 문제는 혈장 농도에 영향을 미치는 변수에 대한 조건문이 무시된다는 것입니다. 두 번째 조건이 무엇인지에 관계없이 코드의 일부를 주석 처리하고 그래프를 보면 출력은 초기 조건에만 응답합니다.

ode45 (또는 ode23)에 의해 해석되도록 값을 업데이트하고이 문제에 접근하는 더 좋은 방법이 있습니까?

내 코드는 다음과 같습니다.

function dx = iir_2(t, x) 
dx = [0; 0; 0; 0]; 
a11 = 1; 
a31 = 1; 
a41 = 1; 
a12 = 1; 
a22 = 3; 
a32 = 1.5; 
a42 = 1; 
a23 = 1; 
a33 = 0.5; 

b1 = -1; 
b2 = 1; 
b3 = 1; 
b4 = -1; 

u1 = 0; 
u2 = 0; 
u3 = 0; 
u4 = 0; 

tau = 0.1; 

if x(4) >= 0.5 
    a21x4 = 0; 
else 
    a21x4 = cos(pi * x(4)); 
end 

if x(4) > 1 
    x(4) = 1; 
end 

dx(1) = (a11 - a12 * x(3)) * x(1) + b1 * u1; 
dx(2) = a21x4 * a22 * x(1) * (t - tau) * x(3) * (t - tau) - a23 * (x(2) - 2) + b2 * u2; 
dx(3) = a31 * x(2) - (a32 + a33 * x(1)) * x(3) + b3 * u3; 
dx(4) = a41 * x(1) - a42 * x(4) + b4 * u4; 

위의 함수가 호출이 두 번째 스크립트 ..

Case = ['Subclinical' 'Clinical' 'Chronic' 'Lethal']; 
x1_0 = [1.5 2 2.57 3]; 
x2_0 = 2; 
x3_0 = (1*x2_0)/1.5; 
x4_0 = 0; 

for i = 1:4 
    if i == 1 
     state = 'Subclinical'; 
    elseif i == 2 
     state = 'Clinical'; 
    elseif i == 3 
     state = 'Chronic'; 
    else 
     state = 'Lethal'; 
    end 

    options = odeset('RelTol', 1e-3, 'NonNegative', [1 2 3 4]); 
    [t,x] = ode45(@iir_2, [0 10], [x1_0(i) x2_0 x3_0 x4_0], options); % use ode23??? 
    figure 
    plot(t,x); 
    str = sprintf('Case Number = %d\nx(0) = %d\n%s', i, x1_0(i), state); 
    title(str); 
    axis([0,10,0,10]) 
    legend('Pathogen', 'Plasma Cell', 'Antibody', 'Organ'); 
end 

조건식은 기능에 있습니다.

if x(4) >= 0.5 
     a21x4 = 0; 
    else 
     a21x4 = cos(pi * x(4)); 
    end 

    if x(4) > 1 
     x(4) = 1; 
    end 

답변

2

이렇게하는 데는 두 가지 방법이 있습니다.

어느 쪽이든 x (4)를 1로 설정하고 모든 것이 최선의 결과를 내기를 바랍니다. Matlab은 x (4)의 이전 값이 메모리에 모두 저장되어 있기 때문에 신경 쓰지 않습니다. 또한, x (4)의 향후 값은 dx (4)와 이전에 보유 된 x (4) 값 (사용자가 설정할 수 없음)에 의해 결정됩니다.

당신은 당신의 문제에 대한 두 가지 해결책이 있습니다

A) 세트 DX을 (4) = 0 당신의 상태가 완벽한 X 발생하지 않습니다하지만

if x(4) > 1 
    dx(4) = 0; 
else 
    dx(4) = a41 * x(1) - a42 * x(4) + b4 * u4; 
end 

을 충족 될 때 (4) = 1, 약간의 오류가있을 것입니다.

B) herehere에 대해 알 수있는 이벤트 함수를 호출합니다.이 경우 새 함수 호출을 위해 인터럽트 할 수 있습니다.

+0

@Rasman. 나는 아직까지 동일한 결과물을 얻지 못했다. 그러나 당신의 반응은 나를 올바른 방향으로 이끌었다. – cer