zoukankan      html  css  js  c++  java
  • 李代数E8 的根系 python绘图

    安装前需安装依赖:(针对Centos7)

    yum install -y cairo
    pip install cairocffi

    源代码:https://github.com/neozhaoliang/pywonderland/blob/master/src/misc/e8.py

    # -*- coding: utf-8 -*-
    """
    ~~~~~~~~~~~~~~
    The E8 picture
    ~~~~~~~~~~~~~~
    This script draws the picture of E8 projected to its Coxeter plane.
    For a detailed discussion of the math see Humphreys's book
        "Reflection Groups and Coxeter Groups", section 17, chapter 3.
    """
    from itertools import product, combinations
    import cairocffi as cairo
    import numpy as np
    
    
    COLORS = [(0.894, 0.102, 0.11),
              (0.216, 0.494, 0.72),
              (0.302, 0.686, 0.29),
              (0.596, 0.306, 0.639),
              (1.0, 0.5, 0),
              (1.0, 1.0, 0.2),
              (0.65, 0.337, 0.157),
              (0.97, 0.506, 0.75)]
    
    
    # --- step one: compute all roots and edges ---
    
    # There are 240 roots in the root system,
    # mutiply them by a factor 2 to be handy for computations.
    roots = []
    
    # Roots of the form (+-1, +-1, 0, 0, 0, 0, 0, 0),
    # signs can be chosen independently and the two non-zeros can be anywhere.
    for i, j in combinations(range(8), 2):
        for x, y in product([-2, 2], repeat=2):
            v = np.zeros(8)
            v[i] = x
            v[j] = y
            roots.append(v)
    
    # Roots of the form 1/2 * (+-1, +-1, ..., +-1), signs can be chosen
    # indenpendently except that there must be an even numer of -1s.
    for v in product([-1, 1], repeat=8):
        if sum(v) % 4 == 0:
            roots.append(v)
    roots = np.array(roots).astype(np.int)
    
    
    # Connect a root to its nearest neighbors,
    # two roots are connected if and only if they form an angle of pi/3.
    edges = []
    for i, r in enumerate(roots):
        for j, s in enumerate(roots[i+1:], i+1):
            if np.sum((r - s)**2) == 8:
                edges.append([i, j])
    
    
    # --- Step two: compute a basis of the Coxeter plane ---
    
    # A set of simple roots listed by rows of 'delta'
    delta = np.array([[1, -1, 0, 0, 0, 0, 0, 0],
                      [0, 1, -1, 0, 0, 0, 0, 0],
                      [0, 0, 1, -1, 0, 0, 0, 0],
                      [0, 0, 0, 1, -1, 0, 0, 0],
                      [0, 0, 0, 0, 1, -1, 0, 0],
                      [0, 0, 0, 0, 0, 1, 1, 0],
                      [-.5, -.5, -.5, -.5, -.5, -.5, -.5, -.5],
                      [0, 0, 0, 0, 0, 1, -1, 0]])
    
    # Dynkin diagram of E8:
    # 1---2---3---4---5---6---7
    #                 |
    #                 8
    # where vertex i is the i-th simple root.
    
    # The cartan matrix:
    cartan = np.dot(delta, delta.transpose())
    
    # Now we split the simple roots into two disjoint sets I and J
    # such that the simple roots in each set are pairwise orthogonal.
    # It's obvious to see how to find such a partition given the
    # Dynkin graph above: I = [1, 3, 5, 7] and J = [2, 4, 6, 8],
    # since roots are not connected by an edge if and only if they are orthogonal.
    # Then a basis of the Coxeter plane is given by
    # u = sum (c[i] * delta[i]) for i in I,
    # v = sum (c[j] * delta[j]) for j in J,
    # where c is an eigenvector for the minimal
    # eigenvalue of the Cartan matrix.
    eigenvals, eigenvecs = np.linalg.eigh(cartan)
    
    # The eigenvalues returned by eigh() are in ascending order
    # and the eigenvectors are listed by columns.
    c = eigenvecs[:, 0]
    u = np.sum([c[i] * delta[i] for i in [0, 2, 4, 6]], axis=0)
    v = np.sum([c[j] * delta[j] for j in [1, 3, 5, 7]], axis=0)
    
    # Gram-Schimdt u, v and normalize them to unit vectors.
    u /= np.linalg.norm(u)
    v = v - np.dot(u, v) * u
    v /= np.linalg.norm(v)
    
    
    # --- step three: project to the Coxeter plane ---
    roots_2d = [(np.dot(u, x), np.dot(v, x)) for x in roots]
    
    # Sort these projected vertices by their modulus in the coxter plane,
    # every successive 30 vertices form one ring in the resulting pattern,
    # assign these 30 vertices a same color.
    vertex_colors = np.zeros((len(roots), 3))
    modulus = np.linalg.norm(roots_2d, axis=1)
    ind_array = modulus.argsort()
    for i in range(8):
        for j in ind_array[30*i: 30*(i+1)]:
            vertex_colors[j] = COLORS[i]
    
    
    # --- step four: render to png image ---
    image_size = 600
    # The axis lie between [-extent, extent] x [-extent, extent]
    extent = 2.4
    linewidth = 0.0018
    markersize = 0.05
    
    surface = cairo.ImageSurface(cairo.FORMAT_RGB24, image_size, image_size)
    ctx = cairo.Context(surface)
    ctx.scale(image_size/(extent*2.0), -image_size/(extent*2.0))
    ctx.translate(extent, -extent)
    ctx.set_source_rgb(1, 1, 1)
    ctx.paint()
    
    for i, j in edges:
        x1, y1 = roots_2d[i]
        x2, y2 = roots_2d[j]
        ctx.set_source_rgb(0.2, 0.2, 0.2)
        ctx.set_line_width(linewidth)
        ctx.move_to(x1, y1)
        ctx.line_to(x2, y2)
        ctx.stroke()
    
    for i in range(len(roots)):
        x, y = roots_2d[i]
        color = vertex_colors[i]
        grad = cairo.RadialGradient(x, y, 0.0001, x, y, markersize)
        grad.add_color_stop_rgb(0, *color)
        grad.add_color_stop_rgb(1, *color/2)
        ctx.set_source(grad)
        ctx.arc(x, y, markersize, 0, 2*np.pi)
        ctx.fill()
    
    surface.write_to_png('E8.png')

    结果:

  • 相关阅读:
    开放就像死亡访问之后就能回头——Leo鉴书84
    将博客搬至CSDN
    将博客搬至CSDN
    滚动条
    Perl Pack写的一个数据报表程序
    利用hash 数组打印标题
    Linux显示只显示目录文件
    Linux显示按文件名降序文件
    Linux显示以时间生升序显示文件
    Linux显示按文件大小降序排列
  • 原文地址:https://www.cnblogs.com/bugutian/p/11098325.html
Copyright © 2011-2022 走看看