解決一個延遲微分方程（DDE）系統約束非負的解決方案

Question

在MATLAB中，ode45有一個參數NonNegative約束的解決方案是非負。 They even wrote a paper about how this method works以及它是如何愚蠢的，只要它變成負數就將y_i設置爲0，因爲這通常不起作用。

現在，MATLAB也有dde23用於求解延遲微分方程，但這個積分器沒有等效的參數NonNegative。

不幸的是，我的任務是爲現有ODE系統添加延遲，該系統使用ode45解決，啓用NonNegative。

任何想法我應該如何進行？

編輯：

我不知道這是有幫助的，但是......

我的系統的DDE部分基本上看起來像：

dx = 1/(1+y*z) - x; 
dy = (y*z)^2/(1+(y*z)^2) - y; 
dz = X - z; 

其中X（首都在第三個等式中的字母變量）是x的延遲版本。然後，我通過在x和z的公式中加入一對項來將這個DDE系統連接到現有的（和更大的）ODE系統，然後將所有組合系統集成到一起。

Answer 1

你有一個棘手的問題，我不確定是否有一個一步的解決方案。我會更樂意爲願意提供替代答案的任何人提供榮譽。

根據延遲的長度，一個選項是多次運行該方程，每次迭代都將舊的x值傳遞到最新的更新。

例如，假設您的延遲時間爲一小時。在第一個小時內，運行標記爲NonNegative的ode45。將值與Time參數一起存儲到新矩陣中，然後再次運行該算法。這個時候請確保您添加兩個輸入參數：您的舊解決方案矩陣和舊時代的矩陣

dx = 1/(1+y*z) - x; 
dy = (y*z)^2/(1+(y*z)^2) - y; 
tindex = find(told>t,1) -1 % find the upper index which best approximates t 
X = xold(tindex) + (xold(tindex+1)-xold(tindex))*(t-told(tindex))/(told(tindex+1)-told(tindex)) % or interpolation method of your choosing 
dz = X - z; 

現在清洗，沖洗，然後重複。請注意，X現在是一個準時間相關項，如示例3中的ode45所示。

Answer 2

我最近有這個問題，我的一些代碼。在「簡單」的解決方案是做到以下幾點，首先，一旦溶液達到0，你必須保持在0， 從而改變

dx = 1/(1+y*z) - x; 

來（注意一下這裏x == 0情況評價）

if x > 0 
    dx = 1/(1+y*z) - x; 
else % if x <= 0 
    dx = 0; 
end 

或也許向（取決於爲什麼它可能永遠不會爲0）

dxTmp = 1/(1+y*z) - x; 
if x > 0 
    dx = dxTmp; 
elseif dxTmp > 0 
    % x becomes positive again 
    dx = dxTmp; 
else 
    dx = 0; 
end 

注意然而，這會在第一衍生物的跳躍不連續，WHI當DDE求解器到達接近這個點的位置時，它不能很好地解決它，除非它知道這個不連續點的確切位置（通常你使用一個額外的選項來告訴Matlab它在哪裏，但是如果你按照下面的步驟，那麼將不需要）。

要確定此中斷的地方，你需要使用的events（DDE選項向下滾動到「活動地點屬性」，你也可以看看這些examples，其中的一個例子其實有類似的系統，其中負交易值在ODE中是不允許的 - ODE和DDE的事件幾乎是相同的）。基本上，一個事件是帶有向量輸出的Matlab函數，向量的每個條目都是對變量的一些或其他評估;在每一步Matlab都檢查它們中的一個是否爲eqauls 0，當其中一個等於0時，DDE停止並返回到該點的解決方案，從中必須重新啓動具有該部分解決方案的DDE作爲歷史記錄，即不是運行

sol = dde23(ddefun, lags, history, [t0 tEnd], options); 

運行（注意sol和t0改變）

sol = dde23(ddefun, lags, sol, [tCurrent tEnd], options); 

在這種情況下，矢量的條目之一將是x（只要你想的DDE停止時x等於0）。此外，代碼行elseif dxTmp <= 0會產生另一個不連續性，因此您需要一個事件，當dxTmp變爲0時，即1/(1+y*z) - x將成爲矢量輸出的另一個組件。

現在，當您重新啓動ODE時，Matlab會自動假定在該點存在不連續性，因此您不必擔心告訴Matlab該處有一個。

接下來的問題是，Matlab的從未完全達到其權，x，y，z和X值會略有下降。因此，如果要創建一個問題，你會希望在計算衍生品之前

if x < 0 
    x = 0; 
end 

糾正x值（和類似的其它值）。但這隻會在本地更改x的值。因此，您可能需要將最終解決方案中x的所有負值更改爲0。我建議你在輸入sol = dde23(ddefun, lags, sol, [tCurrent tEnd], options);之前不要試圖更改sol，因爲我在這方面做了很多嘗試，但無法使其工作。

Answer 3

還有一個多更簡單的回答，使用createOptimProblem設置優化。您必須使用此方法爲每個參數包含邊界，但強制參數保持正值變得微不足道。

詳情點擊這裏：http://www.mathworks.com/help/gads/createoptimproblem.html