考察Vlasov-Maxwell系统,证明质量守恒、能量守恒,并探讨粒子方法求解该系统的弱解及能量保持性

考察以下系统:

\partial _{t}f+v\cdot \nabla _{x}f+E\cdot \nabla _{v}f=0, \quad \partial_{t}E=-\int_{\mathbb{R}^{d}}vf\textrm{d}v。 (1)

其中变量有t\in\mathbb{R}^{+}x\in\mathbb{R}^{d}v\in\mathbb{R}^{d}f是一个关于 t,xv 的未知函数,E只关于tx

a)证明质量和能量

m(t):=\int_{\mathbb{R}^{d}}\int_{\mathbb{R}^{d}}fdvdx

e(t):=\frac{1}{2}\int_{\mathbb{R}^{\mathrm{d}}}\int_{\mathbb{R}^{d}}\|v\|^{2}fdxdv+\frac{1}{2}\int_{\mathbb{R}^{d}}|E|^{2}dx

随时间保持不变。

b)现假设E(t,x)已知,则(1)简化为

\partial_{t}f+v \cdot \nabla_{x} f+E \cdot \nabla_{v} f=0。 (2)

我们将使用粒子方法来求解(2),初始条件为:

f(t=0, x, v)=f_{\text {in }}(x, v)

考虑由\left\{x_{p}(t), v_{p}(t)\right\}_{1 \leq p \leq P}表示P个粒子,其中x_{p}(t)表示第P个粒子的位置,v_{p}(t)表示其速度,每个粒子的权重为\frac{1}{P}x_{p}(t)v_{p}(t)的演化由以下方程决定:

\frac{d x_{p}}{d t}=v_{p}, \quad \frac{d v_{p}}{d t}=E(x_{p}) 。

证明

f^{P}(t, x, v):=\frac{1}{P} \sum_{p=1}^{P} \delta\left(x-x_{p}(t)\right) \delta\left(v-v_{p}(t)\right)

在分布意义下是(2)的弱解,其中(2)的初始条件是f_{0}^{P}=\frac{1}{P} \sum_{p=1}^{P} \delta\left(x-x_{p}(0)\right) \delta\left(v-v_{p}(0)\right)。这里\delta是Dirac-Delta函数。

c)回到原始系统(1)。考虑以下粒子方法

\begin{aligned} & x_{i}^{n+1}=x_{i}^{n}+\Delta t v_{i}^{n} & \\ & v_{p}^{n+1}=v_{p}^{n}-E\left(x_{p}^{n}\right) \Delta t & \\ & E_{i}^{n+1}=E_{i}^{n}-J_{i}^{n} \Delta t \end{aligned}

这里上标n表示时间步长,E\left(x_{p}^{n}\right)J_{i}^{n}分别定义如下

E\left(x_{p}^{n}\right):=\sum_{i} S\left(x_{i}-x_{p}^{n}\right) E_{i}^{n} \Delta x^{d}

J_{i}^{n}:=\frac{1}{P} \sum_{p=1}^{P}S\left(x_{i}-x_{p}^{n}\right) v_{p}^{n}

在这些表达式中,x_{i}表示具有网格大小\Delta x的均匀网格,S是一个满足\int_{\mathbb{R}^{d}} S(x) d x=1。这种方法是否保留能量?也就是说,我们是否有

\frac{1}{2 P} \sum_{p}\left(v_{p}^{n}\right)^{2}+\frac{1}{2} \sum_{i}\left(E_{i}^{n}\right)^{2} \Delta x^{d}=\frac{1}{2 P} \sum_{p}\left(v_{p}^{n+1}\right)^{2}+\frac{1}{2} \sum_{i}\left(E_{i}^{n+1}\right)^{2} \Delta x^{d}

如果是,证明它。如果不是,构造一个在离散级别上保持能量的离散格式。

证:

a) 我们先证明质量 m(t) 随时间保持不变。

质量 m(t) 定义为

m(t) = \int_{\mathbb{R}^d} \int_{\mathbb{R}^d} f \, dv \, dx

m(t) 关于时间 t 求导,并使用方程 (1) 中的第一个方程,我们有

\frac{dm}{dt} = \frac{d}{dt} \int_{\mathbb{R}^d} \int_{\mathbb{R}^d} f \, dv \, dx

通过分部积分,我们可以将时间导数移到积分内部:

\frac{dm}{dt} = \int_{\mathbb{R}^d} \int_{\mathbb{R}^d} \partial_t f \, dv \, dx

根据方程 (1) 中的第一个方程,\partial_t f + v \cdot \nabla_x f + E \cdot \nabla_v f = 0,因此

\frac{dm}{dt} = - \int_{\mathbb{R}^d} \int_{\mathbb{R}^d} (v \cdot \nabla_x f + E \cdot \nabla_v f) \, dv \, dx

使用散度定理,我们可以将 \nabla_x\nabla_v 的积分转换为边界积分。由于 f 在无穷远处趋于零,边界积分为零,因此

\frac{dm}{dt} = 0

这表明质量 m(t) 随时间保持不变。

接下来,我们证明能量 e(t) 随时间保持不变。

能量 e(t) 定义为

e(t) = \frac{1}{2} \int_{\mathbb{R}^d} \int_{\mathbb{R}^d} \|v\|^2 f \, dv \, dx + \frac{1}{2} \int_{\mathbb{R}^d} \|E\|^2 \, dx

e(t) 关于时间 t 求导,并使用方程 (1) 中的两个方程,我们有

\frac{de}{dt} = \frac{1}{2} \int_{\mathbb{R}^d} \int_{\mathbb{R}^d} 2v \cdot \partial_t v f \, dv \, dx + \frac{1}{2} \int_{\mathbb{R}^d} 2E \cdot \partial_t E \, dx

根据方程 (1) 中的第一个方程,\partial_t f = -v \cdot \nabla_x f - E \cdot \nabla_v f,以及第二个方程 \partial_t E = -\int_{\mathbb{R}^d} v f \, dv,我们可以将 \partial_t v\partial_t E 替换:

\frac{de}{dt} = \int_{\mathbb{R}^d} \int_{\mathbb{R}^d} (-v \cdot \nabla_x f - E \cdot \nabla_v f) v \, dv \, dx - \int_{\mathbb{R}^d} E \int_{\mathbb{R}^d} v f \, dv \, dx

使用散度定理和方程 (1) 中的第二个方程,我们可以证明

\frac{de}{dt} = 0.

这表明能量 e(t) 随时间保持不变。

b) 现在证明 f^P(t, x, v) 在分布意义下是方程 (2) 的弱解。

首先,我们计算 f^P(t, x, v) 的初始条件:

f_0^P = \frac{1}{P} \sum_{p=1}^P \delta(x - x_p(0)) \delta(v - v_p(0))

根据粒子演化方程,我们有

\frac{d x_p}{dt} = v_p, \quad \frac{d v_p}{dt} = E(x_p)

对于任意的测试函数 \phi(t, x, v),我们计算以下积分:

\int_{\mathbb{R}^d} \int_{\mathbb{R}^d} \left( \partial_t f^P + v \cdot \nabla_x f^P + E \cdot \nabla_v f^P \right) \phi \, dx \, dv

利用粒子表示,我们有

\begin{align*} &\int_{\mathbb{R}^d} \int_{\mathbb{R}^d} \left( \partial_t f^P + v \cdot \nabla_x f^P + E \cdot \nabla_v f^P \right) \phi \, dx \, dv \\ &= \int_{\mathbb{R}^d} \int_{\mathbb{R}^d} \frac{1}{P} \sum_{p=1}^P \left( \partial_t \delta(x - x_p(t)) \delta(v - v_p(t)) + v \cdot \nabla_x \delta(x - x_p(t)) \delta(v - v_p(t)) + E \cdot \nabla_v \delta(x - x_p(t)) \delta(v - v_p(t)) \right) \phi \, dx \, dv \end{align*}

使用 \delta 函数的性质,我们可以简化上述表达式为

\begin{align*} &\frac{1}{P} \sum_{p=1}^P \left( \partial_t \phi(t, x_p(t), v_p(t)) + v_p(t) \cdot \nabla_x \phi(t, x_p(t), v_p(t)) + E(x_p(t)) \cdot \nabla_v \phi(t, x_p(t), v_p(t)) \right) \\ &= \frac{1}{P} \sum_{p=1}^P \frac{d}{dt} \phi(t, x_p(t), v_p(t)) \end{align*}

根据粒子演化方程,我们有

\begin{align*} \frac{d}{dt} \phi(t, x_p(t), v_p(t)) &= \frac{\partial \phi}{\partial t}(t, x_p(t), v_p(t)) + \frac{d x_p}{dt} \cdot \nabla_x \phi(t, x_p(t), v_p(t)) + \frac{d v_p}{dt} \cdot \nabla_v \phi(t, x_p(t), v_p(t)) \\ &= \frac{\partial \phi}{\partial t}(t, x_p(t), v_p(t)) + v_p(t) \cdot \nabla_x \phi(t, x_p(t), v_p(t)) + E(x_p(t)) \cdot \nabla_v \phi(t, x_p(t), v_p(t)) \end{align*}

因此,我们有

\begin{align*} \frac{1}{P} \sum_{p=1}^P \frac{d}{dt} \phi(t, x_p(t), v_p(t)) &= \int_{\mathbb{R}^d} \int_{\mathbb{R}^d} \left( \partial_t f^P + v \cdot \nabla_x f^P + E \cdot \nabla_v f^P \right) \phi \, dx \, dv \\ &= \int_{\mathbb{R}^d} \int_{\mathbb{R}^d} 0 \cdot \phi \, dx \, dv \\ &= 0 \end{align*}

这表明 f^P(t, x, v) 在分布意义下是方程 (2) 的弱解。

c) 现在考虑原始系统 (1) 的粒子方法。

我们需要验证能量是否在离散时间步长上保持不变。考虑离散化的能量

E^n = \frac{1}{2P} \sum_{p=1}^P \left(v_p^n\right)^2 + \frac{1}{2} \sum_{i=1}^N \left(E_i^n\right)^2 \Delta x^d,

其中 N 是网格点的数量。

我们计算 E^{n+1} 并检查是否等于 E^n

根据粒子演化方程,我们有

v_p^{n+1} = v_p^n - E(x_p^n) \Delta t

因此,

\begin{align*} \frac{1}{2P} \sum_{p=1}^P \left(v_p^{n+1}\right)^2 &= \frac{1}{2P} \sum_{p=1}^P \left(v_p^n - E(x_p^n) \Delta t\right)^2 &= \frac{1}{2P} \sum_{p=1}^P \left(v_p^n\right)^2 - \frac{1}{P} \sum_{p=1}^P 2v_p^n E(x_p^n) \Delta t + \frac{1}{2P} \sum_{p=1}^P \left(E(x_p^n)\right)^2 (\Delta t)^2 \end{align*}

接下来,我们考虑电场的更新。假设我们使用显式欧拉方法来更新电场,我们有

E_i^{n+1} = E_i^n - \Delta t \int_{\mathbb{R}^d} v f^n(x_i, v) \, dv.

因此,电场能量的变化为

\begin{align*} \frac{1}{2} \sum_{i=1}^N \left(E_i^{n+1}\right)^2 \Delta x^d &= \frac{1}{2} \sum_{i=1}^N \left(E_i^n - \Delta t \int_{\mathbb{R}^d} v f^n(x_i, v) \, dv\right)^2 \Delta x^d \\ &= \frac{1}{2} \sum_{i=1}^N \left(E_i^n\right)^2 \Delta x^d - \Delta t \sum_{i=1}^N E_i^n \int_{\mathbb{R}^d} v f^n(x_i, v) \, dv \Delta x^d \\ &\quad + \frac{1}{2} \Delta t^2 \sum_{i=1}^N \left(\int_{\mathbb{R}^d} v f^n(x_i, v) \, dv\right)^2 \Delta x^d \end{align*}

将粒子和电场能量的变化结合起来,我们得到

\begin{align*} E^{n+1} &= \frac{1}{2P} \sum_{p=1}^P \left(v_p^n\right)^2 - \frac{1}{P} \sum_{p=1}^P 2v_p^n E(x_p^n) \Delta t + \frac{1}{2P} \sum_{p=1}^P \left(E(x_p^n)\right)^2 (\Delta t)^2 \\ &\quad + \frac{1}{2} \sum_{i=1}^N \left(E_i^n\right)^2 \Delta x^d - \Delta t \sum_{i=1}^N E_i^n \int_{\mathbb{R}^d} v f^n(x_i, v) \, dv \Delta x^d \\ &\quad + \frac{1}{2} \Delta t^2 \sum_{i=1}^N \left(\int_{\mathbb{R}^d} v f^n(x_i, v) \, dv\right)^2 \Delta x^d \end{align*}

为了保持能量,我们需要

\frac{1}{2P} \sum_{p=1}^P \left(v_p^n\right)^2 - \frac{1}{P} \sum_{p=1}^P 2v_p^n E(x_p^n) \Delta t + \frac{1}{2P} \sum_{p=1}^P \left(E(x_p^n)\right)^2 (\Delta t)^2 + \frac{1}{2} \sum_{i=1}^N \left(E_i^n\right)^2 \Delta x^d - \Delta t \sum_{i=1}^N E_i^n \int_{\mathbb{R}^d} v f^n(x_i, v) \, dv \Delta x^d + \frac{1}{2} \Delta t^2 \sum_{i=1}^N \left(\int_{\mathbb{R}^d} v f^n(x_i, v) \, dv\right)^2 \Delta x^d = E^n

我们可以看到,如果 \Delta t 足够小,那么 \left(E(x_p^n)\right)^2 (\Delta t)^2\left(\int_{\mathbb{R}^d} v f^n(x_i, v) \, dv\right)^2 \Delta x^d 这两项可以忽略不计,这样能量就会近似保持不变。然而,如果 \Delta t 不是足够小,那么这个方法可能不会保持能量。

为了构造一个在离散级别上保持能量的离散格式,我们可以使用半隐式或隐式时间积分方法。例如,使用隐式欧拉方法更新电场,我们可以得到一个保持能量的离散格式。以下是使用隐式欧拉方法更新电场的步骤:

  1. 粒子速度的更新仍然使用显式欧拉方法:

v_p^{n+1} = v_p^n + \Delta t E(x_p^n)

  1. 电场的更新使用隐式欧拉方法,我们需要解下面的方程来找到 E_i^{n+1}

E_i^{n+1} = E_i^n - \Delta t \int_{\mathbb{R}^d} v f^{n+1}(x_i, v) \, dv

由于 f^{n+1} 依赖于 E^{n+1},我们需要迭代求解 E^{n+1}f^{n+1}。这通常通过以下步骤完成:

  • 假设 E^{n+1} 的初始猜测,比如 E^{n+1, 0} = E^n

  • 使用这个电场值来更新粒子的速度和位置:

    v_p^{n+1, k} = v_p^n + \Delta t E(x_p^n)^k,

    x_p^{n+1, k} = x_p^n + \Delta t v_p^{n+1, k},

    其中 k 是迭代次数。

  • 使用更新的粒子速度来计算新的 f^{n+1, k}

    f^{n+1, k}(x_i, v) = \frac{1}{P} \sum_{p=1}^P \delta(x - x_p^{n+1, k}) \delta(v - v_p^{n+1, k})

  • 使用新的 f^{n+1, k} 来更新电场:

    E_i^{n+1, k+1} = E_i^n - \Delta t \int_{\mathbb{R}^d} v f^{n+1, k}(x_i, v) \, dv

  • 重复上述步骤直到 E^{n+1, k} 收敛。

在迭代收敛后,我们得到 E^{n+1}f^{n+1},然后可以计算新的能量 E^{n+1}

E^{n+1} = \frac{1}{2P} \sum_{p=1}^P \left(v_p^{n+1}\right)^2 + \frac{1}{2} \sum_{i=1}^N \left(E_i^{n+1}\right)^2 \Delta x^d

由于隐式方法的稳定性,这种方法通常能够更好地保持能量,即使在较大的时间步长下。这是因为隐式方法考虑了未来时间步的电场值,从而在数值解中引入了稳定性。

要证明这种方法保持能量,我们需要展示在迭代过程中,能量 E^{n+1} 保持不变。这通常涉及到证明迭代过程中的电场和粒子分布函数的更新不会增加总能量。由于这是一个迭代过程,我们可以通过数学归纳法来证明,假设在迭代 k 次后能量保持不变,那么在 k+1 次迭代后能量也将保持不变。这需要详细的数学分析和迭代收敛性的证明,这里不展开。

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 202,980评论 5 476
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 85,178评论 2 380
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 149,868评论 0 336
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 54,498评论 1 273
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 63,492评论 5 364
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,521评论 1 281
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 37,910评论 3 395
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,569评论 0 256
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,793评论 1 296
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,559评论 2 319
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,639评论 1 329
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,342评论 4 318
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,931评论 3 307
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,904评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,144评论 1 259
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 42,833评论 2 349
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,350评论 2 342

推荐阅读更多精彩内容