美しい角が作りたい！クロソイド曲線で曲率連続な曲線を作る。

今作ってるソレ、角の面取りは円で満足ですか？

嫌だと思う人は読み進めてください。

# はじめに

この記事ではコードから3Dモデルが作成できる[OpenSCAD](https://www.openscad.org/)を使います。
しかし、一般的な内容なのでほかのモデリングソフトでも同様のことができると思います。
Fusion360などのソフトには今から実装する内容が元から実装されているため、読む必要はありません。

## この記事に書いていること

- 円弧とクロソイド曲線の繋げ方
- クロソイド曲線で使われている積分をコンピュータで近似する方法

## この記事に書いていないこと

- OpenSCADの使い方
- 円弧のみの角を作ったときと等しい位置に円弧が来る曲線の作り方（僕も知りたいです）

## できるもの

![比較.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/271592/9a0bded7-86e0-be19-4cc1-f1a2c4c7ce80.png)

右の四角のほうが角の丸い部分と直線部分の境目が滑らかに繋がっているように見えませんか？
**注）** 先述の通り、左右では円弧の位置が異なるので(右側のほうが角が内側です)、公平な比較でないことをご容赦ください。

# 前提となる知識

## クロソイド曲線とは

詳しくは[Wikipedia クロソイド曲線](https://ja.wikipedia.org/wiki/%E3%82%AF%E3%83%AD%E3%82%BD%E3%82%A4%E3%83%89%E6%9B%B2%E7%B7%9A)に丸投げします。（そもそも僕も理解してません。フレネル積分って何？）

ここで知っておくべき性質は以下の3つです。

- 曲線の長さに比例して曲率が大きくなる $\Leftrightarrow$ 曲率半径が小さくなる。
- 曲率半径 $R$ と、始点からの距離 $L$ の積は一定。$RL=A^2\quad(Aは定数でクロソイドパラメータと呼ばれる)$
- 曲線の座標は媒介変数 $\theta$ と始点からの曲線長 $l$ を使って以下のように表せる。

```math
x(l) = \int_{0}^{l}\cos\frac{\theta^2}{2}\mathrm{d}\theta\\
y(l) = \int_{0}^{l}\sin\frac{\theta^2}{2}\mathrm{d}\theta
```

これらの性質を使うことで滑らかに直線部と角の円弧部を繋ぐことができます。

## 積分を近似する方法

**注）** 実は、僕もこれで良いのか分かってません。積分の定義通りなのでそこまで間違ったことはしてないと思いますが（結果も近いし）、数学的に良くないことをしていたら教えてください。

高校で習った区分求積法っぽいものを使います。ただし「限りなく小さい」をコンピュータで再現できないので「そこそこ小さい」を使って近似しています。

```math
s>0\:を「そこそこ小さい数」,\; 集合\,X=\{\, a+0,\: a+s,\: a+2s,\: \cdots,\: b \, \}\:とする.\\
\int_{a}^{b}f(x)\mathrm{d}x \;\simeq\; \sum_{x\in X}f(x)\times s
```

つまり、積分をグラフ上の細い長方形の面積の和として考えます。

# 計算

今回の最終目標である美しい角を作るためには、

1. 曲線（位置）が連続であること
- 傾きが連続であること
- 曲率が連続であること

の3つの条件を満たさなくてはなりません（G2連続）。

これらの条件を 傾き → 曲率 → 座標 の順に満たしていきます。

## 傾き連続

クロソイド曲線において、始点から $l$ 進んだ点の傾きは

```math
\begin{align}
\frac{y'(l)}{x'(l)}
&= \frac{\frac{\mathrm{d}}{\mathrm{d}\theta}\int_{0}^{l}\sin\frac{\theta^2}{2}\mathrm{d}\theta}{\frac{\mathrm{d}}{\mathrm{d}\theta}\int_{0}^{l}\cos\frac{\theta^2}{2}\mathrm{d}\theta}\\
&= \frac{\sin\frac{l^2}{2}}{\cos\frac{l^2}{2}}\\
&= \tan\frac{l^2}{2}
\end{align}
```

となる。これが円弧との境目の傾きと等しくなる。

**例）** 30度が円弧、左右30度がクロソイド曲線の合計90度の角。
繋ぎ目の角度が30度$=\pi/6$でないといけないので、

```math
\begin{align}
\tan\frac{l^2}{2}
&= \tan30^\circ\\
&= \tan\frac{\pi}{6}\\
\frac{l^2}{2}
&= \frac{\pi}{6}\\
\therefore \quad l &=\sqrt{\frac{\pi}{3}}
\end{align}
```

いろいろ同値でないのは承知していますが、$l$ は正が都合が良いなどで、ちゃんとうまくいきます。

|![接線連続.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/271592/76f954a5-6891-a1ab-3be8-0339ff3a58cd.png)|
|:-:|
|青が円弧、緑がクロソイド曲線。繋ぎ目の傾きが等しい。|

さて、円との境目で傾きが等しくなるクロソイド曲線が作れました。始点から $l'$ の場所の傾きが等しくなったとします。

## 曲率（曲率半径）連続

先ほどの曲線は、「始点から $l$ の場所で傾きが円と等しい」のでした。次は曲率を合わせたいわけです。ここで $RL=A^2$ を使います。
円の半径を $r$ とすると $r \times l' = A'^2$ となり、この $A'^2$ を相似比として先ほどの曲線を $A^2$ 倍します。

|![曲率連続.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/271592/026afb64-6f8f-256c-3223-b81841f06590.png)|
|:-:|
|先ほどの曲線を $A'^2$ 倍した。傾きに加えて曲率も等しい。|

## 位置連続

これは簡単で、どちらかを移動させて円弧の始点とクロソイド曲線の終点の座標を一致させるだけです。今までのサンプル画像は見やすいように既に位置連続になってました。

両側にクロソイド曲線がないと角全体としてG2連続になりませんから、反対側にも同じ曲線を繋げておきましょう。

|![反対側にも.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/271592/04548349-277a-5496-52ca-a01a3357e338.png)|
|:-:|
|反対側にも同じ曲線を付けた。G2連続の角。|

# 実装

では、OpenSCADで実際にこの角丸を実装します。
具体的な値は以下の通りとします。

- 円弧の半径 = 8 [mm]
- 円弧の角度 = 30$^\circ$
- 角全体の角度 = 90$^\circ$

計算の章で載せたサンプル画像と同じものです。

まずコードをすべて載せてしまいましょう。

```c
// 配列を逆順にする関数
function reverse(v) = [for (i = [0 : len(v)-1]) v[len(v)-1-i]];

function degree2radian(d) = d * PI / 180;

// ラジアンを入力にする三角関数。組み込みの三角関数の入力が度なので。
function sin_radian(d) = sin(180 * d / PI);
function cos_radian(d) = cos(180 * d / PI);

// クロソイド曲線の積分の内側部分
function clothoid(t) = [cos_radian(t * t / 2), sin_radian(t * t / 2)];

// 累積和。積分代わり
function cum_sum(v, _i=0,  _a=[]) =
    _i == len(v) ? _a :
    cum_sum(
        v, _i+1,
        concat(
            _a,
            [_i == 0 ? v[_i] : _a[len(_a)-1]+v[_i]]
        )
    );

// 角丸の本体
module continuous_corner(r, a = 90, circle_a=30) {
    n = 100; // 100分割
    clothoid_a = (a - circle_a) / 2;
    l = sqrt(degree2radian(clothoid_a * 2));

    _arc = [for (i = [0 : n]) let (t = clothoid_a + circle_a * i / n) [cos(t), sin(t)]];
    arc = _arc * r;

    dt = l / n;

    __c1 = [for (t = [0 : dt : l]) clothoid(t)] * dt;
    _c1 = cum_sum(concat([[0, 0]], reverse(__c1)));
    c1 = _c1 * l * r;

    d = arc[0][0] + (c1[len(c1)-1] - c1[0])[1];
    translate([-d, -d]) {
        color("#2196F3") polygon(arc);
        color("#4CAF50") {
            translate(arc[0]) rotate(-90) polygon(c1);
            mirror([-1, 1, 0]) translate(arc[0]) rotate(-90) polygon(c1);
        }
    }
}

continuous_corner(8);
```

コメントで粗方説明しているので、分かりにくいコードになっている部分をかいつまんで説明します。

`arc` 変数は円弧の部分の座標の配列です。今回は90度ある角の中央30度分なので、30度～60度の部分となります。

`c1` 変数はクロソイド曲線の座標の配列です。

```c
__c1 = [for (t = [0 : dt : l]) clothoid(t)] * dt;
```

この部分は、積分するために関数の値を配列にしています。さらに dt(そこそこ小さな値)を掛けることで、細長い長方形の面積が配列になっています。

```c
_c1 = cum_sum(concat([[0, 0]], reverse(__c1)));
```

次に、細長い長方形の面積を順に足し合わせていきます。この累積和が積分の近似になっているわけです。
途中 `reverse()` で配列の順序を逆にしていますが、本質的ではなく後で円と繋げるときに少し楽だからです。

```c
c1 = _c1 * l * r;
```

$RL=A^2$ で $A^2$ を相似比として全体をサイズ変更しています。ここでクロソイド曲線の端の曲率が円弧のそれと等しくなります。

# さいごに

$A$ は長さの次元を持つ量らしく、$A^2$ を使ってサイズ変更するところがなぜなのかよく分かっていません。正直、ごにょごにょしてたら上手くいったという感じです。

今回作ったような美しい曲線はAppleの製品に使われていることでも有名です。しかし、Appleのは角の出っ張り具合が円弧だけの場合と変わらない一段とすごいものなので、今後はそれが実装できるようになりたいと思います。

数学面・実装面ともに正確さに不安があるので、皆様からの愛あるマサカリをいつでも歓迎しています。何かありましたら是非コメントをよろしくお願いします。