有关在Python中为Laplacian Solver正确实现Green函数的建议？

如何解决有关在Python中为Laplacian Solver正确实现Green函数的建议？

我发现this paper详细介绍了一种快速解决拉普拉斯/泊松问题（尤其是图像）的方法。作者介绍了一种非常简单但功能强大的算法，该算法比其他更复杂的算法可以更快地解决这些问题。本文的主要见解是易于计算的离散格林函数的表述。

尽管论文（第8页）中明确指出了该算法的伪代码，但我在计算Green函数时遇到了麻烦。我在下面复制了一个示例问题。

本文使用卷积关系A * B = Finv（F（A）⊙F（B））来找到内核，其中⊙是Hadamard积。应用这种关系，并注意到拉普拉斯内核K_L与格林函数V_mono的卷积等于狄拉克增量，格林函数V_mono可以从其他两个已知量计算得出。

但是在推导Green函数时，我总是遇到错误。例如，在下面的代码中，狄拉克增量的傅立叶变换中没有零项；但是拉普拉斯内核的傅立叶变换中只有一个零项。这使得这些量的Hadamard划分无法确定Green的函数。但是，将拉普拉斯算子应用于已知图像，并在相反的过程中使用它来查找格林函数时，结果几乎相等。 （唯一的功能差异是有问题的[0,0]条目在反向计算的条目中具有10 ^ 18条目，而经过分析计算的条目在其中具有无限值。）

由于这是本文的重点，因此我无法应用重构算法。我该如何解决算法问题？

#Import Test Image
import numpy as np
from scipy import misc
from scipy import fftpack
from scipy import signal
import matplotlib.pyplot as plt

img = misc.face(gray=True)#misc.face()[:,:,0]
pad=4
def pad_with(vector,pad_width,iaxis,kwargs):
    pad_value = kwargs.get('padder',10)
    vector[:pad_width[0]] = pad_value
    vector[-pad_width[1]:] = pad_value

img=np.pad(img,pad,pad_with,padder=0)#img is padded to allow for proper integration later


#Establish Dirac Kernel
#   Dirac Kernel Specified by Paper
kernelPaper=np.zeros(img.shape)
kernelPaper[1,1]=1
#   Dirac Kernel Specified by Link
kernelDirac=np.asarray([[0,0],[0,1,0]],dtype="float")
#kernel = np.outer(signal.gaussian(2,3),signal.gaussian(2,3))
sz = (img.shape[0] - kernelDirac.shape[0],img.shape[1] - kernelDirac.shape[1])  # total amount of padding
kernelDirac = np.pad(kernelDirac,(((sz[0]+1)//2,sz[0]//2),((sz[1]+1)//2,sz[1]//2)),'constant')
kernelDirac = fftpack.ifftshift(kernelDirac)

#Establish Laplace Kernel
#   Laplace Kernel Specified by Paper
kernelLaplacePaper=np.zeros(img.shape)
kernelLaplacePaper[np.tile(range(3),(3,1)),np.transpose(np.tile(range(3),1)))]=np.asarray([[0,-1,[-1,4,-1],dtype="float")

#   Laplace Kernel Specified by Link
kernelLaplace=np.asarray([[0,3))
sz = (img.shape[0] - kernelLaplace.shape[0],img.shape[1] - kernelLaplace.shape[1])  # total amount of padding
kernelLaplace = np.pad(kernelLaplace,'constant')
kernelLaplace = fftpack.ifftshift(kernelLaplace)


#Establish Test Case: Apply Laplacian to Known Image
filtered = np.real(fftpack.ifft2(fftpack.fft2(img) * fftpack.fft2(kernelLaplace)))
##plt.imshow(filtered,vmin=0,vmax=255)
##plt.show()


# ************Desired Result: Green's Function *************
#   Green's Function Specified by Paper
green_F_Paper = fftpack.fft2(kernelDirac)/fftpack.fft2(kernelLaplace) #This is actually F(Green's Function)
#   Green's Function Specified by Inverse Comparison
green_F_ComputeThroughReversal = fftpack.fft2(img)/fftpack.fft2(filtered)
print("Difference between paper Green's Function kernel and Green's Function computed using original image: ")
print(np.sum(abs(green_F_ComputeThroughReversal-green_F_Paper)))
print("Difference between paper Green's Function kernel and Green's Function computed using original image without [0,0] entry: ")
greenDiff=green_F_ComputeThroughReversal-green_F_Paper
greenDiff[0,0]=0
print(np.sum(abs(greenDiff)))

#   Laplace Kernel Specified by Inverse Comparison
kernelLaplace_ComputeThroughReversal = fftpack.ifft2(fftpack.fft2(kernelDirac)/green_F_ComputeThroughReversal)
print("Difference between paper Laplacian kernel and formulation found online: ")
print(np.sum(abs(kernelLaplace_ComputeThroughReversal-kernelLaplace)))



#Desired Result: Solving Laplacian
Ic=(fftpack.ifft2(fftpack.fft2(filtered) * fftpack.fft2(green_F_ComputeThroughReversal))).real #checking to see if reverse operation is preserved
Ic=Ic-Ic[0,0]
#Ic=Ic[pad:-pad,pad:-pad]
print("Difference between solved Laplacian and original image using kernels computed using original image: ")
print(np.sum(abs(Ic-img)))

Ic_Paper=(fftpack.ifft2(fftpack.fft2(filtered) * fftpack.fft2(green_F_Paper))).real #checking to see if reverse operation is preserved
Ic_Paper=Ic_Paper-Ic_Paper[0,pad:-pad]
print("Difference between reconstructed image and original image using paper formulation: ")
print(np.sum(abs(Ic_Paper-img)))

P.S。出于引用目的，该论文被Beaini等人命名为“用于使用Green函数卷积进行梯度域图像编辑的快速和最佳Laplacian解算器”。