Python：有什么方法可以在没有填充数据的情况下获得周期性3D信号/场的n阶抗素？

发布于 2025-01-23 06:25:40 字数 578 浏览 0 评论 0 原文

如标题中所述，我想获得3D场（例如带有形状的数组（1024,1024,1024））的n-th-th（例如4-th）订单抗动力，每一侧l。 如果我需要超过周期L的抗体值，则必须在求解抗体之前定期添加数据（“包装”）。 ，这种方法是高度的存储器消耗，特别是对于1024^3阵列。我的Python代码片段如下所示：

import numpy as np
from scipy.interpolate import make_interp_spline
data = np.load("xxx.npy") # shape: (1024,1024,1024)
data_pad = np.pad(data, ( (1024,1024), (1024,1024), (1024,1024) ), "wrap")

# 4th antiderivative of the padded data along x
interp_x = make_interp_spline(x_pad, data_pad, k=3, axis=0)
interp_x4 = interp_x.antiderivative(4)

是否有其他方法可以实现相同的目标而无需填充数据？

原文

As stated in the title, I want to obtain the n-th (e.g. 4-th) order antiderivative of a 3D field (e.g. array with shape (1024,1024,1024) ) with period L on each side. If I need the antiderivative values beyond period L, I have to pad the data periodically ("wrap") before solving the antiderivative. However, this approach is highly memory consumption, particularly for the 1024^3 array. My Python code snippet is shown below:

import numpy as np
from scipy.interpolate import make_interp_spline
data = np.load("xxx.npy") # shape: (1024,1024,1024)
data_pad = np.pad(data, ( (1024,1024), (1024,1024), (1024,1024) ), "wrap")

# 4th antiderivative of the padded data along x
interp_x = make_interp_spline(x_pad, data_pad, k=3, axis=0)
interp_x4 = interp_x.antiderivative(4)

Is there other way to achieve the same goal without padding data?

分享到QQ

分享到微博

如果你对这篇内容有疑问，欢迎到本站社区发帖提问参与讨论，获取更多帮助，或者扫码二维码加入 Web 技术交流群。

发布评论

需要登录才能够评论，你可以免费注册一个本站的账号。

裂开嘴轻声笑有多痛 2025-01-30 06:25:40

正如我在定期信号的评论中提到的那样，您可以从傅立叶转变中获得抗动力。

如果您的函数为（n+1） times可区分，则可以获取 n 衍生品，此后，您将遇到与 gibbs netomenon

可悲的是，我们在这里没有方程式支持

，但使用，AS d（Omega） * F（Omega），如果 d（Omega）到处都是零，则可以颠倒派生。

让我们从一个1D示例开始，

import numpy as np

y = np.linspace(-4, 4, 256)
Y = np.fft.fft(y)
dt = np.exp(-np.fft.fftfreq(len(Y)) * 2j * np.pi)
dyOp = (dt**-0.5 - dt**0.5) * len(y) / (4 - (-4))
y_derivative = np.fft.ifft(dyOp * Y).real 
y_a1 = np.fft.ifft(np.divide(Y, dyOp, where=dyOp!=0)).real
y_a2 = np.fft.ifft(np.divide(Y, dyOp**2, where=dyOp!=0)).real
y_a3 = np.fft.ifft(np.divide(Y, dyOp**3, where=dyOp!=0)).real

plt.plot(y, y)
plt.plot(y, y_a1)
plt.plot(y, y_a2)
plt.plot(y, y_a3)
plt.ylim([-8, 8])
plt.legend(['y=t/4', '$d^{-1}ydt
如果您的数据具有多个维度，则可以为每个轴个人应用此示例
，假设您的尺寸为 0＆lt; = x，y，z＆lt; 1 在周期1的情况
下，您可以写作
period_x = 1
period_y = 1
period_z = 1

# used 1j pi, to avoid taking the square root down there
dx = np.exp(- 1j * np.pi * np.fft.fftfreq(data.shape[0])).reshape(-1, 1, 1)
dy = np.exp(- 1j * np.pi * np.fft.fftfreq(data.shape[1])).reshape(1, -1, 1)
dz = np.exp(- 1j * np.pi * np.fft.fftfreq(data.shape[2])).reshape(1, 1, -1)

data_dx = (dx.conj() - dx) * data.shape[0] / period_x
data_dy = (dx.conj() - dx) * data.shape[1] / period_y
data_dz = (dx.conj() - dx) * data.shape[2] / period_z

# maybe you will want to set some tolerance arround zero
data_dx[data_dx == 0] = 1
data_dy[data_dy == 0] = 1
data_dz[data_dz == 0] = 1

假设要求解毒药方程，您可以通过
D = np.fft.nfftn(data)
divergent_op = (data_dx**2 + data_dy**2 + data_dz**2)
data_divergent = np.fft.ifftn(D / divergent_op).real

, '$d^{-2}y dt^2
如果您的数据具有多个维度，则可以为每个轴个人应用此示例
，假设您的尺寸为 0＆lt; = x，y，z＆lt; 1 在周期1的情况
下，您可以写作

假设要求解毒药方程，您可以通过

, '$d^{-3}y dt^3
如果您的数据具有多个维度，则可以为每个轴个人应用此示例
，假设您的尺寸为 0＆lt; = x，y，z＆lt; 1 在周期1的情况
下，您可以写作

假设要求解毒药方程，您可以通过

], ncol=2)

如果您的数据具有多个维度，则可以为每个轴个人应用此示例

，假设您的尺寸为 0＆lt; = x，y，z＆lt; 1 在周期1的情况

下，您可以写作

假设要求解毒药方程，您可以通过

As I mentioned in the comments for periodic signals you can get antiderivatives from the Fourier transfrom.

If your function is (n+1) times differentiable, you can get n derivatives, after that you will have problems related to the Gibbs phenomenon

It is sad that we don't have equations support here

but using the translation identities you can get e.g. (f(x+a) - f(x-a)) / (2*a), as D(omega) * F(omega), and if D(omega) is non-zero everywhere, you can invert the derivation.

Let's start with a 1D example

import numpy as np

y = np.linspace(-4, 4, 256)
Y = np.fft.fft(y)
dt = np.exp(-np.fft.fftfreq(len(Y)) * 2j * np.pi)
dyOp = (dt**-0.5 - dt**0.5) * len(y) / (4 - (-4))
y_derivative = np.fft.ifft(dyOp * Y).real 
y_a1 = np.fft.ifft(np.divide(Y, dyOp, where=dyOp!=0)).real
y_a2 = np.fft.ifft(np.divide(Y, dyOp**2, where=dyOp!=0)).real
y_a3 = np.fft.ifft(np.divide(Y, dyOp**3, where=dyOp!=0)).real

plt.plot(y, y)
plt.plot(y, y_a1)
plt.plot(y, y_a2)
plt.plot(y, y_a3)
plt.ylim([-8, 8])
plt.legend(['y=t/4', '$d^{-1}ydt
If your data has multiple dimensions you can apply this for for each axis individually
Assuming that your dimensions are 0 <= x,y,z < 1 with period 1
then you could write
period_x = 1
period_y = 1
period_z = 1

# used 1j pi, to avoid taking the square root down there
dx = np.exp(- 1j * np.pi * np.fft.fftfreq(data.shape[0])).reshape(-1, 1, 1)
dy = np.exp(- 1j * np.pi * np.fft.fftfreq(data.shape[1])).reshape(1, -1, 1)
dz = np.exp(- 1j * np.pi * np.fft.fftfreq(data.shape[2])).reshape(1, 1, -1)

data_dx = (dx.conj() - dx) * data.shape[0] / period_x
data_dy = (dx.conj() - dx) * data.shape[1] / period_y
data_dz = (dx.conj() - dx) * data.shape[2] / period_z

# maybe you will want to set some tolerance arround zero
data_dx[data_dx == 0] = 1
data_dy[data_dy == 0] = 1
data_dz[data_dz == 0] = 1

Suppose you want to solve the Poison's equation, you can do it by
D = np.fft.nfftn(data)
divergent_op = (data_dx**2 + data_dy**2 + data_dz**2)
data_divergent = np.fft.ifftn(D / divergent_op).real

, '$d^{-2}y dt^2
If your data has multiple dimensions you can apply this for for each axis individually
Assuming that your dimensions are 0 <= x,y,z < 1 with period 1
then you could write

Suppose you want to solve the Poison's equation, you can do it by

, '$d^{-3}y dt^3
If your data has multiple dimensions you can apply this for for each axis individually
Assuming that your dimensions are 0 <= x,y,z < 1 with period 1
then you could write

Suppose you want to solve the Poison's equation, you can do it by

], ncol=2)