Оптимизация моделирования n-тела с использованием Numpy и Glumpy - PullRequest
0 голосов
/ 21 ноября 2018

Я делаю процесс перехода от MatLab к Python и пытаюсь симулировать N-body в режиме реального времени как способ освоить язык.Я закончил сценарий, который использует клейкий для визуализации и векторизованный код для числовой интеграции, но он замедляет waaaaaay вниз более чем на 150 пунктов.По сути, я ищу советы о том, как оптимизировать код, чтобы он мог работать в режиме реального времени, по крайней мере, на порядок больше очков.В частности, есть ли способ реализовать интеграцию в шейдерный код с помощью графического процессора?

Мне известны ограничения с использованием метода Эйлера, который я использовал (и масштабирование ускорения), и я ищу информацию по оптимизации реализации программного обеспечения, в отличие от самого численного метода.

import numpy as np
from glumpy import app, gloo, gl

#VERTEX SHADER CODE
vertex = """
    attribute vec2 position;
    void main() {
         gl_PointSize = 5.0;
         gl_Position = vec4(position, 0.0, 1.0);

    } """

#FRAGMENT SHADER CODE
fragment = """
    void main() {
        gl_FragColor = vec4(vec3(0.5), 1.0);
    } """

#CREATE ARRAY OF POSITIONS
n_points = int(300)                         #number of points
p = np.random.uniform(-1,1,(n_points,2))    #array of positions (x,y)
v = np.zeros((n_points,2))                  #array of velocities (xdot,ydot)
m = np.ones((n_points,1))                   #array of masses
a = np.zeros((n_points,2))                  #array of accelerations
index = np.linspace(0,n_points-1,n_points)

#OPEN WINDOW AND CREATE SHADER PROGRAM
window = app.Window(512,512, color = (1,1,1,1))
points = gloo.Program(vertex,fragment,count = n_points)
points["position"] = p

@window.event
def on_draw(dt):
    global v, p, a, m, theta

    for i in range(n_points):    #for each point
        delta_x = p[index != i, 0] - p[i, 0]    #distance in x to each other point
        delta_y = p[index != i, 1] - p[i, 1]    #distance in y to each other point
        mass = m[index != i,:]    #mass of other points
        theta = np.arctan2(delta_y, delta_x)  #angle to other points
        r = (delta_x**2 + delta_y**2)*0.5 + 0.01  #distance to other points
        a[i,0] = np.sum(mass * np.cos(theta) / r**2)  #acceleration in x
        a[i,1] = np.sum(mass * np.sin(theta) / r**2)  #acceleration in y

    print(dt)
    a = a*0.00000001
    v = v + a*dt
    p = p + v*dt

    window.clear()
    points["position"] = p
    points.draw(gl.GL_POINTS)

app.run()
...