fortran66のブログ

fortran について書きます。

八元数の虚部と例外群 G2

例外群 G_2

コンパクト Lie 群の例外群 G_2八元数スカラーのノルムを保存する変換とみなせます。特に虚部のみをみると7次元空間での制限された回転になっています。これは四元数の虚部の変換が、3次元空間での回転に対応しているのと類比されます。

SO(3) では f 電子と関係する角運動量 3 に関わる 6j symbol で選択則的に許容なのに値が 0 になるものがある事と関係しています。
\begin{Bmatrix}
5\, 5\, 3\\
3\, 3\, 3
\end{Bmatrix}=0
rank 3 のテンソル演算子が消えるので 21-(2*3+1)=14。

定義式の数値的確認

G_2=\left\{ g\in\mathrm{GL(O)}| g(uv)=g(u)g(v), \forall u, v \in\mathrm{O} \right\}

定義式から虚部を比較することにします。7次元空間での制限された回転行列の具体形は、以下の卒論から借りてくることにします。

https://www.math.hmc.edu/seniorthesis/archives/2005/rarenas/rarenas-2005-thesis.pdf

SO(7) の場合、本来 7*6/2 = 21 自由度ありますが、ここから 7 自由度減って 14 自由度の変換になっています。

実行結果

左右両辺の計算結果が一致しています。

g(u)*g(v)   0.000  -0.325   0.199   0.179   0.400  -0.041   0.540   0.204
  g(uv)     0.000  -0.325   0.199   0.179   0.400  -0.041   0.540   0.204
g(u)*g(v)   0.000   0.005   0.000   0.178   0.564  -0.298   0.434   0.214
  g(uv)     0.000   0.005   0.000   0.178   0.564  -0.298   0.434   0.214
g(u)*g(v)   0.000   0.388   0.341   0.023   0.084  -0.022   0.090   0.624
  g(uv)     0.000   0.388   0.341   0.023   0.084  -0.022   0.090   0.624
g(u)*g(v)   0.000   0.102  -0.223  -0.446   0.485  -0.349   0.199   0.129
  g(uv)     0.000   0.102  -0.223  -0.446   0.485  -0.349   0.199   0.129
g(u)*g(v)   0.000  -0.066   0.295  -0.243  -0.666  -0.171  -0.001   0.222
  g(uv)     0.000  -0.066   0.295  -0.243  -0.666  -0.171  -0.001   0.222
g(u)*g(v)   0.000  -0.467  -0.309  -0.505   0.072  -0.074  -0.121   0.279
  g(uv)     0.000  -0.467  -0.309  -0.505   0.072  -0.074  -0.121   0.279
g(u)*g(v)   0.000  -0.052  -0.089  -0.331  -0.413   0.047  -0.570  -0.232
  g(uv)     0.000  -0.052  -0.089  -0.331  -0.413   0.047  -0.570  -0.232
g(u)*g(v)   0.000   0.674  -0.099   0.068   0.331   0.019  -0.292   0.093
  g(uv)     0.000   0.674  -0.099   0.068   0.331   0.019  -0.292   0.093
g(u)*g(v)   0.000  -0.097  -0.025   0.571   0.007  -0.283   0.499  -0.086
  g(uv)     0.000  -0.097  -0.025   0.571   0.007  -0.283   0.499  -0.086
g(u)*g(v)   0.000  -0.333   0.168   0.182   0.654  -0.005   0.079   0.255
  g(uv)     0.000  -0.333   0.168   0.182   0.654  -0.005   0.079   0.255
g(u)*g(v)   0.000   0.111  -0.037  -0.200  -0.560  -0.089  -0.499   0.217
  g(uv)     0.000   0.111  -0.037  -0.200  -0.560  -0.089  -0.499   0.217
続行するには何かキーを押してください . . .

ソース・プログラム

rot_1~7 それぞれ2変換角パラメータが合計14自由度を出します。ここでは面倒なので、rot_1..rot_7 の固定順番で、乱数で決めた角度で変換して、左右両辺の結果を比較しています。

八元数の積は、非結合演算なので演算順序に注意が必要です。

    module m_comp
      implicit none
      integer, parameter :: kd = kind(0.0d0)
      
      interface operator(+) 
        procedure :: add_r, add_c, add_h, add_o
      end interface
      interface operator(*) 
        procedure :: mul_r, mul_c, mul_h, mul_o
      end interface
      interface operator(/) 
        procedure :: scale_r, scale_c, scale_h, scale_o
        procedure :: div_r, div_c, div_h, div_o
      end interface
      interface operator(-) 
        procedure :: minus_r, minus_c, minus_h, minus_o
        procedure :: sub_r, sub_c, sub_h, sub_o
      end interface
      interface operator(.cc.) 
        procedure :: conj_r, conj_c, conj_h, conj_o
      end interface
      interface conj 
        procedure :: conj_r, conj_c, conj_h, conj_o
      end interface
      interface dot 
        procedure :: dot_r, dot_c, dot_h, dot_o
      end interface
      interface cross 
        procedure :: cross_r, cross_c, cross_h, cross_o
      end interface
      interface im 
        procedure :: im_r, im_c, im_h, im_o
      end interface

      type :: t_r
        real(kd) :: re
      end type t_r    

      type :: t_c
        type(t_r) :: re, im
      end type t_c    

      type :: t_h
        type(t_c) :: re, im
      end type t_h    

      type :: t_o
        type(t_h) :: re, im
      end type t_o    

    contains  
    
      type(t_r) pure function conj_r(a) result(c)
        type(t_r), intent(in) :: a
        c%re = a%re 
      end function conj_r
      type(t_r) pure function im_r(a) result(c)
        type(t_r), intent(in) :: a
        c%re = 0.0_kd
      end function im_r
      type(t_r) pure function minus_r(a) result(c)
        type(t_r), intent(in) :: a
        c%re = -a%re 
      end function minus_r
      type(t_r) pure function add_r(a, b) result(c)
        type(t_r), intent(in) :: a, b
        c%re = a%re + b%re
      end function add_r
      type(t_r) pure function sub_r(a, b) result(c)
        type(t_r), intent(in) :: a, b
        c%re = a%re - b%re
      end function sub_r
      type(t_r) pure function div_r(a, b) result(c)
        type (t_r), intent(in) :: a, b
        c = a * conj(b) / dot(b, b)
      end function div_r
      type(t_r) pure function mul_r(a, b) result(c)
        type (t_r), intent(in) :: a, b
        c%re = a%re * b%re
      end function mul_r
      type(t_r) pure function scale_r(a, b) result(c)
        type (t_r), intent(in) :: a
        real(kd), intent(in) :: b
        c%re = a%re / b
      end function scale_r
      real(kd) pure function dot_r(a, b) result(c)
        type (t_r), intent(in) :: a, b
        c = a%re * b%re
      end function dot_r
      type(t_r) pure function cross_r(a, b) result(c)
        type (t_r), intent(in) :: a, b
        c = (conj(b) * a - conj(a) * b) / 2.0_kd
      end function cross_r

      type(t_c) pure function conj_c(a) result(c)
        type(t_c), intent(in) :: a
        c%re =  a%re 
        c%im = -a%im 
      end function conj_c
      type(t_c) pure function im_c(a) result(c)
        type(t_c), intent(in) :: a
        c%re = im(a%re)
        c%im = a%im
      end function im_c
      type(t_c) pure function minus_c(a) result(c)
        type(t_c), intent(in) :: a
        c%re = -a%re 
        c%im = -a%im 
      end function minus_c
      type(t_c) pure function add_c(a, b) result(c)
        type (t_c), intent(in) :: a, b
        c%re = a%re + b%re
        c%im = a%im + b%im
      end function add_c
      type(t_c) pure function sub_c(a, b) result(c)
        type (t_c), intent(in) :: a, b
        c%re = a%re - b%re
        c%im = a%im - b%im
      end function sub_c
      type(t_c) pure function mul_c(a, b) result(c)
        type (t_c), intent(in) :: a, b
        c%re = a%re * b%re + (-b%im * a%im)
        c%im = b%im * a%re +   a%im * b%re
      end function mul_c
      type(t_c) pure function div_c(a, b) result(c)
        type (t_c), intent(in) :: a, b
        c = a * conj(b) / dot(b, b)
      end function div_c
      type(t_c) pure function to_c(r0, r1)
        real(kd), intent(in) :: r0, r1
        to_c = t_c(t_r(r0), t_r(r1))
      end function to_c  
      type(t_c) pure function scale_c(a, b) result(c)
        type (t_c), intent(in) :: a
        real(kd), intent(in) :: b
        c%re = a%re / b
        c%im = a%im / b
      end function scale_c
      real(kd) pure function dot_c(a, b) result(c)
        type (t_c), intent(in) :: a, b
        c = dot_r(a%re, b%re) + dot_r(a%im, b%im)
      end function dot_c
      type(t_c) pure function cross_c(a, b) result(c)
        type (t_c), intent(in) :: a, b
        c = (conj(b) * a - conj(a) * b) / 2.0_kd
      end function cross_c
      
      type(t_h) pure function conj_h(a) result(c)
        type(t_h), intent(in) :: a
        c%re = conj(a%re) 
        c%im =     -a%im 
      end function conj_h
      type(t_h) pure function im_h(a) result(c)
        type(t_h), intent(in) :: a
        c%re = im(a%re)
        c%im = a%im
      end function im_h
      type(t_h) pure function minus_h(a) result(c)
        type(t_h), intent(in) :: a
        c%re = -a%re 
        c%im = -a%im 
      end function minus_h
      type(t_h) pure function add_h(a, b) result(c)
        type (t_h), intent(in) :: a, b
        c%re = a%re + b%re
        c%im = a%im + b%im
      end function add_h
      type(t_h) pure function sub_h(a, b) result(c)
        type (t_h), intent(in) :: a, b
        c%re = a%re - b%re
        c%im = a%im - b%im
      end function sub_h
      type(t_h) pure function mul_h(a, b) result(c)
        type (t_h), intent(in) :: a, b
        c%re = a%re * b%re - conj(b%im) * a%im
        c%im = b%im * a%re + a%im * conj(b%re)
      end function mul_h
      type(t_h) pure function div_h(a, b) result(c)
        type (t_h), intent(in) :: a, b
        c = a * conj(b) / dot(b, b)
      end function div_h
      type(t_h) pure function to_h(r0, r1, r2, r3)
        real(kd), intent(in) :: r0, r1, r2, r3
        to_h = t_h(to_c(r0, r1), to_c(r2, r3))
      end function to_h
      type(t_h) pure function scale_h(a, b) result(c)
        type (t_h), intent(in) :: a
        real(kd), intent(in) :: b
        c%re = a%re / b
        c%im = a%im / b
      end function scale_h
      real(kd) pure function dot_h(a, b) result(c)
        type (t_h), intent(in) :: a, b
        c = dot_c(a%re, b%re) + dot_c(a%im, b%im)
      end function dot_h
      type(t_h) pure function cross_h(a, b) result(c)
        type (t_h), intent(in) :: a, b
        c = (conj(b) * a - conj(a) * b) / 2.0_kd
      end function cross_h
      type(t_h) pure function to_h_arr(a) result(c)
        real(kd), intent(in) :: a(4)
        c = to_h(a(1), a(2), a(3), a(4))
      end function to_h_arr
      
      type(t_o) pure function conj_o(a) result(c)
        type(t_o), intent(in) :: a
        c%re = conj(a%re) 
        c%im =     -a%im 
      end function conj_o
      type(t_o) pure function im_o(a) result(c)
        type(t_o), intent(in) :: a
        c%re = im(a%re)
        c%im = a%im
      end function im_o
      type(t_o) pure function minus_o(a) result(c)
        type(t_o), intent(in) :: a
        c%re = -a%re 
        c%im = -a%im 
      end function minus_o
      type(t_o) pure function add_o(a, b) result(c)
        type (t_o), intent(in) :: a, b
        c%re = a%re + b%re
        c%im = a%im + b%im
      end function add_o
      type(t_o) pure function sub_o(a, b) result(c)
        type (t_o), intent(in) :: a, b
        c%re = a%re - b%re
        c%im = a%im - b%im
      end function sub_o
      type(t_o) pure function mul_o(a, b) result(c)
        type (t_o), intent(in) :: a, b
        c%re = a%re * b%re - conj(b%im) * a%im
        c%im = b%im * a%re + a%im * conj(b%re)
      end function mul_o
      type(t_o) pure function div_o(a, b) result(c)
        type (t_o), intent(in) :: a, b
        c = a * conj(b) / dot(b, b)
      end function div_o
      type(t_o) pure function scale_o(a, b) result(c)
        type (t_o), intent(in) :: a
        real(kd), intent(in) :: b
        c%re = a%re / b
        c%im = a%im / b
      end function scale_o
      type(t_o) pure function cross_o(a, b) result(c)
        type (t_o), intent(in) :: a, b
        c = (conj(b) * a - conj(a) * b) / 2.0_kd
      end function cross_o
          
      type(t_o) pure function arr8_to_o(r)
        real(kd), intent(in) :: r(8)
        arr8_to_o = t_o(to_h(r(1), r(2), r(3), r(4)), to_h(r(5), r(6), r(7), r(8)))
      end function arr8_to_o

      type(t_o) pure function to_o(r0, r1, r2, r3, r4, r5, r6, r7)
        real(kd), intent(in) :: r0, r1, r2, r3, r4, r5, r6, r7
        to_o = t_o(to_h(r0, r1, r2, r3), to_h(r4, r5, r6, r7))
      end function to_o  
      type(t_o) pure function r7_to_o(a) result(b)
        real(kd), intent(in) :: a(7)
        b = t_o(to_h(0.0_kd, a(1), a(2), a(3)), to_h(a(4), a(5), a(6), a(7)))
      end function r7_to_o  
      pure function o_to_r7(a) result(b)
        type(t_o), intent(in) :: a
        real(kd) :: b(7)
        associate(hrcr=>a%re%re, hrci=>a%re%im, hicr=>a%im%re, hici=>a%im%im)
          b = [hrcr%im%re, hrci%re%re, hrci%im%re, hicr%re%re, hicr%im%re, hici%re%re, hici%im%re]
        end associate
      end function o_to_r7  
      real(kd) pure function dot_o(a, b) result(c)
        type (t_o), intent(in) :: a, b
        c = dot_h(a%re, b%re) + dot_h(a%im, b%im)
      end function dot_o
      
      type (t_o) pure function associator_o(a, b, c) result(d)
        type (t_o), intent(in) :: a, b, c
        d = (a * b) * c - a * (b * c)  
      end function associator_o

    end module m_comp

    
    module m_g2
      use m_comp
      implicit none
    contains
      pure function so7_1(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        real(kd), intent(in) :: a(7)
        real(kd) :: b(7), r(7, 7)
        r = 0.0
        r(1, 1) = 1.0
        r(2, 2) =  cos(x); r(2, 3) =  sin(x)
        r(3, 2) = -sin(x); r(3, 3) =  cos(x)
        r(4, 4) =  cos(x - y); r(4, 5) = -sin(x - y)
        r(5, 4) =  sin(x - y); r(5, 5) =  cos(x - y)
        r(6, 6) =  cos(y); r(6, 7) =  sin(y)
        r(7, 6) = -sin(y); r(7, 7) =  cos(y)
        b = matmul(r, a)
      end function so7_1 
      pure function so7_2(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        real(kd), intent(in) :: a(7)
        real(kd) :: b(7), r(7, 7)
        r = 0.0
        r(1, 1) =  cos(x); r(1, 3) = -sin(x)
        r(2, 2) = 1.0
        r(3, 1) =  sin(x); r(3, 3) =  cos(x)
        r(4, 4) =  cos(- x + y); r(4, 6) = sin(- x + y)
        r(5, 5) =  cos(y); r(5, 7) = -sin(y)
        r(6, 4) = -sin(- x + y); r(6, 6) = cos(- x + y)
        r(7, 5) =  sin(y); r(7, 7) =  cos(y)
        b = matmul(r, a)
      end function so7_2 
      pure function so7_3(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        real(kd), intent(in) :: a(7)
        real(kd) :: b(7), r(7, 7)
        r = 0.0
        r(1, 1) =  cos(x); r(1, 2) = sin(x)
        r(2, 1) = -sin(x); r(2, 2) = cos(x)
        r(3, 3) = 1.0
        r(4, 4) =  cos(- x + y); r(4, 7) = sin(- x + y)
        r(5, 5) =  cos(y); r(5, 6) = sin(y)
        r(6, 5) = -sin(y); r(6, 6) = cos(y)
        r(7, 4) = -sin(- x + y); r(7, 7) = cos(- x + y)
        b = matmul(r, a)
      end function so7_3 
      pure function so7_4(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        real(kd), intent(in) :: a(7)
        real(kd) :: b(7), r(7, 7)
        r = 0.0
        r(1, 1) =  cos(x); r(1, 5) = -sin(x)
        r(2, 2) =  cos(- x + y); r(2, 6) = -sin(- x + y)
        r(3, 3) =  cos(y); r(3, 7) =  sin(y)
        r(4, 4) = 1.0
        r(5, 1) =  sin(x); r(5, 5) =  cos(x)
        r(6, 2) =  sin(- x + y); r(6, 6) =  cos(- x + y)
        r(7, 3) = -sin(y); r(7, 7) =  cos(y)
        b = matmul(r, a)
      end function so7_4 
      pure function so7_5(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        real(kd), intent(in) :: a(7)
        real(kd) :: b(7), r(7, 7)
        r = 0.0
        r(1, 1) =  cos(x); r(1, 4) =  sin(x)
        r(2, 2) =  cos(x + y); r(2, 7) = sin(x + y)
        r(3, 3) =  cos(y); r(3, 6) =  sin(y)
        r(4, 1) = -sin(x); r(4, 4) =  cos(x)
        r(5, 5) = 1.0
        r(6, 3) = -sin(y); r(6, 6) =  cos(y)
        r(7, 2) = -sin(x + y); r(7, 7) =  cos(x + y)
        b = matmul(r, a)
      end function so7_5 
      pure function so7_6(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        real(kd), intent(in) :: a(7)
        real(kd) :: b(7), r(7, 7)
        r = 0.0
        r(1, 1) =  cos(x - y); r(1, 7) = -sin(x - y)
        r(2, 2) =  cos(x); r(2, 4) =  sin(x)
        r(3, 3) =  cos(y); r(3, 5) =  sin(y)
        r(4, 2) = -sin(x); r(4, 4) =  cos(x)
        r(5, 3) = -sin(y); r(5, 5) =  cos(y)
        r(6, 6) = 1.0
        r(7, 1) =  sin(x - y); r(7, 7) =  cos(x - y)
        b = matmul(r, a)
      end function so7_6 
      pure function so7_7(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        real(kd), intent(in) :: a(7)
        real(kd) :: b(7), r(7, 7)
        r = 0.0
        r(1, 1) =  cos(x); r(1, 6) = -sin(x)
        r(2, 2) =  cos(x - y); r(2, 5) = -sin(x - y)
        r(3, 3) =  cos(y); r(3, 4) = -sin(y)
        r(4, 3) =  sin(y); r(4, 4) =  cos(y)
        r(5, 2) =  sin(x - y); r(5, 5) =  cos(x - y)
        r(6, 1) =  sin(x); r(6, 6) =  cos(x)
        r(7, 7) = 1.0
        b = matmul(r, a)
      end function so7_7 
!      
!      
      type(t_o) pure function rot_1(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        type(t_o), intent(in) :: a
        b = r7_to_o(so7_1(x, y, o_to_r7(a)))
      end function rot_1
      type(t_o) pure function rot_2(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        type(t_o), intent(in) :: a
        b = r7_to_o(so7_2(x, y, o_to_r7(a)))
      end function rot_2
      type(t_o) pure function rot_3(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        type(t_o), intent(in) :: a
        b = r7_to_o(so7_3(x, y, o_to_r7(a)))
      end function rot_3
      type(t_o) pure function rot_4(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        type(t_o), intent(in) :: a
        b = r7_to_o(so7_4(x, y, o_to_r7(a)))
      end function rot_4
      type(t_o) pure function rot_5(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        type(t_o), intent(in) :: a
        b = r7_to_o(so7_5(x, y, o_to_r7(a)))
      end function rot_5
      type(t_o) pure function rot_6(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        type(t_o), intent(in) :: a
        b = r7_to_o(so7_6(x, y, o_to_r7(a)))
      end function rot_6
      type(t_o) pure function rot_7(x, y, a) result(b)
        real(kd), intent(in) :: x, y
        type(t_o), intent(in) :: a
        b = r7_to_o(so7_7(x, y, o_to_r7(a)))
      end function rot_7
    
    end module m_g2

    
    program supercomplex
      use m_comp
      use m_g2
      implicit none
      real(kd), parameter :: pi = 4 * atan(1.0)
      real(kd) :: theta, tx(7), ty(7)
      integer :: i
      type(t_o) :: ou, ov, uv, u0, v0, uv0
      ou = to_o(0.0_kd,  2.0_kd,-1.0_kd,  1.0_kd, -2.0_kd,  2.0_kd, -1.0_kd, -2.0_kd)
      ov = to_o(0.0_kd, -1.0_kd, 1.0_kd, -2.0_kd,  1.0_kd, -1.0_kd,  2.0_kd, -1.0_kd)
      u0 = ou / sqrt(dot(ou, ou))
      v0 = ov / sqrt(dot(ov, ov))
      uv0 = u0*v0 
      do i = 1, 100
        call random_number(tx)
        call random_number(ty)
        tx = tx * pi
        ty = ty * pi
        u0 = rot_1(tx(1), ty(1), u0)
        v0 = rot_1(tx(1), ty(1), v0)
        u0 = rot_2(tx(2), ty(2), u0)
        v0 = rot_2(tx(2), ty(2), v0)
        u0 = rot_3(tx(3), ty(3), u0)
        v0 = rot_3(tx(3), ty(3), v0)
        u0 = rot_4(tx(4), ty(4), u0)
        v0 = rot_4(tx(4), ty(4), v0)
        u0 = rot_5(tx(5), ty(5), u0)
        v0 = rot_5(tx(5), ty(5), v0)
        u0 = rot_6(tx(6), ty(6), u0)
        v0 = rot_6(tx(6), ty(6), v0)
        u0 = rot_7(tx(7), ty(7), u0)
        v0 = rot_7(tx(7), ty(7), v0)
        uv0 = rot_1(tx(1), ty(1), uv0)
        uv0 = rot_2(tx(2), ty(2), uv0)
        uv0 = rot_3(tx(3), ty(3), uv0)
        uv0 = rot_4(tx(4), ty(4), uv0)
        uv0 = rot_5(tx(5), ty(5), uv0)
        uv0 = rot_6(tx(6), ty(6), uv0)
        uv0 = rot_7(tx(7), ty(7), uv0)
        print '(a, 8f8.3)', 'g(u)*g(v)', u0 * v0 + to_o(dot(u0, v0), 0.0_kd, 0.0_kd, 0.0_kd, 0.0_kd, 0.0_kd, 0.0_kd, 0.0_kd)            
        print '(a, 8f8.3)', '  g(uv)  ', uv0            
      end do  

    end program supercomplex