矢量化/并行化GEKKO方程

Question

使用Phython的GEKKO库，并按照以下示例来求解波动方程，their original paper，需要创建很多方程。

m.Equations([u[i].dt() == v[i]for i in range(npx)])
m.Equations([v[i].dt() == c**2*(u[mod(i+1,npx)]-2.0*u[mod(i,npx)]+u[mod(i-1,npx)])/dx**2 for i in arange(0,npx,1)])

（（我已经通过mod将其简化为两行，并结合了边界条件，并让我从0到99，但是基本上与示例中的相同。）

它仍在使用for循环，即v[i]for i in range(npx)]。可以避免这种情况，从而像通过使用numpy数组一样对它们进行并行化/矢量化处理，从而使代码更具可读性或加快仿真速度吗？]

当然，对于此问题，我的笔记本电脑上的计算仅需<18s，但较大的问题可能会非常缓慢。

文档中有一些关于数组的内容，但是我不理解，也不知道这是否正是我想要的。

我使用上面的代码here附加了MWE。

Answer 1

如果要提高速度，请尝试使用m.options.IMODE=7进行顺序解决。对于仿真，它通常比IMODE=4（同时解决方案）快得多。

from gekko import GEKKO
import numpy as np
import time

start = time.time()
pi = np.pi

#Initialise model
m = GEKKO(remote = False)

#Discretisations(time and space)
npt = 200
npx = 100

m.time = np.linspace(0, 1, npt)
xpos = np.linspace(0, 2*pi, npx)
dx = xpos[1]-xpos[0]

#Define Variables
c = m.Const(value = 10)
u = [m.Var(value = np.cos(xpos[i])) for i in range(npx)]
v = [m.Var(value = np.sin(2*xpos[i])) for i in range(npx)]

#Automatic discretisation in time and manual discretisation in space
m.Equations([u[i].dt() == v[i] for i in range(npx)])
m.Equations([v[i].dt() == c**2*(u[np.mod(i+1,npx)]-2.0*u[np.mod(i,npx)]\
                                 +u[np.mod(i-1,npx)])/dx**2 \
                          for i in np.arange(0,npx,1)])

#Set options
m.options.imode = 7
m.options.solver = 1
m.options.nodes = 3
m.solve(disp = False, GUI = False)
print('Solve Time: ' + str(m.options.SOLVETIME))
print('Total Time: ' + str(time.time()-start))

使用IMODE=4，花费更多的时间来解决问题，而花费更少的时间来解决它。

Solve Time: 2.7666
Total Time: 25.547

使用IMODE=7，花费更少的时间来解决问题，但是花费更多的时间来解决它。

Solve Time: 8.2271
Total Time: 11.819

当您遇到更大的问题时，IMODE=7将比IMODE=4花费更少的时间进行仿真。

对于并行选项，您可以并行构建模型，但是该模型的解决方案除了IMODE=7之外，没有很多其他选项可以实现更快的解决方案。您可以通过以下方式创建并行的Gekko应用程序：

• 在IPOPT中将并行线性求解器与ma77，ma97等一起使用。与我对大规模问题所做的某些测试相比，这通常仅会使速度提高20-60％。这些选项在公开发布的IPOPT版本中不可用，因为求解器需要许可证。线性求解器MUMPS与Gekko一起分发，但不包括并行支持（尽管稍后可能会提供）。问题在于，求解器只是解决方案的一部分，即使求解器无限快速，自动微分，目标求值和方程式求值仍然需要大约50％的CPU时间。

• 独立运行。通常将其称为“大规模并行”，因为可以将进程拆分为单独的线程，然后在所有子进程完成后再次合并代码。此instructional material on multi-threading和How to do parallel Python Gekko?显示了如何并行计算。您的问题未针对多线程设置，因为方程式一起求解。