2010-12-07 58 views
2

我想有這樣的事情,如何使初始條件成爲ndsolve中的一個變量?

w[w1_] := 
NDSolve[{y''[x] + y[x] == 2, y[0] == w1, y'[0] == 0}, y, {x, 0, 30}] 

這似乎是它的作品更好,但我想我缺少SMTN

w := NDSolve[{y''[x] + y[x] == 2, y[0] == w1, y'[0] == 0}, 
    y, {x, 0, 30}] 
w2 = Table[y[x] /. w, {w1, 0.0, 1.0, 0.5}] 

,因爲當我嘗試做一個表,它不「T工作:

Table[Evaluate[y[x] /. w2], {x, 10, 30, 10}] 

我得到一個錯誤:

ReplaceAll::reps: {<<1>>[x]} is neither a list of replacement rules nor a valid dispatch table, and so cannot be used for replacing. >> 

ps:有沒有更好的地方來提問這樣的問題? mathematica沒有支持的論壇,只有mathGroup電子郵件列表。這將是很好,如果stackoverflow將具有更多具體的數學標籤,如簡化,ndsolve,劇情操作

+0

我覺得有沒有足夠的MMA的用戶參與SO開分標籤,遺憾的是 – 2010-12-07 21:44:10

回答

4

有很多方法可以做到這一點。其一是:

w[w1_] := NDSolve[{y''[x] + y[x] == 2, 
        y'[0] == 0},  y[0] == w1, 
         y[x], {x, 0, 30}]; 

Table[Table[{w1,x,y[x] /. w[w1]}, {w1, 0., 1.0, 0.5}]/. x -> u, {u, 10, 30, 10}] 

輸出:

{{{0., 10, {3.67814}}, {0.5, 10, {3.25861}}, {1.,10, {2.83907}}}, 
{{0., 20, {1.18384}}, {0.5, 20, {1.38788}}, {1.,20, {1.59192}}}, 
{{0., 30, {1.6915}}, {0.5, 30, {1.76862}}, {1.,30, {1.84575}}}} 
1

我看你已經選擇了一個答案,但我想折騰這個解決方案爲線性方程組起來的家庭。具體來說,這是模擬Lotka-Volterra上的一個小變化。

(*Put everything in a module to scope x and y correctly.*) 
Module[{x, y}, 

(*Build a function to wrap NDSolve, and pass it 
       the initial conditions and range.*) 
soln[iCond_, tRange_, scenario_] := 
    NDSolve[{ 
    x'[t] == -scenario[[1]] x[t] + scenario[[2]] x[t]*y[t], 
    y'[t] == (scenario[[3]] - scenario[[4]]*y[t]) - 
     scenario[[5]] x[t]*y[t], 
    x[0] == iCond[[1]], 
    y[0] == iCond[[2]] 
    }, 
    {x[t], y[t]}, 
    {t, tRange[[1]], tRange[[2]]} 
    ]; 

(*Build a plot generator*) 
GeneratePlot[{iCond_, tRange_, scen_, 
    window_}] := 
    (*Find a way to catch errors and perturb iCond*)  
    ParametricPlot[ 
    Evaluate[{x[t], y[t]} /. soln[iCond, tRange, scen]], 
    {t, tRange[[1]], tRange[[2]]}, 
    PlotRange -> window, 
    PlotStyle -> Thin, LabelStyle -> Medium 
    ]; 

(*Call the plot generator with different starting conditions*) 
graph[scenario_, tRange_, window_, points_] := 
    {plots = {}; 
    istep = (window[[1, 2]] - window[[1, 1]])/(points[[1]]+1); 
    jstep = (window[[2, 2]] - window[[2, 1]])/(points[[2]]+1); 
    Do[Do[ 
    AppendTo[plots, {{i, j}, tRange, scenario, window}] 
    , {j, window[[2, 1]] + jstep, window[[2, 2]] - jstep, jstep} 
    ], {i, window[[1, 1]] + istep, window[[1, 2]] - istep, istep}]; 
    Map[GeneratePlot, plots] 
    } 
] 
] 

然後我們可以使用動畫(或桌子上,但動畫是真棒)

tRange = {0, 4}; 
window = {{0, 8}, {0, 6}}; 
points = {5, 5} 
Animate[Show[graph[{3, 1, 8, 2, 0.5}, 
     {0, t}, window, points]], {t, 0.01, 5}, 
     AnimationRunning -> False] 
相關問題