我正在考虑使用petsc4py代替scipy.integrate.odeint（它是 Fortran 求解器的包装器）来解决涉及 ODE 系统求解的问题。这个问题有可能变得僵硬。写下它的雅可比是非常困难的。

到目前为止，我已经能够通过用“类似 C”（使用numba或Cython）编写 RHS 函数来产生合理的速度增益。我想获得更多的性能，因此我考虑了 PETSc。

介绍问题的玩具版本：变形虫

考虑下图，它定义了变形虫的基础：

所以，变形虫由 $N$ 节点和边，受限于二维空间。描述变形虫的节点的位置随时间而变化。

让位置 $n$ th ( $n \in \{0, 1, ..., N-1\}$ ) 节点被描述为一个位置向量， $\textbf{x}_n$ ：

x_{n} = [\begin{array}{cc} x_{0, n} \\ x_{1, n} \end{array}]

$\textbf{x}_n = \left[\begin{array}{cc} x_{0,n} \\ x_{1,n} \end{array}\right]$

为简化起见，假设变形虫生活在一个具有一阶运动方程的世界中，因此给定一个作用力 $\textbf{F}$ 在节点上，节点速度 $\dot{\textbf{x}}$ 被生产：

F = c \dot{x}

$\textbf{F} = c \dot{\textbf{x}}$

在哪里 $c$ 是我选择的常数。

我们对跟踪每个节点的位置非常感兴趣，所以我们立即生成 $2n$ 感兴趣的 ODE（每个空间分量一个 ODE）（假设边是弹性弹簧）：

\frac{d p_{n}}{d t} = elastic edge forces + F_{n, a c t i v e}

$\frac{\textrm{d}\textbf{p}_n}{dt} = \text{elastic edge forces} + \textbf{F}_{n, active}$

$\textbf{F}_{n, active}$ 是变形虫在节点处产生的力 $n$ . 可以想象可以设计出一系列规则来确定如何 $\textbf{F}_{n, active}$ 是确定的，我建议我们考虑以下规则：

1）变形虫里面有红球和蓝球。^ 2) 红球和蓝球都有“开启”（象征 $r^*$ / $b^*$ ）和“关闭”（象征 $r$ / $b$ ）状态。

3）红球和蓝球都可以在变形虫的身体中找到，它们总是被关闭，或者在变形虫的节点处，它们可以被打开或关闭。

4）无论状态位置如何，红球的总数都是一个常数 $R$ 同样对于蓝球，一个常数 $B$

5）当红球或蓝球在体内时，它们可以自由地跳到他们想要的任何节点上，但是当它们在节点上时，由于某种“扩散”，它们中的一些可能会被诱导移动到相邻节点

6）在节点处打开红球会导致试图减少变形虫面积的力

7）在节点处打开蓝色球会导致试图增加变形虫面积的力

8）打开的红球往往会关闭打开的蓝球，反之亦然，打开的蓝球往往会关闭打开的红球

9) ...

好吧，你明白我的意思了。我可以提出一个模型，该模型在各种感兴趣的变量如何相互作用方面相当复杂。

由于涉及大量方程，考虑写下雅可比已经很乏味了。但是，由于我们可以使用计算机，因此我们可能可以编写管理特定交互的函数，这些函数没有简洁的分析形式（更不用说它们的导数是否具有简洁的分析形式），所以我们可能会有一堆乱七八糟的东西。如果我们要继续尝试产生雅可比行列式，则需要分段函数来逼近它们......

我看到的所有关于 PETSc 时间步进问题的玩具示例都定义了雅可比矩阵，所以我想知道我是否会从切换到它中获得速度增益，如果我的计算成本高的原因之一可能是因为没有被能够提供雅可比函数吗？