For n1n 1 lkn0nkn0n1 bkn0nsinpixn1 lknnnknnn1 bknnn0

This preview shows page 16 - 19 out of 19 pages.

for n=1:N-1 L(k(n,0,N),k(n,0,N))=1; B(k(n,0,N))=sin(pi*X(n+1)); L(k(n,N,N),k(n,N,N))=1; B(k(n,N,N))=0; end; Finally we define the parts of the matrix L and vector B that correspond to the equations for the interior points. for i=1:N-1 for j=1:N-1 L(k(i,j,N),k(i,j,N))=-4; L(k(i,j,N),k(i+1,j,N))=1; L(k(i,j,N),k(i-1,j,N))=1; L(k(i,j,N),k(i,j+1,N))=1; L(k(i,j,N),k(i,j-1,N))=1; end end Now we can solve the equation for F and plot the result. To do this we have to put the F values in a two dimensional grid FF and use the mesh command to do the 3-d plot. If X and Y are vectors of length n and Z is an nxn matrix then mesh(X,Y,Z) plots the points [X(j),Y(i),Z(i,j)] . F=L\B; FF=zeros(N+1,N+1); for i=0:N for j=0:N FF(j+1,i+1)=F(k(i,j,N)); end end mesh(X,X,FF); 9
We can print out the resulting graph using print laplace1.pdf (or print laplace1.jpg or print laplace1.eps ). This will produce a pdf file (or jpg or eps file) containing the graph. Here is the result: 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 Run the file laplaceeqn.m to produce this picture. Then modify the code to solve Laplace’s equation with boundary conditions: f ( x, 0) = sin( π x ) 0 x 1 f (0 , y ) = 0 0 y 1 f ( x, 1) = sin( π x ) 0 x 1 f (1 , y ) = 0 0 y 1 Say what code you modified, and hand in the resulting picture. Finally modify the code to solve Laplace’s equation with boundary conditions: f ( x, 0) = sin( π x ) 0 x 1 f (0 , y ) = 0 0 y 1 f y ( x, 1) = 0 0 x 1 f (1 , y ) = 0 0 y 1 The third boundary condition is called a Neumann boundary condition, and corresponds to detaching the rubber membrane from the wire along the top boundary. Again, say what code you modified, and hand in the resulting picture. 10
The first modification is to change the code defining the top and bottom boundary conditions to for n=1:N-1 L(k(n,0,N),k(n,0,N))=1; B(k(n,0,N))=sin(pi*X(n+1)); L(k(n,N,N),k(n,N,N))=1; B(k(n,N,N))=sin(pi*X(n+1)); end; The resulting plot looks like 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 The second modification is to change the code defining the top and bottom boundary conditions to for n=1:N-1 L(k(n,0,N),k(n,0,N))=1; B(k(n,0,N))=sin(pi*X(n+1)); L(k(n,N,N),k(n,N,N))=-1; L(k(n,N,N),k(n,N-1,N))=1; B(k(n,N,N))=0; end; The resulting plot looks like 11
0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 12

  • Left Quote Icon

    Student Picture

  • Left Quote Icon

    Student Picture

  • Left Quote Icon

    Student Picture