CODEWITHSUNDEEP
matlaboctavetrigonometry
Ziezi
Ziezi

Reputation: 6467

Why does plotting this equation not generate the correct curve?

The intention is to plot the following equation: P*sin(x)/x + cos(x), for P = 1.6 and x in [0, 5], which, ignoring the green-filled areas, should look someting like:

enter image description here

However, from the following code:

x = 0 : 0.01 : 5;     % ka/pi, where k-wavevector, a-lattice spacing.

P = 1.6;              % 2*m*U_0 / hbar^2.
rhs =  P * sinc(x*pi) + cos(x*pi);
rhs2 = P * ( sin(x*pi) / x*pi) + cos(x*pi);

plot(x, rhs, '--b', x, rhs2, 'b', x, -1*ones(size(x)), 'r', x, 1*ones(size(x)), 'r')
axis([0 5 -3 3])
xlabel('ka/pi')
legend('P*sinc(x) + cos(x)', '(2mU_0b)/(hbar^2) * sin(ka)/ka + cos(ka)', 'y = -1', 'y = 1')

what I currently get is:

enter image description here

What am I doing wrong here?

I am on Windows 10, Octave-4.2.1

Upvotes: 4

Views: 101

Answers (2)

Adriaan
Adriaan

Reputation: 18177

The MATLAB definition of the sinc is sinc(t) = sin(pi t)/(pi t), i.e. you must not multiply by pi in the rhs definition:

x = 0 : 0.01 : 5;     % ka/pi, where k-wavevector, a-lattice spacing.

P = 1.6;              % 2*m*U_0 / hbar^2.
rhs =  P * sinc(x)+ cos(x*pi);
rhs2 = P * (sin(x*pi) / x*pi) + cos(x*pi);

plot(x, rhs, 'b', x, rhs2, '--b', x, -1*ones(size(x)), 'r', x, 1*ones(size(x)), 'r')
axis([0 5 -3 3])
xlabel('ka/\pi')
legend('P*sinc(x) + cos(x)', '(2mU_0b)/(hbar^2) * sin(ka)/ka + cos(ka)', 'y = -1', 'y = 1')

enter image description here

Also note that for t=0 sinc(t)=1, whereas your rhs2 has sin(x pi)/(x pi), which for x=0 returns NaN, hence the difference in the two signals, as the second is a pure cosine.

I missed the element wise division and lack of brackets in the OP's sinc implementation, see am304's answer for that. Note that even when using element wise division and brackets you'll still miss the point at x=0, since that'll result to NaN.

Upvotes: 5

am304
am304

Reputation: 13876

There are two mistakes in your code:

  • the first one, as Adriaan pointed out, is incorrect use of the sinc function.
  • the second one, is you forgot to do the element-wise division when implementing your own version of sinc

Here's the corrected version - note the ./:

x = 0 : 0.01 : 5;     % ka/pi, where k-wavevector, a-lattice spacing.

P = 1.6;              % 2*m*U_0 / hbar^2.
rhs =  P * sinc(x)+ cos(x*pi);
rhs2 = P * (sin(x*pi) ./ (x*pi)) + cos(x*pi);

plot(x, rhs, 'b', x, rhs2, '--b', x, -1*ones(size(x)), 'r', x, 1*ones(size(x)), 'r')
axis([0 5 -3 3])
xlabel('ka/pi')
legend('P*sinc(x) + cos(x)', '(2mU_0b)/(hbar^2) * sin(ka)/ka + cos(ka)', 'y = -1', 'y = 1')

And the resulting plot is:

enter image description here

Upvotes: 5

Related Questions