- Notifications
You must be signed in to change notification settings - Fork1
📈 poissonpy is a Python Poisson Equation library for scientific computing, image and video processing, and computer graphics.
License
bchao1/poissonpy
Folders and files
| Name | Name | Last commit message | Last commit date | |
|---|---|---|---|---|
Repository files navigation
Plug-and-play standalone library for solving 2D Poisson equations. Useful tool in scientific computing prototyping, image and video processing, computer graphics.
- Solves the Poisson equation on sqaure or non-square rectangular grids.
- Solves the Poisson equation on regions with arbitrary shape.
- Supports arbitrary boundary and interior conditions using
sympyfunction experssions ornumpyarrays. - Supports Dirichlet, Neumann, or mixed boundary conditions.
This package is only used to solve 2D Poisson equations. If you are looking for a general purpose and optimized PDE library, you might want to checkout theFEniCSx project.
Import necessary libraries.poissonpy utilizesnumpy andsympy greatly, so its best to import both:
importnumpyasnpfromsympyimportsin,cosfromsympy.abcimportx,yfrompoissonpyimportfunctional,utils,sovlers
In the following examples, we use a ground truth function to create a mock Poisson equation and compare the solver's solution with the analytical solution.
Define functions usingsympy function expressions ornumpy arrays:
f_expr=sin(x)+cos(y)# create sympy function expressionlaplacian_expr=functional.get_sp_laplacian_expr(f_expr)# create sympy laplacian function expressionf=functional.get_sp_function(f_expr)# create sympy functionlaplacian=functional.get_sp_function(laplacian_expr)# create sympy function
Define interior and Dirichlet boundary conditions:
interior=laplacianboundary= {"left": (f,"dirichlet"),"right": (f,"dirichlet"),"top": (f,"dirichlet"),"bottom": (f,"dirichlet")}
Initialize solver and solve Poisson equation:
solver=Poisson2DRectangle(((-2*np.pi,-2*np.pi), (2*np.pi,2*np.pi)),interior,boundary,X=100,Y=100)solution=solver.solve()
Plot solution and ground truth:
poissonpy.plot_3d(solver.x_grid,solver.y_grid,solution)poissonpy.plot_3d(solver.x_grid,solver.y_grid,f(solver.x_grid,solver.y_grid))
| Solution | Ground truth | Error |
|---|---|---|
 |  |  |
You can also define Neumann boundary conditions by specifyingneumann_x andneumann_y in the boundary condition parameter.
x_derivative_expr=functional.get_sp_derivative_expr(f_expr,x)y_derivative_expr=functional.get_sp_derivative_expr(f_expr,y)interior=laplacianboundary= {"left": (f,"dirichlet"),"right": (functional.get_sp_function(x_derivative_expr),"neumann_x"),"top": (f,"dirichlet"),"bottom": (functional.get_sp_function(y_derivative_expr),"neumann_y")}
| Solution | Ground truth | Error |
|---|---|---|
 |  |  |
If the boundary condition is purely Neumann, then the solution is not unique. Naively solving the Poisson equation gives bad results. In this case, you can set thezero_mean paramter toTrue, such that the solver finds a zero-mean solution.
solver=solvers.Poisson2DRectangle( ((-2*np.pi,-2*np.pi), (2*np.pi,2*np.pi)),interior,boundary,X=100,Y=100,zero_mean=True)
zero_mean=False | zero_mean=True | Ground truth |
|---|---|---|
 |  | |
It's also straightforward to define a Laplace equation -we simply set the interior laplacian value to 0. In the following example, we set the boundary values to be spatially-varying periodic functions.
interior=0# laplace equation formleft=poissonpy.get_2d_sympy_function(sin(y))right=poissonpy.get_2d_sympy_function(sin(y))top=poissonpy.get_2d_sympy_function(sin(x))bottom=poissonpy.get_2d_sympy_function(sin(x))boundary= {"left": (left,"dirichlet"),"right": (right,"dirichlet"),"top": (top,"dirichlet"),"bottom": (bottom,"dirichlet")}
Solve the Laplace equation:
solver=Poisson2DRectangle( ((-2*np.pi,-2*np.pi), (2*np.pi,2*np.pi)),interior,boundary,100,100)solution=solver.solve()poissonpy.plot_3d(solver.x_grid,solver.y_grid,solution,"solution")poissonpy.plot_2d(solution,"solution")
| 3D surface plot | 2D heatmap |
|---|---|
 |  |
Use thePoisson2DRegion class to solve the Poisson eqaution on a arbitrary-shaped function domain.poissonpy can be seamlessly integrated in gradient-domain image processing algorithms.
The following is an example wherepoissonpy is used to implement the image cloning algorithm proposed inPoisson Image Editing by Perez et al., 2003. Seeexamples/poisson_image_editing.py for more details.
# compute laplacian of interpolation functionGx_src,Gy_src=functional.get_np_gradient(source)Gx_target,Gy_target=functional.get_np_gradient(target)G_src_mag= (Gx_src**2+Gy_src**2)**0.5G_target_mag= (Gx_target**2+Gy_target**2)**0.5Gx=np.where(G_src_mag>G_target_mag,Gx_src,Gx_target)Gy=np.where(G_src_mag>G_target_mag,Gy_src,Gy_target)Gxx,_=functional.get_np_gradient(Gx,forward=False)_,Gyy=functional.get_np_gradient(Gy,forward=False)laplacian=Gxx+Gyy# solve interpolation functionsolver=solvers.Poisson2DRegion(mask,laplacian,target)solution=solver.solve()# alpha-blend interpolation and target functionblended=mask*solution+ (1-mask)*target
Another example of usingpoissonpy to implement flash artifacts and reflection removal, using the algorithm proposed inRemoving Photography Artifacts using Gradient Projection and Flash-Exposure Sampling by Agrawal et al. 2005. Seeexamples/flash_noflash.py for more details.
Gx_a,Gy_a=functional.get_np_gradient(ambient)Gx_f,Gy_f=functional.get_np_gradient(flash)# gradient projectiont= (Gx_a*Gx_f+Gy_a*Gy_f)/ (Gx_a**2+Gy_a**2+1e-8)Gx_f_proj=t*Gx_aGy_f_proj=t*Gy_a# compute laplacian (div of gradient)lap=functional.get_np_div(Gx_f_proj,Gy_f_proj)# integrate laplacian fieldsolver=solvers.Poisson2DRegion(mask,lap,flash)res=solver.solve()
