【MATLAB】Cody5:Hard なるものにチャレンジしてみた

# はじめに
研究室で [MATLAB](https://jp.mathworks.com/products/matlab.html) に触れることが多くなったので, MATLABer 内輪で行われているコードゲーム[^1] である [Cody](https://www.mathworks.com/matlabcentral/cody/) に挑戦してみました。
# [MATLAB](https://jp.mathworks.com/products/matlab.html) とは
MATLAB は, **直観的な方法での数値の解析に長けているプログラミング言語**になります。**豊富な組み込み関数, 広範な学術分野の拡張機能**が特徴です。

一時期 NASA も使ってるという宣伝文句が使われてましたが, [NASA の Open Source](https://code.nasa.gov/) を見てもかなりの割合で MATLAB ベースのものがあるので, これは現在についても本当かと。あと, [神経科学の分野でもデファクトだという記事](https://special.nikkeibp.co.jp/atclh/NXT/18/mathworks1206/)もありました。

MATLAB のドキュメンテーションはコーダーではなく, 一般の科学者向けの説明になっているため, 使う人によっては用語に違和感があるかと思われます。短所を挙げるなら, インタラクティブな機能を前提としてないためか, 描画機能は低速と感じます。（あと, ライセンスが有料。私は大学がライセンスを与えてくれたので使えています。）
# [Cody](https://www.mathworks.com/matlabcentral/cody/) とは
MATLABer が内輪でやってるコードゲームです。~~需要は正直わかりかねますが, 研究員が研究の暇つぶしにやってたりするんですかね。~~ほぼすべての大問が英語になっており, 題意の関数を書ければ AC という無難なつくりです。
# [Cody5:Hard](https://www.mathworks.com/matlabcentral/cody/groups/35) の第１問を解いてみる
[Cody5:Hard](https://www.mathworks.com/matlabcentral/cody/groups/35) は Cody の 5 周年を記念して作られた問題集になります。MATLABer は理数系が得意な人が集まる印象なので, これまた難しいんでしょうか。

>It's Cody's 5th birthday, and you've been tasked with putting the candles on the cake. Your goal is to maximize the distance between the 5 required candles.

>Given a rectangular cake with specified length and width, return the highest possible minimum distance between any two candles.

>**Important notes**

>You may assume that a candle can be placed directly on the edge of the cake (even though that would be physically challenging!).
Required tolerance: <0.01

要は, **「長方形の２辺の長さが与えられるので, なるべく各点の距離が離れるように 5 つの点を打て」**ってことですね。返すのは 2 点の距離の最小値になります。回答は下図を見れば理解できるでしょう。

$\mathrm{mn}$ が短辺, $\mathrm{mx}$ が長辺になります。
![birthday2.jpg](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/577202/4ae27431-b15e-c7ea-a1d8-b5e3d84cd21e.jpeg)

以下が回答のコードです。geometry は正直苦手なので, 上の図が書けるまでちょっと時間がかかりました。

```matlab:birthdaycandles.m
function ans = birthdaycandles(len,wid)
    mx = max(len,wid);
    mn = min(len,wid);
    l1 = mn * sqrt(3);
    l2 = l1 * 4/3;
    if mx <= l1
        hypot(mn,mx)/2
    elseif mx <= l2
        2*sqrt(mn^2 + (mx/2)^2 - mn*(mx/2)*sqrt(3))
    else
        hypot(mn,mx/4)
    end
end
```
## 所望の値のふるまい
折角なので, 所望の値がどのようにふるまうか, プロットして確認してみましょう。

```matlab:birthdayplot1.m
dist = zeros(1000);
parfor i = 1:1000
    for j = 1:1000
        mx = max(i,j);
        mn = min(i,j);
        l1 = mn * sqrt(3);
        l2 = l1 * 4/3;
        if mx <= l1
            dist(i,j) = hypot(mn,mx)/2;
        elseif mx <= l2
            dist(i,j) = 2*sqrt(mn^2 + (mx/2)^2 - mn*(mx/2)*sqrt(3));
        else
            dist(i,j) = hypot(mn,mx/4);
        end
    end
end
imagesc(dist);
hold on
contour(dist,20,'--black');
```
`parfor` に関しては [Parallel Computing Toolbox](https://jp.mathworks.com/help/parallel-computing/index.html?s_tid=CRUX_lftnav) をインストールしている場合に, 自動的にタスクを分割し CPU 上で並列ループ計算を実行できます。

結果は下図で, 左上の頂点が座標と共通の長方形を書いたら, 右下の頂点が所望の値に位置するようになってます。
![birthdayplot1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/577202/91ce4ae1-57a8-d71b-0983-9770ec5f2097.png)
ちょっと情報量が多い気がするので, 長辺を固定してプロットしてみましょう。

```matlab:birthdayplot2.m
hold off
plot(dist(1000,:))
title('long side = 1000')
xticks(0:200:1000)
xticklabels(0:0.2:1)
```
結果は下図。
![birthdayplot2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/577202/66fb8648-0cc9-b521-940c-c16c1ceae75b.png)
図から, ケーキは $1:\sqrt3$の割合だと一番効率よく 5 本のロウソクが置けるのでは（？）と推測されました（無理やりな結論）。

# おわりに
この拙文で、なんとなく Cody のノリが理解いただけたら幸いに思います。MATLAB を普段使っていらしてる方は力試しになるのではないでしょうか。（AtCoder でも使えたらいいなとか思いつつ）

なんだかんだこのアカウントで記事を書くのも初めてなので, 生存報告代わりに更新していこうかな:sleeping:


[^1]: AtCoder などは競技としてチートなどができぬようしっかり設計されていますが、Cody はそんなにガチではなく、ほぼすべての問題でチートがトップを占めてます。やるなら気にしない姿勢が大事。