2012-08-16 66 views
0

我試圖在下面的代碼中繪製函數hh的圖形(請跳過代碼的最後一行)。我已經設置了PlotPoints->2MaxRecursion->0,但代碼仍在運行,運行了大約8個小時。函數hh非常複雜,涉及大量迭代。有什麼辦法讓代碼運行得更快嗎?當功能非常複雜時,如何更快繪製曲面?

s0[a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, x_] := 
{a, b, u, c, d, v, p, q, z, s, t, w, 1, 0, (a + b) x + u} 

s1[a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, n_, x_] := 
Which[0 <= x <= c, {a, b, u, c, d, v, p, q, z, s, t, w, 2, n - 1, x/c*q + p}, 
    c <= x <= c + d, {a, b, u, c, d, v, p, q, z, s, t, w, 2, n, (x - c)/d*p}, 
    c + d <= x <= 1, {a, b, u, c, d, v, p, q, z, s, t, w, 1, n + 1, (x - (c + d))/v*u}] 

s2[a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, n_, x_] := 
Which[0 <= x <= w, {a, b, u, c, d, v, p, q, z, s, t, w, 2, n - 1, x/w*z + p + q}, 
    w <= x <= 1 - s, {a, b, u, c, d, v, p, q, z, s, t, w, 3, n - 1, (x - w)/t*v + c + d}, 
    1 - s <= x <= 1, {a, b, u, c, d, v, p, q, z, s, t, w, 3, n, (x - (1 - s))/s*d + c}] 

s3[a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, n_, x_] := 
Which[0 <= x <= u, {a, b, u, c, d, v, p, q, z, s, t, w, 4, n - 1, x/u*t + w}, 
    u <= x <= 1 - a, {a, b, u, c, d, v, p, q, z, s, t, w, 4, n, (x - u)/b*w}, 
    1 - a <= x <= 1, {a, b, u, c, d, v, p, q, z, s, t, w, 3,n + 1, (x - (1 - a))/a*c}] 

s4[a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, n_, x_] := 
Which[0 <= x <= p, {a, b, u, c, d, v, p, q, z, s, t, w, 4, n - 1, x/p*s + 1 - s}, 
    p <= x <= p + q, {a, b, u, c, d, v, p, q, z, s, t, w, 5, n - 1, (x - p)/q*a/(a+ b)+ b/(a + b)}, 
    p + q <= x <= 1, {a, b, u, c, d, v, p, q, z, s, t, w, 5, n, (x - (p + q))/z*b/(a + b)}] 

f[{a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, k_, n_, x_}] := 
Which[k == 0, s0[a, b, u, c, d, v, p, q, z, s, t, w, x], 
    k == 1, s1[a, b, u, c, d, v, p, q, z, s, t, w, n, x], 
    k == 2, s2[a, b, u, c, d, v, p, q, z, s, t, w, n, x], 
    k == 3, s3[a, b, u, c, d, v, p, q, z, s, t, w, n, x], 
    k == 4, s4[a, b, u, c, d, v, p, q, z, s, t, w, n, x]] 

g[a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, x_] := 
NestWhile[f, {a, b, u, c, d, v, p, q, z, s, t, w, 0, 0, x}, Function[e, Extract[e, {13}] != 5]] 

h[a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, x_] := 
Extract[g[a, b, u, c, d, v, p, q, z, s, t, w, x], {15}] + 
    Extract[g[a, b, u, c, d, v, p, q, z, s, t, w, x], {14}] 

ff[{a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, x_}] := {a, b, u, c, d, v, p, q, z, s, t, w, h[a, b, u, c, d, v, p, q, z, s, t, w, x - Floor[x]] + Floor[x]} 

gg[a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_] := 
N[Extract[Nest[ff, N[{a, b, u, c, d, v, p, q, z, s, t, w, 0}], 10^3], {13}]/10^3] 

hh[x_, y_] := 
gg[x, y, 1 - x - y, x, y, 1 - x - y, x, y, 1 - x - y, x, y, 1 - x - y] 

Plot3D[hh[x, y], {x, 0, 1}, {y, 0, 1}, RegionFunction -> Function[{x, y, z}, x + y <= 1], PlotPoints -> 2, MaxRecursion -> 0] 
+0

某處有問題,你的函數:評估,例如,'HH [1,1]',打印信息和運行無限:'提取:: partw:「 f [s0 [1.,1., - 1.,1., - 1.,1.,1.,1.,1., - 1.,0。]的第13部分]不存在。「' – 2012-08-16 10:27:35

+0

函數'hh'不允許輸入'(1,1)'。參數'x,y'必須滿足'x + y <1'。感謝評論。 – 2012-08-16 19:34:04

回答

2

我覺得你的函數有一些問題。 但是,這裏有一些想法;我修改了幾個功能:

s0[a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, n_, x_] = 
{a, b, u, c, d, v, p, q, z, s, t, w, 1, 0, (a + b) x + u}  

(* so it matches the others *) 

f[{a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, k_, n_, x_}] := 
{s0[a, b, u, c, d, v, p, q, z, s, t, w, n, x],       
    s1[a, b, u, c, d, v, p, q, z, s, t, w, n, x], 
    s2[a, b, u, c, d, v, p, q, z, s, t, w, n, x], 
    s3[a, b, u, c, d, v, p, q, z, s, t, w, n, x], 
    s4[a, b, u, c, d, v, p, q, z, s, t, w, n, x]}[[k + 1]] 

g[a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, x_] := 
    NestWhile[f, {a, b, u, c, d, v, p, q, z, s, t, w, 0, 0, x}, #[[13]] != 5 &] 


h[a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_, x_] := 
    Total[g[a, b, u, c, d, v, p, q, z, s, t, w, x][[{14, 15}]]] 

(* in order to avoid calculating the same quantity twice *) 

gg[a_, b_, u_, c_, d_, v_, p_, q_, z_, s_, t_, w_] := 
    Nest[ff, {a, b, u, c, d, v, p, q, z, s, t, w, 0}, 10^3][[13]]/10^3 

{elapsed, data} = Outer[If[#1 + #2 <= 1, {#1, #2, hh[#1, #2]}, {#1, #2, "bad"}] &, 
    Range[0.1, 1, 0.25], Range[0.1, 1, 0.25]] // AbsoluteTiming 

(* {6.820061, {{{0.1, 0.1, -1.99961}, {0.1, 0.35, -1.99961}, {0.1, 0.6, -2.00009}, 
    {0.1, 0.85, -2.00001}}, {{0.35, 0.1, -1.99993}, {0.35, 0.35, -2.00004}, 
    {0.35, 0.6, -2.00017}, {0.35, 0.85, "bad"}}, {{0.6, 0.1, -1.99996}, 
    {0.6,0.35, -2.00024}, {0.6, 0.6, "bad"}, {0.6, 0.85, "bad"}}, 
    {{0.85, 0.1, -1.99967}, {0.85, 0.35, "bad"}, {0.85, 0.6, "bad"}, {0.85, 0.85, "bad"}}}} *) 

ListPlot3D[Select[Flatten[data, 1], NumericQ[#[[3]]] &]] 

example

+0

@gatessucks謝謝你的回答。修改後的代碼現在只在我的機器上運行8秒鐘。我確實有一個問題。爲什麼需要使用「選擇」? 'data'中的點似乎都是數字。 – 2012-08-16 19:59:03

+1

@MichaelC不,我用了一個快速的方法來得到一些要點。我包含了不符合約束「x + y <= 1」的點(因此z值爲「bad」),因此我需要在繪圖之前刪除這些點。 – 2012-08-17 08:17:28