121: 1D Poisson with point charge
Solve a Poisson equation
\[- \Delta u = 0\]
in $\Omega=(-1,1)$ with a point charge $Q$ at $x=0$.
module Example121_PoissonPointCharge1D
using Printf
using VoronoiFVM
using ExtendableGrids
function main(;nref=0,Plotter=nothing, verbose=false, unknown_storage=:sparse, brea=false)
# Create grid in (-1,1) refined around 0
hmax=0.2/2.0^nref
hmin=0.05/2.0^nref
X1=geomspace(-1.0,0.0, hmax,hmin)
X2=geomspace(0.0,1.0, hmin,hmax)
X=glue(X1,X2)
grid=VoronoiFVM.Grid(X)
# Edit default region numbers:
# additional boundary region 3 at 0.0
bfacemask!(grid, [0.0],[0.0],3)
# Material 1 left of 0
cellmask!(grid, [-1.0],[0.0],1)
# Material 2 right of 0
cellmask!(grid, [0.0],[1.0],2)
Q=0.0
function flux!(f,u0,edge)
u=unknowns(edge,u0)
f[1]=u[1,1]-u[1,2]
end
function storage!(f,u,node)
f[1]=u[1]
end
# Define boundary reaction defining charge
# Note that the term is written on the left hand side, therefore the - sign
function breaction!(f,u,node)
if node.region==3
f[1]=-Q
end
end
# Create physics
physics=VoronoiFVM.Physics(
flux=flux!,
storage=storage!,
breaction=breaction!
)
# Create system
sys=VoronoiFVM.System(grid,physics,unknown_storage=:dense)
# put potential into both regions
enable_species!(sys,1,[1,2])
# Set boundary conditions
boundary_dirichlet!(sys,1,1,1.0)
boundary_dirichlet!(sys,1,2,0.0)
# Create a solution array
inival=unknowns(sys)
U=unknowns(sys)
inival.=0
# Create solver control info
control=VoronoiFVM.NewtonControl()
control.verbose=verbose
if ispyplot(Plotter)
Plotter.clf()
end
# Solve and plot for several values of charge
for q in [0.1,0.2,0.4,0.8,1.6]
if brea
# Charge in reaction term
Q=q
else
# Charge as boundary condition
sys.boundary_values[1,3]=q
end
solve!(U,inival,sys, control=control)
# Plot data
if ispyplot(Plotter)
Plotter.grid()
coord=grid[Coordinates]
Plotter.plot(coord[1,:],U[1,:],label=@sprintf("Q=%.2f",q))
Plotter.xlabel("x")
Plotter.ylabel("\$\\phi\$")
Plotter.legend(loc="upper right")
Plotter.pause(1.0e-10)
end
end
return sum(U)
end
function test()
main()≈20.254591679579015
end
end
This page was generated using Literate.jl.