index out of bounds because numel(y)=0 error

I'm trying to create a code which uses the trapezoid rule to find the integral of a function over a given bound. Here's my code so far...

H=.5.^[2:1:13];
nH=length(H);
a=1;
b=2;
n=.5.^[2:1:13]
for k=1:nH
    h=(b-a)/(n(k)-1);
    x=[1:h:nH];
    y=log(1+x);
    I=(y(1) + y(nH) + 2.*sum(y(2:end-1))).*h.*.5 
end

It seems to work perfectly if I set n=[2:13] however, if I set it to what it is now ( .5.^[2:13] ) i get the error:

??? Attempted to access y(1); index out of bounds because
numel(y)=0.

Error in ==> trapezoid_rule at 10
    I=(y(1) + y(nH) + 2.*sum(y(2:end-1))).*h.*.5

Any help would be appreciated. Thanks!

[Matt Ogden]

Ok I'm understanding the responses. I think I messed up writing the code for h=

Since h is the step size for x, I changed x=[1,H,2] which should create a certain number of points depending on .5^[2:13].

H=.5.^[2:13];
nH=length(H);
a=1;
b=2;
for k=1:nH
  x=[1:H:2];
  y=log(1+x);
  I=(y(1) + y(nH) + 2.*sum(y(2:end-1))).*h.*.5 
end

My only error now is...

Attempted to access y(12); index out of bounds because
numel(y)=5.

This means that when k=5, H=.03125. So there will be 64 points between 1 and 2 correct? I don't understand why it isn't able to produce a y value for every 64 points. Any clarification would be appreciated. Thanks!

[Walter Roberson]

You cannot use a vector as the "step size" argument for the colon operator.

Perhaps you wanted 1:H(k):2 ?

And perhaps you want to assign to I(k) instead of to I ?

[Matt Ogden]

Those are both good ideas and I will implement them; however, I still receive the same error.

H=.5.^[2:13];
nH=length(H);
a=1;
b=2;
for k=1:nH
    x=1:H(k):2;
    y=log(1+x);
    I(k)=(y(1) + y(nH) + 2.*sum(y(2:end-1))).*h.*.5 
end

I would like my final output to be all of the I(k)'s calculated at various step sizes. The step size changing from .5^[2:13]. Why am I receiving an error once k=5? Should there just be 64 x-values?

Answer 1

[Walter Roberson]

When you have 0.5.^[2:1:13], the first element will be 0.5^2 which is 0.25 .

You start with a=1 and b=2. Then in the expression (b-a)/(n(k)-1) that would be (2-1)/(0.25-1) which is going to be negative. 1 : a negative value : a positive value is going to be the empty vector, so your x will be empty, so your y will be empty.

Are you sure you want n(k)-1 and not n(k)-n(k-1) or something similar ?

[Walter Roberson]

x=1:H(k):2 is creating a variable number of x values, depending on the step size, but y(1) + y(nH) assumes that there are at least as many values in x as length(H) which is a strange assumption.

When k = 1, then H(k) is H(1) = 0.5^2 = 0.25, and 1:0.25:2 will have the elements 1, 1.25, 1.5, 1.75, 2 -- a total of only 5 elements.

Perhaps instead of y(nH) you want y(end) ?

[Matt Ogden]

Changing this code to y(end) works perfectly. Thanks for your help!

Answer 2

[Jan]

You determine the stepsize h by:

h = (b-a) / (n(k)-1)

For n(1) = 0.5^2 = 0.25, you get:

h = -1.3333

Therefore x = 1:h:nH is empty and y is in consequence also.

It looks strange, that H and n have the same values.

Btw. "1:h:nH" is a vector already and including it in square brackets wastes time only.