問 :
三次元空間上に与えられた、点x1と、点x2と、それぞれの法線n1、n2から、
x1、x2を滑らかに繋いだ曲線を求め、その曲線上の点pを求めることを考える。
x1、x2は、二次補間するものとする。


答 :
問題を簡単にするため、まず二次元で考える。
二次方程式の式を、下のように定義する。

  y = -ax^2 + b                   ( 上に凸な放物線 )
    = -a ( x - x0 )( x + x0 )

  y'= -2ax

とすると、点( x0 , 0 ) ( -x0 , 0 )での法線ベクトルは、
傾きがそれぞれ、( 1 , -2ax0 ) ( 1 , 2ax0 )であることから、
( 2ax0 , 1 ) ( -2ax0 , 1 )となる。

この法線を、N1 = ( 2ax0 , 1 )、N2 = ( -2ax0 , 1 )とすると、

  cos = ( N1 * N2 )/sqrt( N1^2 * N2^2 )
      = ( -4 a^2 x0^2 + 1)/( 4 a^2 x0^2 + 1 )

  sin = sqrt( 1 - cos^2 )
      = sqrt( |N1|^2 * |N2|^2 - ( N1 * N2 )^2 )/sqrt( |N1|^2 * |N2|^2 )
      = 4a/( 4 a^2 x0^2 + 1 )

簡略化のため、x0 = 1にスケーリングすると、

  cos = ( -4 a'^2 + 1)/( 4 a'^2 + 1 )

  sin = 4a'/( 4 a'^2 + 1 )

  cos/sin = ( -4 a'^2 + 1)/4a'

  z = cos/sin と置くと

  z = ( -4 a'^2 + 1)/4a'
  4a'z = -4 a'^2 + 1
  4 a'^2 + 4a'z - 1 = 0

  a'= ( -2z +- sqrt( 4z^2 + 4 ) )/4
    = ( -z +- sqrt( z^2 + 1 ) )/2

  z = cos/sin
    = ( N1 . N2 )/sqrt( |N1|^2 * |N2|^2 - ( N1 . N2 )^2 )

  zをa'に代入して、

  a'= ( -z +- sqrt( z^2 + 1 ) )/2
    = ( -( N1 . N2 )/sqrt( |N1|^2 * |N2|^2 - ( N1 . N2 )^2 )
        +- sqrt( ( N1 . N2 )^2/( |N1|^2 * |N2|^2 - ( N1 . N2 )^2 ) + 1 ) )/2
    = ( -( N1 . N2 )/sqrt( |N1|^2 * |N2|^2 - ( N1 . N2 )^2 )
        +- sqrt( ( |N1|^2 * |N2|^2 )/( |N1|^2 * |N2|^2 - ( N1 . N2 )^2 ) ) )/2

    = ( -( N1 . N2 ) +- ( |N1|^2 * |N2|^2 ) )/ 2 * sqrt( |N1|^2 * |N2|^2 - ( N1 . N2 )^2 )

  a'>=0であるとすると、

    = ( -( N1 . N2 ) + ( |N1|^2 * |N2|^2 ) )/ 2 * sqrt( |N1|^2 * |N2|^2 - ( N1 . N2 )^2 )
    = sqrt( |N1|^2 * |N2|^2 - ( N1 . N2 )^2 )/2

    = sqrt( 1 - ( N1 * N2 )^2 )/2  ( N1,N2を長さ1とする )

ここで、三次元でこれを考え直すと、Ｘ軸を辺の方向、Ｙ軸を辺の法線方向と考えることが
できることから、
  
  P = ( X1 + X2 )/2 + b' * Np  ( Npは、(N1 + N2)/2を正規化したもの )

スケーリングを元に戻して、b'を算出すると、
  y/x0 = -a'x^2 + b'
       = -a'( x'/x0 - 1 )( x'/x0 + 1 )
       = -a'(x'/x0)^2 + a
  y    = -a'/x0 * x'^2 + a'x0
  b'   = a'x0

  P = ( X1 + X2 )/2 + a' * l/2 * Np  ( lはX1からX2までの長さ
                                       NpはPでの長さ1の法線ベクトル )

a'は前述したとおりであるが、それは二次元平面上を想定しているので、それを考慮し、
N1、N2を、NpとX1 -> X2を軸とする平面に展開してみる。

  N1' = ( (N1 . L)/l , N1 . Np )
  N2' = ( (N2 . L)/l , N2 . Np )   ( L = X2 - X1 , l = |L|とする )

これをa'に代入する。

  a' = sqrt( |N1'|^2 * |N2'|^2 - ( N1' . N2' )^2 )/2
     = sqrt( ( ( N1 . L )^2 / l^2 + ( N1 . Np )^2 )( ( N2 . L )^2 / l^2 + ( N2 . Np )^2 )
             - ( ( N1 . L )( N2 . L ) / l^2 + ( N1 . Np )( N2 . Np ) )^2 )/2
     = sqrt( ( ( N1 . L )^2 ( N2 . Np )^2 ) / l^2 + ( ( N1 . Np )^2 ( N2 . L )^2 ) / l^2
             - 2( N1 . L )( N2 . L )( N1 . Np )( N2 . Np ) / l^2 )/2
     = sqrt( ( ( N1 . L )( N2 . Np ) - ( N2 . L )( N1 . Np ) )^2 )/2l
     = ( ( N1 . L )( N2 . Np ) - ( N2 . L )( N1 . Np ) )/2l

ここで、a'はプラスしか想定してないので、Npにその方向を含ませるように
できないかを考える。方向を有効にするには、外積は有効であるので、
N1'とN2'の外積を考えると、

  N1' x N2' = ( ( N1 . L )( N2 . Np ) - ( N2 . L )( N1 . Np ) )/l

となり、

  a' = ( N1' x N2' )/2 となる。

N1' x N2'は、最初の放物線の定義より、N1'とN2'が上に開いたようになっていると－、
閉じたようになっていると＋になるので、

  P = ( X1 + X2 )/2 + a' * l/2 * Np
    = ( X1 + X2 )/2 - ( N1' x N2' )/2 * l/2 * Np
    = ( X1 + X2 )/2 - ( ( N1 . L )( N2 . Np ) - ( N2 . L )( N1 . Np ) )/l * l/4 * Np
    = ( X1 + X2 )/2 - ( ( N1 . L )( N2 . Np ) - ( N2 . L )( N1 . Np ) )/4 * Np
    = ( X1 + X2 )/2 - ( ( N1 . L )( N2 . Np' ) - ( N2 . L )( N1 . Np' ) )/4 * Np' / |Np'|^2

  ( L = X2 -X1、Np'は、Npを正規化しなかったもの )