在 Pygame 中使用 numpy 进行更高效的风洞模拟
more efficient wind tunnel simulation in Pygame, using numpy
我是一名航空学生，正在为我们的 python 编程课程开展一个学校项目。作业是仅使用 Pygame 和 numpy 创建一个程序。我决定创建一个风洞模拟来模拟二维机翼上的气流。我想知道从编程的角度来看是否有更有效的计算方法。我来解释程序：
我在这里附上了一张图片：
（稳态）流场使用涡面板法建模。基本上，我使用的是 Nx 乘以 Ny 点的网格，其中每个点都有一个速度 (u,v) 向量。然后使用 Pygame 我将这些网格点映射为圆圈，这样它们就像一个影响区域。网格点是下图中的灰色圆圈：
我创建了 N 个粒子并通过如下迭代确定它们的速度：
创建粒子列表。
创建网格列表。
对于网格列表中的每个网格点：
对于粒子列表中的每个粒子：
若粒子A在网格点n(xn,yn)的影响范围内：
粒子A其速度=网格点n处的速度
可视化 Pygame 中的所有内容。
这种基本方法是我能想到的在 Pygame 中可视化流程的唯一方法。模拟效果很好，但是如果我增加网格点的数量（增加流场的精度），性能就会下降。我的问题是是否有更有效的方法来仅使用 pygame 和 numpy?
我在这里附上了代码：
import pygame,random,sys,numpy
from Flow import Compute
from pygame.locals import *
import random, math, sys
#from PIL import Image
pygame.init()
Surface = pygame.display.set_mode((1000,600))
#read the airfoil geometry from a dat file
with open ('./resources/naca0012.dat') as file_name:
x, y = numpy.loadtxt(file_name, dtype=float, delimiter='\t', unpack=True)
#parameters used to describe the flow
Nx=30# 30 column grid
Ny=10#10 row grid
N=20#number of panels
alpha=0#angle of attack
u_inf=1#freestream velocity
#compute the flow field
u,v,X,Y= Compute(x,y,N,alpha,u_inf,Nx,Ny)
#The lists used for iteration
Circles = []
Particles= []
Velocities=[]
#Scaling factors used to properly map the potential flow datapoints into Pygame
magnitude=400
vmag=30
umag=30
panel_x= numpy.multiply(x,magnitude)+315
panel_y= numpy.multiply(-y,magnitude)+308
#build the grid suited for Pygame
grid_x= numpy.multiply(X,magnitude)+300
grid_y= numpy.multiply(Y,-1*magnitude)+300
grid_u =numpy.multiply(u,umag)
grid_v =numpy.multiply(v,-vmag)
panelcoordinates= zip(panel_x, panel_y)
# a grid area
class Circle:
def __init__(self,xpos,ypos,vx,vy):
self.radius=16
self.x = xpos
self.y = ypos
self.speedx = 0
self.speedy = 0
#create the grid list
for i in range(Ny):
for s in range(Nx):
Circles.append(Circle(int(grid_x[i][s]),int(grid_y[i][s]),grid_u[i][s],grid_v[i][s]))
Velocities.append((grid_u[i][s],grid_v[i][s]))
#a particle
class Particle:
def __init__(self,xpos,ypos,vx,vy):
self.image = pygame.Surface([10, 10])
self.image.fill((150,0,0))
self.rect = self.image.get_rect()
self.width=4
self.height=4
self.radius =2
self.x = xpos
self.y = ypos
self.speedx = 30
self.speedy = 0
#change particle velocity if collision with grid point
def CircleCollide(Circle,Particle):
Particle.speedx = int(Velocities[Circles.index((Circle))][0])
Particle.speedy = int(Velocities[Circles.index((Circle))][1])
#movement of particles
def Move():
for Particle in Particles:
Particle.x += Particle.speedx
Particle.y += Particle.speedy
#create particle streak
def Spawn(number_of_particles):
for i in range(number_of_particles):
i=i*(300/number_of_particles)
Particles.append(Particle(0, 160+i,1,0))
#create particles again if particles are out of wake
def Respawn(number_of_particles):
for Particle in Particles:
if Particle.x >1100:
Particles.remove(Particle)
if Particles==[]:
Spawn(number_of_particles)
#Collsion detection using pythagoras and distance formula
def CollisionDetect():
for Circle in Circles:
for Particle in Particles:
if Particle.y >430 or Particle.y<160:
Particles.remove(Particle)
if math.sqrt( ((Circle.x-Particle.x)**2) + ((Circle.y-Particle.y)**2) ) <= (Circle.radius+Particle.radius):
CircleCollide(Circle,Particle)
#draw everything
def Draw():
Surface.fill((255,255,255))
#Surface.blit(bg,(-300,-83))
for Circle in Circles:
pygame.draw.circle(Surface,(245,245,245),(Circle.x,Circle.y),Circle.radius)
for Particle in Particles:
pygame.draw.rect(Surface,(150,0,0),(Particle.x,Particle.y,Particle.width,Particle.height),0)
#pygame.draw.rect(Surface,(245,245,245),(Circle.x,Circle.y,1,16),0)
for i in range(len(panelcoordinates)-1):
pygame.draw.line(Surface,(0,0,0),panelcoordinates[i],panelcoordinates[i+1],3)
pygame.display.flip()
def GetInput():
keystate = pygame.key.get_pressed()
for event in pygame.event.get():
if event.type == QUIT or keystate[K_ESCAPE]:
pygame.quit();sys.exit()
def main():
#bg = pygame.image.load("pressure.png")
#bg = pygame.transform.scale(bg,(1600,800))
#thesize= bg.get_rect()
#bg= bg.convert()
number_of_particles=10
Spawn(number_of_particles)
clock = pygame.time.Clock()
while True:
ticks = clock.tick(60)
GetInput()
CollisionDetect()
Move()
Respawn(number_of_particles)
Draw()
if __name__ == '__main__': main()
该代码需要另一个脚本来计算流场本身。它还从文本文件中读取数据点以获得机翼的几何形状。
我没有提供这两个文件，但如果需要我可以添加它们。提前谢谢你。
我整理了您的代码并进行了一些更改，即为您的 类 添加范围并引入更多内容。没有 Flow 的进一步知识，我无法对此进行全面测试，但如果你能回复我，我可以做更多的事情。我在这里假设 'flow field' 可以通过 numpy.meshgrid 函数模拟。
import pygame,numpy,sys
import pygame.locals
import math
class Particle:
def __init__(self,xpos,ypos,vx,vy):
self.size = numpy.array([4,4])
self.radius =2
self.pos = numpy.array([xpos,ypos])
self.speed = numpy.array([30,0])
self.rectangle = numpy.hstack((self.pos,self.size))
def move(self):
self.pos += self.speed
self.rectangle = numpy.hstack((self.pos,self.size))
def distance(self,circle1):
return math.sqrt(numpy.sum((circle1.pos - self.pos)**2))
def collision(self,circle1):
result = False
if self.pos[1] >430 or self.pos[1]<160:
result = True
if self.distance(circle1) <= (circle1.radius+self.radius):
self.speed = circle1.speed
return result
class Particles:
def __init__(self,num_particles):
self.num = num_particles
self.particles =[]
self.spawn()
def spawn(self):
for i in range(self.num):
i=i*(300/self.num)
self.particles.append(Particle(0, 160+i,1,0))
def move(self):
for particle in self.particles:
particle.move()
if particle.pos[0] >1100:
self.particles.remove(particle)
if not self.particles: self.spawn()
def draw(self):
for particle in self.particles:
pygame.draw.rect(Surface,(150,0,0),particle.rectangle,0)
def collisiondetect(self,circle1):
for particle in self.particles:
if particle.collision(circle1):
self.particles.remove(particle)
def GetInput():
keystate = pygame.key.get_pressed()
for event in pygame.event.get():
if event.type == pygame.locals.QUIT or keystate[pygame.locals.K_ESCAPE]:
pygame.quit()
sys.exit()
#draw everything
def Draw(sw,cir):
Surface.fill((255,255,255))
cir.draw()
for i in range(panelcoordinates.shape[1]):
pygame.draw.line(Surface,(0,0,0),panelcoordinates[0,i-1],panelcoordinates[0,i],3)
sw.draw()
pygame.display.flip()
# a grid area
class Circle:
def __init__(self,xpos,ypos,vx,vy):
self.radius=16
self.pos = numpy.array([xpos,ypos])
self.speed = numpy.array([vx,vy])
class Circles:
def __init__(self,columns,rows):
self.list = []
grid_x,grid_y = numpy.meshgrid(numpy.linspace(0,1000,columns),numpy.linspace(200,400,rows))
grid_u,grid_v = numpy.meshgrid(numpy.linspace(20,20,columns),numpy.linspace(-1,1,rows))
for y in range(rows):
for x in range(columns):
c1= Circle(int(grid_x[y,x]),int(grid_y[y,x]),grid_u[y,x],grid_v[y,x])
self.list.append(c1)
def draw(self):
for circle in self.list:
pygame.draw.circle(Surface,(245,245,245),circle.pos,circle.radius)
def detectcollision(self,parts):
for circle in self.list:
parts.collisiondetect(circle)
if __name__ == '__main__':
#initialise variables
number_of_particles=10
Nx=30
Ny=10
#circles and particles
circles1 = Circles(Nx,Ny)
particles1 = Particles(number_of_particles)
#read the airfoil geometry
panel_x = numpy.array([400,425,450,500,600,500,450,425,400])
panel_y = numpy.array([300,325,330,320,300,280,270,275,300])
panelcoordinates= numpy.dstack((panel_x,panel_y))
#initialise PyGame
pygame.init()
clock = pygame.time.Clock()
Surface = pygame.display.set_mode((1000,600))
while True:
ticks = clock.tick(6)
GetInput()
circles1.detectcollision(particles1)
particles1.move()
Draw(particles1,circles1)
我也粗略地画了一个机翼，同样我不了解带有坐标的数据文件。
代码中的一个瓶颈可能是碰撞检测。 CollisionDetect() 计算每个粒子和每个圆之间的距离。然后，如果检测到碰撞，CircleCollide() 会在 Circles 中找到圆的索引（线性搜索），以便可以从 Velocities 中的相同索引中检索速度。显然，改进的时机已经成熟。
首先，Circleclass已经有了speedx/speedy属性中的速度，所以Velocities可以去掉。
其次，因为圆在固定的位置，你可以从粒子的位置计算出哪个圆最接近任何给定的粒子。
# You may already have these values from creating grid_x etc.
# if not, you only need to calculated them once, because the
# circles don't move
circle_spacing_x = Circles[1].x - Circles[0].x
circle_spacing_y = Circles[Nx].y - Circles[0].y
circle_first_x = Circles[0].x - circle_spacing_x / 2
circle_first_y = Circles[0].y - circle_spacing_y / 2
那么CollisionDetect()就变成了：
def CollisionDetect():
for particle in Particles:
if particle.y >430 or particle.y<160:
Particles.remove(particle)
continue
c = (particle.x - circle_first_x) // circle_spacing_x
r = (particle.y - circle_first_y) // circle_spacing_y
circle = Circles[r*Nx + c]
if ((circle.x - particle.x)**2 + (circle.y - particle.y)**2
<= (circle.radius+particle.radius)**2):
particle.speedx = int(circle.speedx)
particle.speedy = int(circle.speedy)
我是一名航空学生，正在为我们的 python 编程课程开展一个学校项目。作业是仅使用 Pygame 和 numpy 创建一个程序。我决定创建一个风洞模拟来模拟二维机翼上的气流。我想知道从编程的角度来看是否有更有效的计算方法。我来解释程序：
我在这里附上了一张图片：
（稳态）流场使用涡面板法建模。基本上，我使用的是 Nx 乘以 Ny 点的网格，其中每个点都有一个速度 (u,v) 向量。然后使用 Pygame 我将这些网格点映射为圆圈，这样它们就像一个影响区域。网格点是下图中的灰色圆圈：
我创建了 N 个粒子并通过如下迭代确定它们的速度：
创建粒子列表。
创建网格列表。
对于网格列表中的每个网格点：
对于粒子列表中的每个粒子：
若粒子A在网格点n(xn,yn)的影响范围内：
粒子A其速度=网格点n处的速度
可视化 Pygame 中的所有内容。
这种基本方法是我能想到的在 Pygame 中可视化流程的唯一方法。模拟效果很好，但是如果我增加网格点的数量（增加流场的精度），性能就会下降。我的问题是是否有更有效的方法来仅使用 pygame 和 numpy?
我在这里附上了代码：
import pygame,random,sys,numpy
from Flow import Compute
from pygame.locals import *
import random, math, sys
#from PIL import Image
pygame.init()
Surface = pygame.display.set_mode((1000,600))
#read the airfoil geometry from a dat file
with open ('./resources/naca0012.dat') as file_name:
x, y = numpy.loadtxt(file_name, dtype=float, delimiter='\t', unpack=True)
#parameters used to describe the flow
Nx=30# 30 column grid
Ny=10#10 row grid
N=20#number of panels
alpha=0#angle of attack
u_inf=1#freestream velocity
#compute the flow field
u,v,X,Y= Compute(x,y,N,alpha,u_inf,Nx,Ny)
#The lists used for iteration
Circles = []
Particles= []
Velocities=[]
#Scaling factors used to properly map the potential flow datapoints into Pygame
magnitude=400
vmag=30
umag=30
panel_x= numpy.multiply(x,magnitude)+315
panel_y= numpy.multiply(-y,magnitude)+308
#build the grid suited for Pygame
grid_x= numpy.multiply(X,magnitude)+300
grid_y= numpy.multiply(Y,-1*magnitude)+300
grid_u =numpy.multiply(u,umag)
grid_v =numpy.multiply(v,-vmag)
panelcoordinates= zip(panel_x, panel_y)
# a grid area
class Circle:
def __init__(self,xpos,ypos,vx,vy):
self.radius=16
self.x = xpos
self.y = ypos
self.speedx = 0
self.speedy = 0
#create the grid list
for i in range(Ny):
for s in range(Nx):
Circles.append(Circle(int(grid_x[i][s]),int(grid_y[i][s]),grid_u[i][s],grid_v[i][s]))
Velocities.append((grid_u[i][s],grid_v[i][s]))
#a particle
class Particle:
def __init__(self,xpos,ypos,vx,vy):
self.image = pygame.Surface([10, 10])
self.image.fill((150,0,0))
self.rect = self.image.get_rect()
self.width=4
self.height=4
self.radius =2
self.x = xpos
self.y = ypos
self.speedx = 30
self.speedy = 0
#change particle velocity if collision with grid point
def CircleCollide(Circle,Particle):
Particle.speedx = int(Velocities[Circles.index((Circle))][0])
Particle.speedy = int(Velocities[Circles.index((Circle))][1])
#movement of particles
def Move():
for Particle in Particles:
Particle.x += Particle.speedx
Particle.y += Particle.speedy
#create particle streak
def Spawn(number_of_particles):
for i in range(number_of_particles):
i=i*(300/number_of_particles)
Particles.append(Particle(0, 160+i,1,0))
#create particles again if particles are out of wake
def Respawn(number_of_particles):
for Particle in Particles:
if Particle.x >1100:
Particles.remove(Particle)
if Particles==[]:
Spawn(number_of_particles)
#Collsion detection using pythagoras and distance formula
def CollisionDetect():
for Circle in Circles:
for Particle in Particles:
if Particle.y >430 or Particle.y<160:
Particles.remove(Particle)
if math.sqrt( ((Circle.x-Particle.x)**2) + ((Circle.y-Particle.y)**2) ) <= (Circle.radius+Particle.radius):
CircleCollide(Circle,Particle)
#draw everything
def Draw():
Surface.fill((255,255,255))
#Surface.blit(bg,(-300,-83))
for Circle in Circles:
pygame.draw.circle(Surface,(245,245,245),(Circle.x,Circle.y),Circle.radius)
for Particle in Particles:
pygame.draw.rect(Surface,(150,0,0),(Particle.x,Particle.y,Particle.width,Particle.height),0)
#pygame.draw.rect(Surface,(245,245,245),(Circle.x,Circle.y,1,16),0)
for i in range(len(panelcoordinates)-1):
pygame.draw.line(Surface,(0,0,0),panelcoordinates[i],panelcoordinates[i+1],3)
pygame.display.flip()
def GetInput():
keystate = pygame.key.get_pressed()
for event in pygame.event.get():
if event.type == QUIT or keystate[K_ESCAPE]:
pygame.quit();sys.exit()
def main():
#bg = pygame.image.load("pressure.png")
#bg = pygame.transform.scale(bg,(1600,800))
#thesize= bg.get_rect()
#bg= bg.convert()
number_of_particles=10
Spawn(number_of_particles)
clock = pygame.time.Clock()
while True:
ticks = clock.tick(60)
GetInput()
CollisionDetect()
Move()
Respawn(number_of_particles)
Draw()
if __name__ == '__main__': main()
该代码需要另一个脚本来计算流场本身。它还从文本文件中读取数据点以获得机翼的几何形状。 我没有提供这两个文件，但如果需要我可以添加它们。提前谢谢你。
我整理了您的代码并进行了一些更改，即为您的 类 添加范围并引入更多内容。没有 Flow 的进一步知识，我无法对此进行全面测试，但如果你能回复我，我可以做更多的事情。我在这里假设 'flow field' 可以通过 numpy.meshgrid 函数模拟。
import pygame,numpy,sys
import pygame.locals
import math
class Particle:
def __init__(self,xpos,ypos,vx,vy):
self.size = numpy.array([4,4])
self.radius =2
self.pos = numpy.array([xpos,ypos])
self.speed = numpy.array([30,0])
self.rectangle = numpy.hstack((self.pos,self.size))
def move(self):
self.pos += self.speed
self.rectangle = numpy.hstack((self.pos,self.size))
def distance(self,circle1):
return math.sqrt(numpy.sum((circle1.pos - self.pos)**2))
def collision(self,circle1):
result = False
if self.pos[1] >430 or self.pos[1]<160:
result = True
if self.distance(circle1) <= (circle1.radius+self.radius):
self.speed = circle1.speed
return result
class Particles:
def __init__(self,num_particles):
self.num = num_particles
self.particles =[]
self.spawn()
def spawn(self):
for i in range(self.num):
i=i*(300/self.num)
self.particles.append(Particle(0, 160+i,1,0))
def move(self):
for particle in self.particles:
particle.move()
if particle.pos[0] >1100:
self.particles.remove(particle)
if not self.particles: self.spawn()
def draw(self):
for particle in self.particles:
pygame.draw.rect(Surface,(150,0,0),particle.rectangle,0)
def collisiondetect(self,circle1):
for particle in self.particles:
if particle.collision(circle1):
self.particles.remove(particle)
def GetInput():
keystate = pygame.key.get_pressed()
for event in pygame.event.get():
if event.type == pygame.locals.QUIT or keystate[pygame.locals.K_ESCAPE]:
pygame.quit()
sys.exit()
#draw everything
def Draw(sw,cir):
Surface.fill((255,255,255))
cir.draw()
for i in range(panelcoordinates.shape[1]):
pygame.draw.line(Surface,(0,0,0),panelcoordinates[0,i-1],panelcoordinates[0,i],3)
sw.draw()
pygame.display.flip()
# a grid area
class Circle:
def __init__(self,xpos,ypos,vx,vy):
self.radius=16
self.pos = numpy.array([xpos,ypos])
self.speed = numpy.array([vx,vy])
class Circles:
def __init__(self,columns,rows):
self.list = []
grid_x,grid_y = numpy.meshgrid(numpy.linspace(0,1000,columns),numpy.linspace(200,400,rows))
grid_u,grid_v = numpy.meshgrid(numpy.linspace(20,20,columns),numpy.linspace(-1,1,rows))
for y in range(rows):
for x in range(columns):
c1= Circle(int(grid_x[y,x]),int(grid_y[y,x]),grid_u[y,x],grid_v[y,x])
self.list.append(c1)
def draw(self):
for circle in self.list:
pygame.draw.circle(Surface,(245,245,245),circle.pos,circle.radius)
def detectcollision(self,parts):
for circle in self.list:
parts.collisiondetect(circle)
if __name__ == '__main__':
#initialise variables
number_of_particles=10
Nx=30
Ny=10
#circles and particles
circles1 = Circles(Nx,Ny)
particles1 = Particles(number_of_particles)
#read the airfoil geometry
panel_x = numpy.array([400,425,450,500,600,500,450,425,400])
panel_y = numpy.array([300,325,330,320,300,280,270,275,300])
panelcoordinates= numpy.dstack((panel_x,panel_y))
#initialise PyGame
pygame.init()
clock = pygame.time.Clock()
Surface = pygame.display.set_mode((1000,600))
while True:
ticks = clock.tick(6)
GetInput()
circles1.detectcollision(particles1)
particles1.move()
Draw(particles1,circles1)
我也粗略地画了一个机翼，同样我不了解带有坐标的数据文件。
代码中的一个瓶颈可能是碰撞检测。 CollisionDetect() 计算每个粒子和每个圆之间的距离。然后，如果检测到碰撞，CircleCollide() 会在 Circles 中找到圆的索引（线性搜索），以便可以从 Velocities 中的相同索引中检索速度。显然，改进的时机已经成熟。
首先，Circleclass已经有了speedx/speedy属性中的速度，所以Velocities可以去掉。
其次，因为圆在固定的位置，你可以从粒子的位置计算出哪个圆最接近任何给定的粒子。
# You may already have these values from creating grid_x etc.
# if not, you only need to calculated them once, because the
# circles don't move
circle_spacing_x = Circles[1].x - Circles[0].x
circle_spacing_y = Circles[Nx].y - Circles[0].y
circle_first_x = Circles[0].x - circle_spacing_x / 2
circle_first_y = Circles[0].y - circle_spacing_y / 2
那么CollisionDetect()就变成了：
def CollisionDetect():
for particle in Particles:
if particle.y >430 or particle.y<160:
Particles.remove(particle)
continue
c = (particle.x - circle_first_x) // circle_spacing_x
r = (particle.y - circle_first_y) // circle_spacing_y
circle = Circles[r*Nx + c]
if ((circle.x - particle.x)**2 + (circle.y - particle.y)**2
<= (circle.radius+particle.radius)**2):
particle.speedx = int(circle.speedx)
particle.speedy = int(circle.speedy)