【问题标题】:Compute divergence with python用python计算分歧
【发布时间】:2023-04-05 19:45:01
【问题描述】:

根据this 的回答,数值向量场的散度可以这样计算:

def divergence(f):
    num_dims = len(f)
    return np.ufunc.reduce(np.add, [np.gradient(f[i], axis=i) for i in range(num_dims)])

但是,我注意到输出似乎很大程度上取决于网格分辨率,所以似乎有问题!

如果我看一个例子:

我们有以下向量场F:

F(x) = cos(x+2y)
F(y) = sin(x-2y)

如果我们计算散度(使用 Mathematica):

Div[{Cos[x + 2*y], Sin[x - 2*y]}, {x, y}]

我们得到:

-2 Cos[x - 2 y] - Sin[x + 2 y]

在y [-2,2]和x [-2,2]的范围内有最大值:

N[Max[Table[-2 Cos[x - 2 y] - Sin[x + 2 y], {x, -2, 2 }, {y, -2, 2}]]] = 2.938

使用此处给出的散度方程,我们得到以下图表,最大值与分辨率的关系(NxN:x 和 y 方向的值数)。这些都没有接近 3。

代码如下:

import numpy as np
import matplotlib.pyplot as plt

# Boundaries
ymin = -2.; ymax = 2.
xmin = -2.; xmax = 2.
# Number of points (NxN)
N = 20

# Divergence function
def divergence(f):
    num_dims = len(f)
    return np.ufunc.reduce(np.add, [np.gradient(f[i], axis=i) for i in range(num_dims)])

# Create Meshgrid
x = np.linspace(xmin,xmax, N)
y = np.linspace(ymin,ymax, N)
xx, yy = np.meshgrid(x, y)


# Define 2D Vector Field
Fx  = np.cos(xx + 2*yy)
Fy  = np.sin(xx - 2*yy)

F = np.array([Fx, Fy])
# Compute Divergence
g = divergence(F)

print("Max: ", np.max(g.flatten()))

plt.imshow(g)
plt.colorbar()

编辑:
创建情节:

# %%
a = []
for N in range(20,100):
    # Number of points (NxN)
    # = 20
    # Boundaries
    ymin = -2.; ymax = 2.
    xmin = -2.; xmax = 2.
    
    
    # Deivergence function
    def divergence(f):
        num_dims = len(f)
        return np.ufunc.reduce(np.add, [np.gradient(f[i], axis=i) for i in range(num_dims)])
    
  
    
    # Create Meshgrid
    x = np.linspace(xmin,xmax, N)
    y = np.linspace(ymin,ymax, N)
    xx, yy = np.meshgrid(x, y)
    
    
    # Define 2D Vector Field
    Fx  = np.cos(xx + 2*yy)
    Fy  = np.sin(xx - 2*yy)
    
    F = np.array([Fx, Fy])
    # Compute Divergence
    g = divergence(F)
    
    print("Max: ", np.max(g.flatten()))
    a.append(np.max(g.flatten()))
plt.plot(a)

【问题讨论】:

    标签:
    python
    matplotlib
    vector