ITADN

Issue solving 1D Poisson Eq.

#337Closedmarkowkes 创建于 2023-05-01
M
markowkescommented
I'm trying to use IterativeSolvers.jl but am struggling to get a very simple example working. This code solves L*P=R where L is a Laplacian operator (discretized with a finite difference) and R is zero. The first and last equations set the boundary condition of P=10 and P=100, respectively. I'm solving the system with conjugate gradient, BiCGStab, IDRS, and the "\\". All the methods are not working except the base Julia operator "\\". Any information on why this is not working would be very helpful. ``` using IterativeSolvers using Plots # Create grid Nx=100 x=range(0,Nx,Nx) dx=x[2]-x[1] # Create Laplacian - eventually I want to do this with LinearMap.jl L = zeros(Nx,Nx) L[ 1, 1] = 1.0 # Left BC: Dirchlet L[Nx,Nx] = 1.0 # Right BC: Dirchlet # Interior points: Laplacian for i=2:Nx-1 L[i,i ] = -2/dx^2 L[i,i-1] = 1/dx^2 L[i,i+1] = 1/dx^2 end # Define RHS R=zeros(Nx) R[ 1] = 10 # Left BC R[Nx] = 100 # Right BC # Solve for p using various methods P1 = cg(L, R; abstol=1e-6,verbose=true) P2 = bicgstabl(L, R; abstol=1e-6,verbose=true) P3 = idrs(L, R, abstol=1e-6,verbose=true) P4 = L\R # Plot results plt1 = plot(x,P1,title="Conjugate Gradient",legend=false) plt2 = plot(x,P2,title="BiCGStab" ,legend=false) plt3 = plot(x,P3,title="IDRS" ,legend=false) plt4 = plot(x,P4,title="Julia" ,legend=false) plot(plt1,plt2,plt3,plt4,layout = 4) ``` ![result](https://user-images.githubusercontent.com/17165240/235509784-a82353a7-f31b-4307-aab6-2743cad9d3e9.svg)
关闭于 2023-05-02 1 条评论