4

I'm going to start off by stating that, yes, this is homework (my first homework question on stackoverflow!). But I don't want you to solve it for me, I just want some guidance!

The equation in question is this: Series in question

I'm told to take N = 50, phi1 = 300, phi2 = 400, 0<=x<=1, and 0<=y<=1, and to let x and y be vectors of 100 equally spaced points, including the end points.

So the first thing I did was set those variables, and used x = linspace(0,1) and y = linspace(0,1) to make the correct vectors.

The question is Write a MATLAB script file called potential.m which calculates phi(x,y) and makes a filled contour plot versus x and y using the built-in function contourf (see the help command in MATLAB for examples). Make sure the figure is labeled properly. (Hint: the top and bottom portions of your domain should be hotter at about 400 degrees versus the left and right sides which should be at 300 degrees).

However, previously, I've calculated phi using either x or y as a constant. How am I supposed to calculate it where both are variables? Do I hold x steady, while running through every number in the vector of y, assigning that to a matrix, incrementing x to the next number in its vector after running through every value of y again and again? And then doing the same process, but slowly incrementing y instead?

If so, I've been using a loop that increments to the next row every time it loops through all 100 values. If I did it that way, I would end up with a massive matrix that has 200 rows and 100 columns. How would I use that in the linspace function?

If that's correct, this is how I'm finding my matrix:

clear
clc
format compact
x = linspace(0,1);
y = linspace(0,1);
N = 50;
phi1 = 300;
phi2 = 400;
phi = 0;
sum = 0;
for j = 1:100
    for i = 1:100
        for n = 1:N
            sum = sum + ((2/(n*pi))*(((phi2-phi1)*(cos(n*pi)-1))/((exp(n*pi))-(exp(-n*pi))))*((1-(exp(-n*pi)))*(exp(n*pi*y(i)))+((exp(n*pi))-1)*(exp(-n*pi*y(i))))*sin(n*pi*x(j)));
        end
        phi(j,i) = phi1 - sum;
    end
end
for j = 1:100
    for i = 1:100
        for n = 1:N
            sum = sum + ((2/(n*pi))*(((phi2-phi1)*(cos(n*pi)-1))/((exp(n*pi))-(exp(-n*pi))))*((1-(exp(-n*pi)))*(exp(n*pi*y(j)))+((exp(n*pi))-1)*(exp(-n*pi*y(j))))*sin(n*pi*x(i)));
        end
        phi(j+100,i) = phi1 - sum;
    end
end

This is the definition of contourf. I think I have to use contourf(X,Y,Z):

contourf(X,Y,Z), contourf(X,Y,Z,n), and contourf(X,Y,Z,v) draw filled contour plots of Z using X and Y to determine the x- and y-axis limits. When X and Y are matrices, they must be the same size as Z and must be monotonically increasing.

Here is the new code:

N = 50;
phi1 = 300;
phi2 = 400;
[x, y, n] = meshgrid(linspace(0,1),linspace(0,1),1:N)
f = phi1-((2./(n.*pi)).*(((phi2-phi1).*(cos(n.*pi)-1))./((exp(n.*pi))-(exp(-n.*pi)))).*((1-(exp(-1.*n.*pi))).*(exp(n.*pi.*y))+((exp(n.*pi))-1).*(exp(-1.*n.*pi.*y))).*sin(n.*pi.*x));
g = sum(f,3);
[x1,y1] = meshgrid(linspace(0,1),linspace(0,1));
contourf(x1,y1,g)
TheTreeMan
  • 913
  • 8
  • 23
  • 38

2 Answers2

3

The reason your code takes so long to calculate the phi matrix is that you didn't pre-allocate the array. The error about size happens because phi is not 100x100. But instead of fixing those things, there's an even better way...

MATLAB is a MATrix LABoratory so this type of equation is pretty easy to compute using matrix operations. Hints:

  1. Instead of looping over the values, rows, or columns of x and y, construct matrices to represent all the possible input combinations. Check out meshgrid for this.

  2. You're still going to need a loop to sum over n = 1:N. But for each value of n, you can evaluate your equation for all x's and y's at once (using the matrices from hint 1). The key to making this work is using element-by-element operators, such as .* and ./.

Using matrix operations like this is The Matlab Way. Learn it and love it. (And get frustrated when using most other languages that don't have them.)

Good luck with your homework!

shoelzer
  • 10,648
  • 2
  • 28
  • 49
  • This is a stupid question, but if I pre-allocate the tray like that, will my matrix stay that huge, even after inserting values? If I try to do the controurf, will it also run through all of those zeroes as well? I'm impressed that it cut down the time to calculate that array by 2/3 though! That's crazy. I'm having a hard time understanding how to use the meshgrid function for this purpose. Do you think you can explain it really fast? It's difficult to read and understand what it's doing, since it's a little beyond my programming level. – TheTreeMan Jan 25 '13 at 04:41
  • 1
    a 200 by 100 element matrix is not huge. it only takes ~160 KB memory... vectorize your code following the advice from @Shoelzer – bla Jan 25 '13 at 04:53
  • 1
    @TheTreeMan natan's answer is a good example of using meshgrid. – shoelzer Jan 25 '13 at 05:14
3

Vectorize the code. For example you can write f(x,y,n) with:

 [x y n] = meshgrid(-1:0.1:1,-1:0.1:1,1:10);
 f=exp(x.^2-y.^2).*n ;

f is a 3D matrix now just sum over the right dimension...

 g=sum(f,3);

in order to use contourf, we'll take only the 2D part of x,y:

 [x1 y1] = meshgrid(-1:0.1:1,-1:0.1:1);    
 contourf(x1,y1,g)

enter image description here

bla
  • 25,846
  • 10
  • 70
  • 101
  • So I understand the meshgrid, and how it makes every combination of the two matrices, and fast. But I'm not sure how to put my sum into the format to use it like this. I have the sum = sum + blah loop, and then have to do the phi1 - sum. The ways you guys are showing me seem so much easier than how I'm doing it, but I'm not sure how to adapt my problem to your method. – TheTreeMan Jan 25 '13 at 05:08
  • So wait, do I do the loop as said, but I'm using the x and y as a result of meshgrid, and using the dot operators to make the operations performed in the equation happen for every single part of the matrix at once? – TheTreeMan Jan 25 '13 at 05:11
  • I added at the end of the original post what I tried out, and it didn't work. Do you think you can tell me where I went wrong? – TheTreeMan Jan 25 '13 at 05:17
  • 1
    see my edited answer, also what didn't work? by the way, don't use the word `sum` for your variable, you are overwriting an important function of matlab... – bla Jan 25 '13 at 05:24
  • I tried to make a program in the format of your previous answer, and I put it in the original post above, at the bottom. It gave me errors because I went about it wrong. So the 3D matrix is x and y values in the first and second dimension, and the actual values I want in the third dimension. So I do phi = sum(f,3)? I looked up how to sum 3D matrices in one dimension, and I found that the sum function can be used as sum(matrix,dimension). Unless I'm misunderstanding everything, which is most probable. Thank you for your help by the way! Things are becoming clear. – TheTreeMan Jan 25 '13 at 05:26
  • Oh my god, IT WORKED. IT ACTUALLY WORKED! And it only used EIGHT LINES. What the hell, man. What the hell. Seriously. – TheTreeMan Jan 25 '13 at 06:00
  • Thank you so much for sticking with me, even when I wasn't getting it! You're the best. – TheTreeMan Jan 25 '13 at 06:12