---
# System prepended metadata

title: 質問部屋

---

# 質問部屋

```julia
using Images 
using FFTW
using NLsolve


#テストするための画像を作る．
u = zeros(128,128) #真っ黒な画像 
# make a view
img_a = @view u[1:128,24:48] #この行列を変化させると大元の行列も変えることができる．
img_a =  fill!(img_a, 1.0) 

img_b = @view u[1:128,72:96] 
img_b =  fill!(img_b, 1.0) 

img_c = @view u[1:128,104:112] 
img_c =  fill!(img_c, 1.0) 

img_d = @view u[40:70,1:128] 
img_d =  fill!(img_d, 0.5) 

 f = u #動かしたくない領域を参照するための行列
 u #動かすための行列．こいつを更新していく


dt = 0.05
e = 0.1
C1 = 300
C2 = 150000
lam = 50000




function nabla1(u)
  fftu = fft(u)
    fftu[1,:] = fftu[1,: ] * 0
    for j in 1:63
        fftu[j+1,: ] =      (2 *pi  *j *1im /128 ) * fftu[j+1 ,: ]
        fftu[129-j ,: ] =  (2 *pi  *(-j) *1im /128 )*fftu[129-j ,: ]
    end #j
real(ifft(fftu))
end

function nabla2(u)
  fftu = fft(u)
    fftu[: ,1 ] = fftu[: ,1 ] * 0
    for k in 1:63
        fftu[: ,k+1 ] =      (2 *pi  *k *1im /128 ) * fftu[: ,k+1 ]
        fftu[: ,129-k ] =  (2 *pi  *(-k) *1im /128 )*fftu[: ,129-k ]
 end #k
real(ifft(fftu))
end

function lap(u) #ラプラシアンの微分演算子．
 nabla1(nabla1(u))  + nabla2(nabla2(u))
end
function lambda(u) #微分方程式に出てくるラムダの項を計算する関数．ok
    lamb =zeros(128,128)
        for i in 1:39
            for j in 1:128
                lamb[i,j] = lam * (f[i,j]-u[i,j])
            end #j
        end #i
        for i in 71:128
            for j in 1:128
                lamb[i,j] = lam * (f[i,j]-u[i,j])
            end #j
        end #i
    lamb
end #function


#NLsolveを用いる
function nls(func, params...; ini = [0.0])
    if typeof(ini) <: Number
        return nlsolve((vout,vin)->vout[1]=func(vin[1],params...), [ini]).zero[1]
    else
        return nlsolve((vout,vin)->vout .= func(vin,params...), ini).zero
    end
end
#nls(f,v,ini = u)

#非線形ソルバーで解くために方程式のスキームを書いた．
 g(v)=(v -u)/dt   + e * lap(lap(v)) - C1 * lap(v) + C2 * v -( lap((1/e) *4*(u^3) -6 * (u^2) + 2*u) + lambda(u) -C1 *(lap(u)) + C2*u)
 
 nls(g,v, ini = u)
 
```

今書いてるやつ貼り付けますね

おｋ

これの最後で数値計算するんですがこけましたね

実行してみる。
どうでしたか？

うん、エラーがでな
実行結果どうなりましたか


```
g (generic function with 1 method)
```

だね。

定義しただけで終わってませんかそれ

たぶん

```
g(v)
```

を実行して終わった。
最後の行のnls()実行すれば計算がone step進められます．

ごめん
```
nls(g,v, ini = u)
```

忘れた。
あ，これだとvの定義忘れてます

dane

vの定義は
v = uです
おｋ
vはu_nに対してu_{n+1}の意味です

現状こんな感じ(めんどうだからgistに上げた)

ごめん、誤っていま消した。
もう一度発行する。

発行しなおした。

https://gist.github.com/muripoLife/500151545e6cb696a25314292d48a478

using NLsolve

が抜けていますね

https://gist.github.com/muripoLife/500151545e6cb696a25314292d48a478

いま、上げなおした。
エラーこれで同じ？

これです．同じです

じゃ、考えるか
エラーコードの
[12] (::var"#kw##nls")(::NamedTuple{(:ini,),Tuple{Array{Float64,2}}}, ::typeof(nls), ::Function, ::Array{Float64,2}) at ./none:0
にこけた原因が書かれてそうだと思いました（細かいところは読めてない．）

```
#NLsolveを用いる
function nls(func, params...; ini = [0.0])
    if typeof(ini) <: Number
        return nlsolve((vout, vin)->vout[1] = func(vin[1], params...), [ini]).zero[1]
    else
        return nlsolve((vout, vin)->vout .= func(vin,params...), ini).zero
    end
end
#nls(f,v,ini = u)
```

Numberって変数どこから出てきてる?
ボスのサイトにあったやつを丸写ししました

おｋ
http://www.cas.cmc.osaka-u.ac.jp/~paoon/misc/julia/post/nlsolve-package/
このサイトです
おｋ

これ結構時間かかる？
とりあえず、動いたみたいだけど、計算かえってこない
非線形方程式のソルバーなので時間はかかると思います．

いま動かしているコード
https://gist.github.com/muripoLife/500151545e6cb696a25314292d48a478

イテレーションごとにログ出ないの動いているかどうかわからんな
これはどこを書き換えましたか？

最期を
```
nls(g, ini=u)
```
にしただけ

更新されているかどうかわからない……
デバックしやすいように書きなおそ
確かに実行はされていそうですね

実行はされて、メモリーが5000mbかかっているところしかみてない。

いわゆるデバッグってどうするんですか

この場合、まずイテレーションごとに値プリントするとかかな？

今回は、gを高階関数で渡しているからかくの少し難しいかもね。

？？？

```
nls(g, ini=u)
```
gって関数でしょ。
関数を引数ように扱うのを高階関数というのよね。
合成関数とかがプログラミングで書きやすい的な。
なるほど

nlsの定義の中にうまく書き込んだらいいのかなと思うのですがどうでしょうか

それを試しているのだがうまく動かないね。

とりあえず長時間走らせてみます

ok

ワイもいろいろ書き換えてみるわ。

あとコード書いてるのjupyter?
そうです
30分やっても終わらないですね

おわらないね。
ちょいコードを考えておく。

```
using Images 
using FFTW
using NLsolve

#種々のパラメーターのセッティング
dt = 1.0
e = 0.8
C1 = 300.0
C2 = 150000.0
lam = 50000.0




#テストするための画像を作る．
u = zeros(32,32) #真っ黒な画像 
# make a view
img_a = @view u[1:32,8:13] #この行列を変化させると大元の行列も変えることができる．
img_a =  fill!(img_a, 1.0) 

img_b = @view u[1:32,20:24] 
img_b =  fill!(img_b, 1.0) 

img_c = @view u[1:32,27:30] 
img_c =  fill!(img_c, 1.0) 

img_d = @view u[16:20,1:32] 
img_d =  fill!(img_d, 0.5) 

 f = u #動かしたくない領域を参照するための行列
 #uを更新していく





function nabla1(u)
  fftu = fft(u)
    fftu[1,:] = fftu[1,: ] * 0
    for j in 1:15
        fftu[j+1,: ] =      (2 *pi  *j *1im /32 ) * fftu[j+1 ,: ]
        fftu[33-j ,: ] =  (2 *pi  *(-j) *1im /32 )*fftu[33-j ,: ]
   　 end #j
　real(ifft(fftu))
end

function nabla2(u)
  fftu = fft(u)
    fftu[: ,1 ] = fftu[: ,1 ] * 0
    for k in 1:15
        fftu[: ,k+1 ] =      (2 *pi  *k *1im /32 ) * fftu[: ,k+1 ]
        fftu[: ,33-k ] =  (2 *pi  *(-k) *1im /32 )*fftu[: ,33-k ]
     end #k
　real(ifft(fftu))
end

function lap(u) #ラプラシアンの微分演算子．
 nabla1(nabla1(u))  + nabla2(nabla2(u))
end
function lambda(u) #微分方程式に出てくるラムダの項を計算する関数．ok
    lamb =zeros(32,32)
        for i in 1:15
            for j in 1:32
                lamb[i,j] = lam * (f[i,j]-u[i,j])
            end #j
        end #i
        for i in 21:32
            for j in 1:32
                lamb[i,j] = lam * (f[i,j]-u[i,j])
            end #j
        end #i
    lamb
end #function


#NLsolveを用いる
function nls(func, params...; ini = [0.0])
    if typeof(ini) <: Number
        return nlsolve((vout,vin)->vout[1]=func(vin[1],params...), [ini]).zero[1]
    else
        return nlsolve((vout,vin)->vout .= func(vin,params...), ini).zero
    end
end
#nls(f,v,ini = u)

#非線形ソルバーで解くために方程式のスキームを書いた．vはuを1sステップ進めたもの
 g(v)=(v -u)/dt   + e * lap(lap(v)) - C1 * lap(v) + C2 * v -( lap((1/e) *4*(u^3) -6 * (u^2) + 2*u) + lambda(u) -C1 *(lap(u)) + C2*u)
 
 for i  in 1:10
    @time u = nls(g, ini = u)
end
 
colorview(Gray,u)
```

ちゃんと動くようにはなりましたが，値にほぼ変化がないです．