Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
Error message using trapz

Subject: Error message using trapz

From: Arie Driga

Date: 29 Sep, 2010 13:30:26

Message: 1 of 19

Hi all,

I would like to integrate two function

y1 = a*t - d cos (w*t); a, d, and w are constant
y2 = sin wt * y1

my t is: [0 0.1 0.2 0.3 0.4 0.5],--> just an example

basically i want to integrate:
0 to 0-1
0 to 0.2 ...
0 to 0.5

What i did:

intgr1=zeros(1,length(t));
intgr2=zeros(1,length(t));

for i=1:size(a,2)
    x=0:i*1/100:t(i);
    y1=exp(a*x-d*cos(w*x));
    y2=sin(w*x).*y1;
    intgr1(i)=trapz(x,y1);
    intgr2(i)=trapz(x,y2);
end

And the error messages are:

??? Subscript indices must either be real positive integers or logicals.

Error in ==> ipermute at 23
inverseorder(order) = 1:numel(order); % Inverse permutation order

Error in ==> trapz at 78
if ~isempty(perm), z = ipermute(z,perm); end

Error in ==> time_response3 at 64
    intgr1(i)=trapz(x,y1);

How to solve this?

Thanks in advance

Subject: Error message using trapz

From: dpb

Date: 29 Sep, 2010 14:12:44

Message: 2 of 19

Arie Driga wrote:
...
> I would like to integrate two function
>
> y1 = a*t - d cos (w*t); a, d, and w are constant
> y2 = sin wt * y1
>
> my t is: [0 0.1 0.2 0.3 0.4 0.5],--> just an example
>
> basically i want to integrate:
> 0 to 0-1
> 0 to 0.2 ...
> 0 to 0.5
> What i did:
>
> intgr1=zeros(1,length(t));
> intgr2=zeros(1,length(t));
>
> for i=1:size(a,2)

a, d, and w are constant --> i==1
> x=0:i*1/100:t(i);

x = 0:1*0.01:t(1) --> x = 0:0.01:0

> y1=exp(a*x-d*cos(w*x));
> y2=sin(w*x).*y1;
> intgr1(i)=trapz(x,y1);
> intgr2(i)=trapz(x,y2);
> end
>
> And the error messages are:
>
> ??? Subscript indices must either be real positive integers or logicals.
...

t(i) is a value, not a vector

trapz() works as an integrator over the argument vector so you'll need
to build a time vector over the integration length for each upper limit
and evaluate the function for that vector and pass the result to trapz()

doc trapz

--

Subject: Error message using trapz

From: dpb

Date: 29 Sep, 2010 14:21:08

Message: 3 of 19

dpb wrote:
...

> t(i) is a value, not a vector
>
> trapz() works as an integrator over the argument vector ...

Sorry, mixed metaphors here...the two-argument version does compute the
integral of y wrt x but it still requires vectors or a vector/array and
you've passed a single point...

--

Subject: Error message using trapz

From: Arie Driga

Date: 29 Sep, 2010 14:48:21

Message: 4 of 19

dpb <none@non.net> wrote in message <i7vi7h$ogs$1@news.eternal-september.org>...
> dpb wrote:
> ...
>
> > t(i) is a value, not a vector
> >
> > trapz() works as an integrator over the argument vector ...
>
> Sorry, mixed metaphors here...the two-argument version does compute the
> integral of y wrt x but it still requires vectors or a vector/array and
> you've passed a single point...
>
> --

Sorry dpb but I don't understand what you are trying to say.
Is it because my initial value of t is 0?

Subject: Error message using trapz

From: dpb

Date: 29 Sep, 2010 16:33:32

Message: 5 of 19

Arie Driga wrote:
> dpb <none@non.net> wrote in message
> <i7vi7h$ogs$1@news.eternal-september.org>...
>> dpb wrote:
>> ...
>>
>> > t(i) is a value, not a vector
>> > > trapz() works as an integrator over the argument vector ...
>>
>> Sorry, mixed metaphors here...the two-argument version does compute
>> the integral of y wrt x but it still requires vectors or a
>> vector/array and you've passed a single point...
>>
>> --
>
> Sorry dpb but I don't understand what you are trying to say.
> Is it because my initial value of t is 0?

No, it's because t(i) is a reference to an element of an array, not a
vector.

Try

size(t(1))

or

t(1)

to demonstrate vis a vis

size(t)

or

t

at the command line...

--

Subject: Error message using trapz

From: Arie Driga

Date: 29 Sep, 2010 17:17:05

Message: 6 of 19

dpb <none@non.net> wrote in message <i7vpvk$8c9$1@news.eternal-september.org>...
> Try
>
> size(t(1))
>
> or
>
> t(1)
>
> to demonstrate vis a vis
>
> size(t)
>
> or
>
> t
>
> at the command line...
>
> --
>

Ok, back to what I want to do:
I will just say
Lets say integral from 0 to 0.3 equal to (integral 0 to 0.1)+(integral 0.1 to 0.2)+(integral 0.2 to 0.3)

Here is what I did:
ic=time_increment; %in my example is 0.1
int_icr=integral_increment;

stp=ic/int_icr;
x = zeros((length(t)-1),(int_icr));
y1 = zeros((length(t)-1),(int_icr));
y2 = zeros((length(t)-1),(int_icr));
int1 = zeros(1,(length(t)-1));
int2 = zeros(1,(length(t)-1));

for i=1:(length(t)-1)
    for j=1:(int_icr+1)
    x(i,j) = t(1,i)+(j-1)*stp;
    y1(i,j) = exp(a*x(i,j)-d*cos(w*x(i,j)));
    y2(i,j) = y1(i,j) + sin(w*x(i,j));
    end
end
x
y1
for i=1:size(x,1)
    int1(1,i) = trapz(x(i,:),y1(i,:));
    int2(1,i) = trapz(x(i,:),y2(i,:));
end

int1_sum=sum(int1);
int2_sum=sum(int2);

I got answer, but I do not know if it is correct or not.

Subject: Error message using trapz

From: AJP

Date: 29 Sep, 2010 17:27:22

Message: 7 of 19

"Arie Driga" <a_driga@yahoo.com> wrote in message <i7vf1i$1s$1@fred.mathworks.com>...
> Hi all,
>
> I would like to integrate two function
>
> y1 = a*t - d cos (w*t); a, d, and w are constant
> y2 = sin wt * y1
>
> my t is: [0 0.1 0.2 0.3 0.4 0.5],--> just an example
>
> basically i want to integrate:
> 0 to 0-1
> 0 to 0.2 ...
> 0 to 0.5
>
> What i did:
>
> intgr1=zeros(1,length(t));
> intgr2=zeros(1,length(t));
>
> for i=1:size(a,2)
> x=0:i*1/100:t(i);
> y1=exp(a*x-d*cos(w*x));
> y2=sin(w*x).*y1;
> intgr1(i)=trapz(x,y1);
> intgr2(i)=trapz(x,y2);
> end
>
> And the error messages are:
>
> ??? Subscript indices must either be real positive integers or logicals.
>
> Error in ==> ipermute at 23
> inverseorder(order) = 1:numel(order); % Inverse permutation order
>
> Error in ==> trapz at 78
> if ~isempty(perm), z = ipermute(z,perm); end
>
> Error in ==> time_response3 at 64
> intgr1(i)=trapz(x,y1);
>
> How to solve this?
>
> Thanks in advance
==================================================

Your two functions are, written out in full:

> y1 = a*t - d cos (w*t)
> y2 = sin wt * a*t - d*cos(w*t)

Since these functions are known and continuous, the best way to get an accurate value of the integtal is to use "quad", not "trapz". Quad is an adaptive quadrature algorithm, which means it will automatically vary the step-spacing in t for maximum accuracy and efficiency.

To use quad you will need to write a function file for each of your functions, something like:

function y1=func1(t,a,d,w)
y1 = a*t - d cos (w*t)
end

and

function y2=func2(t,a,d,w)
y2 = sin w*t. * a*t - d*cos(w*t) % Make sure you want to do an elementwise operation here, i.e. using " .* " not " * ".
end

Then call quad, for example on func1, by using:
a=const.; d=const.; w=const.;
intgr1=quad(@(t)func1(t,a,d,w),0,0.1)

The syntax here basically states that "t" is your independent variable, func1 is the function being integrated, which is a function of t and the constants a, d and w which you passed to the function, and the limits of integration are t=0 and t=0.1

If you want to vary the limits and output "intgr" each time, then an easy way would be to loop over the values in your current t array:
t=[0 0.1 0.2 0.3 0.4 0.5]
On each iteration of the loop, say for n=1:length(t), put t(n) as the upper limit of integration when you call quad.

PS. I may have made a couple of typos, so please ensure the syntax is right by looking at the help files e.g. type "help quad" at the command prompt.

Subject: Error message using trapz

From: AJP

Date: 29 Sep, 2010 17:44:04

Message: 8 of 19

Yeah, I've just seen a mistake there - that sould be t=[0.1,0.2... etc.] since you always want to use 0 as your lower bound you don't need to include it.

Subject: Error message using trapz

From: Arie Driga

Date: 29 Sep, 2010 21:19:04

Message: 9 of 19

"AJP " <adam.noakes07@imperial.ac.uk> wrote in message <i7vtt4$cua$1@fred.mathworks.com>...
> Yeah, I've just seen a mistake there - that sould be t=[0.1,0.2... etc.] since you always want to use 0 as your lower bound you don't need to include it.

%=================================================

Hi AJP

I have tried the quad:

function y1=funct1(t,a,d,w)
y1 = a*t - d*cos(w*t)
end

    function y2=funct2(t,a,d,w)
        y2 = sin(w*t).* a*t - d*cos(w*t)
    end

for i=2:length(t)
    intgr1 = quad(@(t)func1(t,a,d,w),0,t(i));
    intgr2 = quad(@(t)func2(t,a,d,w),0,t(i));
end

intgr1

the above code is written inside another function, will it be a problem to do so?
I tried it, run it and got error:

??? Undefined function or method 'func1' for input arguments of type 'double'.

Error in ==> time_response4>@(t)func1(t,a,d,w)


Error in ==> quad at 77
y = f(x, varargin{:});

Error in ==> time_response4 at 65
    intgr1 = quad(@(t)func1(t,a,d,w),0,t(i));

Could you help me on this?

Subject: Error message using trapz

From: AJP

Date: 1 Oct, 2010 16:39:19

Message: 10 of 19

"Arie Driga" <a_driga@yahoo.com> wrote in message <i80ag8$db7$1@fred.mathworks.com>...
> "AJP " <adam.noakes07@imperial.ac.uk> wrote in message <i7vtt4$cua$1@fred.mathworks.com>...
> > Yeah, I've just seen a mistake there - that sould be t=[0.1,0.2... etc.] since you always want to use 0 as your lower bound you don't need to include it.
>
> %=================================================
>
> Hi AJP
>
> I have tried the quad:
>
> function y1=funct1(t,a,d,w)
> y1 = a*t - d*cos(w*t)
> end
>
> function y2=funct2(t,a,d,w)
> y2 = sin(w*t).* a*t - d*cos(w*t)
> end
>
> for i=2:length(t)
> intgr1 = quad(@(t)func1(t,a,d,w),0,t(i));
> intgr2 = quad(@(t)func2(t,a,d,w),0,t(i));
> end
>
> intgr1
>
> the above code is written inside another function, will it be a problem to do so?
> I tried it, run it and got error:
>
> ??? Undefined function or method 'func1' for input arguments of type 'double'.
>
> Error in ==> time_response4>@(t)func1(t,a,d,w)
>
>
> Error in ==> quad at 77
> y = f(x, varargin{:});
>
> Error in ==> time_response4 at 65
> intgr1 = quad(@(t)func1(t,a,d,w),0,t(i));
>
> Could you help me on this?
===================================================

I politely suggest you read up about using functions - it's pretty fundamental to using MATLAB.

Type "doc function" at the command prompt or look at this webpage (both give the same info): http://www.mathworks.com/help/techdoc/ref/function.html

You should almost always save your function in its own .m file, unless you are using a subfunction, but this is an advanced technique.

Use the editor in matlab, start a new file write each function e.g. for func1:

function y1=funct1(t,a,d,w)
y1 = a*t - d*cos(w*t)
end

then save, ensuring that the name of the function is the same as the filename, "func1". Then do the same for "func2".

Then in a new .m file, i.e. a script file that you use to run your functions, write out your parameters and call the function as above.

Don't forget that if you want to "save" the values of inter1 and intgr2 on each pass of your loop over the values in t, you need to use the following syntax in your script file:

for i=2:length(t)
     intgr1(i) = quad(@(t)func1(t,a,d,w),0,t(i));
     intgr2(i) = quad(@(t)func2(t,a,d,w),0,t(i));
end

Note the inclusion of (i) after intgr1 and intgr2 - this means that on iteration, i , the value of intgr1 and intgr2 will be stored in two separate vectors, i elements in length.

The first element, for the iteration i=2 in your case, will correspond to the value of intgr1 / intgr2 for limits t=0 and t=the 2nd element in t. The next iteration will give you the values for the limits 0 to the 3rd element in t and so on.


Give it a go! Don't give up - you learn a lot by trying it for yourself!

Subject: Error message using trapz

From: AJP

Date: 1 Oct, 2010 18:23:04

Message: 11 of 19

=================================================

OK, I took another look and noticed the following error.

The function files are written like this, yes?

function y1=funct1(t,a,d,w)
y1 = a*t - d*cos(w*t)
end

function y2=funct2(t,a,d,w)
y2 = sin(w*t).* a*t - d*cos(w*t)
end

This defines functions "funct1" and funct2"

However, in your script file, when you're calling the functions in quad, the following is written:

for i=2:length(t)
intgr1(i) = quad(@(t)func1(t,a,d,w),0,t(i));
intgr2(i) = quad(@(t)func2(t,a,d,w),0,t(i));
end

Notice how quad calls the functions "func1" and "func2" - NOT "funct1" and "funct2" as defined above.

You have to make sure you calll the same function name exactly!

I warned you to not to copy and paste my code and to check it for typos before using it - this is the sort of debugging you should do yourself. The error message "undefined function func1" is letting you know that the function called has not been defined.

Subject: Error message using trapz

From: Arie Driga

Date: 2 Oct, 2010 00:18:04

Message: 12 of 19

"AJP " <adam.noakes07@imperial.ac.uk> wrote in message <i858u8$e2n$1@fred.mathworks.com>...

> I warned you to not to copy and paste my code and to check it for typos before using it - this is the sort of debugging you should do yourself. The error message "undefined function func1" is letting you know that the function called has not been defined.

Yes you are right, it was my bad.
I do have doubt now, I have fixed the method using trapz. It works now, but the answer between trapz method and quad method is differ by order of magnitude.

Here is for trapz:
%=================================================
%integration

stp = ic/int_icr;
x = zeros((length(t)-1),(int_icr));
y1 = zeros((length(t)-1),(int_icr));
y2 = zeros((length(t)-1),(int_icr));
int1 = zeros(1,(length(t)-1));
int2 = zeros(1,(length(t)-1));

for i=1:(length(t)-1)
    for j=1:(int_icr+1)
    x(i,j) = t(1,i)+(j-1)*stp;
    y1(i,j) = exp(a*x(i,j)-d*cos(omg*x(i,j)));
    y2(i,j) = exp(a*x(i,j)-d*cos(omg*x(i,j))) * sin(omg*x(i,j));
    end
end

for i=1:size(x,1)
    int1(1,i) = trapz(x(i,:),y1(i,:));
    int2(1,i) = trapz(x(i,:),y2(i,:));
end


int1_sum = zeros(1,length(int1));
int2_sum = zeros(1,length(int1));
for i=1:length(int1)
    int1_sum(1,i) = sum(int1(1:i));
    int2_sum(1,i) = sum(int2(1:i));
end

int1_new=[0 int1_sum];
int2_new=[0 int2_sum];
%==============================================

omg for trapz is the same as w in quad:

And here is the quad method based on your suggestion:

%==============================================

%integration
function y1=funct1(t,a,d,w)
y1 = a*t - d*cos(w*t);
end

    function y2=funct2(y1,t,d,w)
        y2 = y1 - d.*cos(w*t);
    end

intgr1=zeros(1,length(t));
intgr2=zeros(1,length(t));
for i=2:length(t)
    intgr1(i) = quad(@(t)funct1(t,a,d,w),0,t(i));
    intgr2(i) = quad(@(t)funct2(t,a,d,w),0,t(i));
end

%===============================================

I plotted the integral value both for quad and trapz method, at the result (as I mentioned earlier) were much different, and I really do not know how it happens.

Subject: Error message using trapz

From: AJP

Date: 4 Oct, 2010 10:35:20

Message: 13 of 19

"Arie Driga" <a_driga@yahoo.com> wrote in message <i85tns$aks$1@fred.mathworks.com>...
> "AJP " <adam.noakes07@imperial.ac.uk> wrote in message <i858u8$e2n$1@fred.mathworks.com>...
>
> > I warned you to not to copy and paste my code and to check it for typos before using it - this is the sort of debugging you should do yourself. The error message "undefined function func1" is letting you know that the function called has not been defined.
>
> Yes you are right, it was my bad.
> I do have doubt now, I have fixed the method using trapz. It works now, but the answer between trapz method and quad method is differ by order of magnitude.
>
> Here is for trapz:
> %=================================================
> %integration
>
> stp = ic/int_icr;
> x = zeros((length(t)-1),(int_icr));
> y1 = zeros((length(t)-1),(int_icr));
> y2 = zeros((length(t)-1),(int_icr));
> int1 = zeros(1,(length(t)-1));
> int2 = zeros(1,(length(t)-1));
>
> for i=1:(length(t)-1)
> for j=1:(int_icr+1)
> x(i,j) = t(1,i)+(j-1)*stp;
> y1(i,j) = exp(a*x(i,j)-d*cos(omg*x(i,j)));
> y2(i,j) = exp(a*x(i,j)-d*cos(omg*x(i,j))) * sin(omg*x(i,j));
> end
> end
>
> for i=1:size(x,1)
> int1(1,i) = trapz(x(i,:),y1(i,:));
> int2(1,i) = trapz(x(i,:),y2(i,:));
> end
>
>
> int1_sum = zeros(1,length(int1));
> int2_sum = zeros(1,length(int1));
> for i=1:length(int1)
> int1_sum(1,i) = sum(int1(1:i));
> int2_sum(1,i) = sum(int2(1:i));
> end
>
> int1_new=[0 int1_sum];
> int2_new=[0 int2_sum];
> %==============================================
>
> omg for trapz is the same as w in quad:
>
> And here is the quad method based on your suggestion:
>
> %==============================================
>
> %integration
> function y1=funct1(t,a,d,w)
> y1 = a*t - d*cos(w*t);
> end
>
> function y2=funct2(y1,t,d,w)
> y2 = y1 - d.*cos(w*t);
> end
>
> intgr1=zeros(1,length(t));
> intgr2=zeros(1,length(t));
> for i=2:length(t)
> intgr1(i) = quad(@(t)funct1(t,a,d,w),0,t(i));
> intgr2(i) = quad(@(t)funct2(t,a,d,w),0,t(i));
> end

==============================================

OK, I have some general observations, hopefully these will help you.

(1) The functions you have in your trapz method and those in the quad method are totally different.

In the trapz method you have writteN:
> y1(i,j) = exp(a*x(i,j)-d*cos(omg*x(i,j)));
> y2(i,j) = exp(a*x(i,j)-d*cos(omg*x(i,j))) * sin(omg*x(i,j));

Whereas in the functions called by quad are:

> y1 = a*t - d*cos(w*t);
> y2 = y1 - d.*cos(w*t);

Notice there is no exp() in y1 or y2 here, so of course the results will be different. You really must pay attention to what you are writing!

Check "doc quad" at the command line.

(2) Whenever you use trapz, you must always multiply by the step length. You have forgotten to do this, so the trapz function has assumed your step length is 1. This probably explains why the results for quad and trapz are very different.

Look at "doc trapz"

As an example, say I want to find the integtral of sin(x) between 0 and pi. If you do this by hand you will find the real answer is exactly 2.

To do this in MATLAB, I write the following code at the command line or in an .m file:

%====================================
format long % this makes MATLAB display all numbers to their fulll number of places

x=linspace(0,pi,756); % note 756 is just some arbitrary number of elements I chose, it could be any number

dx=x(2)-x(1) % gives the step length

y=sin(x);

intgr=trapz(y)*dx % (=1.99998)
%====================================

If you forget to put " *dx " after "trapz" then the result is 218.99, clearly wrong.


(3) quad is much preferable to trapz in this case as it is more accurate and you don't need to generate any vectors for your y1, y2 and x as in your trapz method.

(4) I said before that you have to save your function files func1 and func2 in separate files, and then a separate script file that calls quad and your own functions. The code you put above suggests you have them all in the same file, which is not right. This may also be why your results are unexpected.

(5) I don't understand why you need the double for loop in your trapz method - if you are only looping over the values in a row vector then you only need a single loop, passing over each element in the row in turn. It is likely this method is either generating incorrect results or is simply inefficient, creating values that arene't needed.

(5) Using loops to generate new vectors or matrices, as in your trapz method, is less efficient than using elementwise operators like .* ./ .^
These operators allow you to perform operations on each element of a matrix or vector.

For example, say I have two 2x2 matrices, x=[1,2;3,4] and y=[4,5;6,7]. If I type x*y, I get ans=[16,19;36,43]. This is like a normal matrix multiplication.

However, if I write x.*y, I get ans=[4,10;18,28] which is like I have multiplied each element in x with the corresponding element in y.

This is much much more efficient than using a loop like this:

for n=1:numel(x)
ans(n)=x(n)*y(n)
end

Search online for "matllab element wise operations" or "matlab vectorized operations" for more info.

Subject: Error message using trapz

From: AJP

Date: 4 Oct, 2010 13:30:30

Message: 14 of 19

=====================================

Further to my previous post I have made some .m files for both the quad and trapz methods, for the functions in your original post.

Download them via this link:

https://fileexchange.imperial.ac.uk/files/4d09b7428ec/trapz_quad_comparison.zip



Hope this helps your understanding. But please read my previous post about debugging your code and making it more efficient.

Subject: Error message using trapz

From: Steven_Lord

Date: 4 Oct, 2010 15:10:28

Message: 15 of 19



"AJP " <adam.noakes07@imperial.ac.uk> wrote in message
news:i8cal8$bsb$1@fred.mathworks.com...
> "Arie Driga" <a_driga@yahoo.com> wrote in message
> <i85tns$aks$1@fred.mathworks.com>...
>> "AJP " <adam.noakes07@imperial.ac.uk> wrote in message
>> <i858u8$e2n$1@fred.mathworks.com>...

*snip*

> (2) Whenever you use trapz, you must always multiply by the step length.
> You have forgotten to do this, so the trapz function has assumed your step
> length is 1. This probably explains why the results for quad and trapz
> are very different.
>
> Look at "doc trapz"
>
> As an example, say I want to find the integtral of sin(x) between 0 and
> pi. If you do this by hand you will find the real answer is exactly 2.
>
> To do this in MATLAB, I write the following code at the command line or in
> an .m file:
>
> %====================================
> format long % this makes MATLAB display all numbers to their fulll number
> of places
>
> x=linspace(0,pi,756); % note 756 is just some arbitrary number of
> elements I chose, it could be any number
>
> dx=x(2)-x(1) % gives the step length
>
> y=sin(x);
>
> intgr=trapz(y)*dx % (=1.99998)

Rather than multiplying by the step length, let TRAPZ do that for you by
using the two input syntax, specifying not just the y data but also the x
data at which that y data was collected.

trapz(x, y)

See Example 1 in the reference page for TRAPZ.

http://www.mathworks.com/help/techdoc/ref/trapz.html

*snip*

--
Steve Lord
slord@mathworks.com
comp.soft-sys.matlab (CSSM) FAQ: http://matlabwiki.mathworks.com/MATLAB_FAQ
To contact Technical Support use the Contact Us link on
http://www.mathworks.com

Subject: Error message using trapz

From: AJP

Date: 4 Oct, 2010 20:27:33

Message: 16 of 19

=============================================

Stephen is right, using the two-argument input into trapz means you don't have to specify the step length yourself. I had only put this down to help our friend here understand of how trapz works and what it's doing.

Good point though - make sure to do this as it's more efficient.

Subject: Error message using trapz

From: AJP

Date: 5 Oct, 2010 12:27:28

Message: 17 of 19

If you have trouble getting the .m files from my previous link, they are now on the File Exchange at the following address:

http://www.mathworks.com/matlabcentral/fileexchange/28938-trapz-and-quad-basic-examples

Subject: Error message using trapz

From: Arie Driga

Date: 12 Oct, 2010 08:58:05

Message: 18 of 19

"AJP " <adam.noakes07@imperial.ac.uk> wrote in message <i8f5jg$iar$1@fred.mathworks.com>...
> If you have trouble getting the .m files from my previous link, they are now on the File Exchange at the following address:
>
> http://www.mathworks.com/matlabcentral/fileexchange/28938-trapz-and-quad-basic-examples

I have tried your code and it appears that not much different between trapz and quad. I have another question:
In the example, I mentioned a,d,w and what if a and d is not constant and the value is based on t so a and d is function of t:

I tried to modify your code:

for n=1:length(t) % Loop over the vector of upper limits, t
    
    intgr1_quad(n) = quad(@(t)funct1(t,a(n),d(n),w,ps),1e-6,t(n));
    % @(t) tells quad that t is the independent variable. The function
    % being integrated is "funct1", which is a function of t,a,d and w,
    % which we defined above. "0" means the lower limit of integration is
    % always zero and t(n) means the upper limit of integration is the
    % current (nth) value from the vector of upper limits, t.
    
    intgr2_quad(n) = quad(@(t)funct2(t,a(n),d(n),w,ps),1e-6,t(n)); % Same for func2.
  
end

It is giving me answer but I don't know if it is correct or not.
How to check my answer?

Thanks

Subject: Error message using trapz

From: AJP

Date: 15 Oct, 2010 18:10:19

Message: 19 of 19

It's usually a good idea to plot your functions to check they "look" right. If all else fails then you can at least estimate the integral graphically an compare this to the numerical integrals you have calculated using quad and trapz. If there is little difference between trapz and quad this is a good indication that your answers are right, but you are right to be careful about checking.

To make a simple plot, try this example:

x=linspace(0,2*pi,200);
y=sin(x);
figure %creates a new figure window
plot(x,y) %plots the graph

Another thing to rememebr is that even if your answers are perfect, you may not have defined the functions to be integrated correctly, so you have the right answers for the wrong function! Check your syntax carefully to ensure you haven't made any silly errors that may affect your answers. If you know the equation of the function then you should be able to have a good idea of the graph and the graph of the integral function.

Another way to check the answer, is to try checking out "doc cumtrapz". Cumtrapz is a function that will allow you to plot out the graph of the integral function. If you want to test it, try using it on, say, sin, and you should get a graph that looks like cos. If you are sure you have got the syntax correct you can then try it on your function, then you can compare the resulting graph to your original function, which is now like the derivative of this new graph. You will be able to see just from looking at the graph if the general shape is right or not, but you can further check by using the "gradient" function on the new graph, which should output some new alues that match your original function.

Alternatively, if the function is simple enough, do the integration by hand and plot the function you find. You can then do the same comparison as above, but with the true integral function, not just a numerical estimate. If the function is too hard to integrate by hand, or if the indefinite integral does not exist (like for the Gaussian function), then you could try out the Wolfram Integrator ( URL: http://integrals.wolfram.com/index.jsp ) which will give you an integral function, even if it only gives you the result in terms of the error function. MATLAB can evaluate this using the syntax erf(x).

If you want to change any of the constants to a function of t then simply write an annonymous function for them, save a function handle as the name of the variable and then pass this to the function in the usual way.

Let me explain...

For example, say at first we have function called "func" that uses some constants "a" and "b", and a vector of (say) time values, "t". Normally the function would be used by first assigning values to a and b, defining the range of the t vector then calling the funtion, say something like:

a=1; b=2; t=1:0.5:10; answer=func(a,b);

where, say, func is defined in its own m-file by something like...

function answer=func(a,b,t)
answer=a.*t+b.*t
end

This I'm sure you understand now.

However, say we want to make a or b into a function of t.

This can be easily done by redefining these variables as a function handle, say something like:

a_new=@(t) t.^2;
b_new=@(t) t.^3;

Now when we call func(a_new,b_new,t), MATLAB will use the new a and b in the function. The result in this case will be:

answer=a_new.*t+b_new.*t = (t.^2).*t + (t.^3).*t = t.^3 + t.^4

This general method is very powerful as it allows you to pass functions (or rather, the values calculated during a function) to other functions, allowing you to create more complicated inputs without having to write a huge argument list or a more complicated function.

Have a go. Also check out the help files on "annonymous function" and "function handle @"


I hope this is helpful.

Adam

=============================================
"Arie Driga" <a_driga@yahoo.com> wrote in message <i917ut$ps0$1@fred.mathworks.com>...
> "AJP " <adam.noakes07@imperial.ac.uk> wrote in message <i8f5jg$iar$1@fred.mathworks.com>...
> > If you have trouble getting the .m files from my previous link, they are now on the File Exchange at the following address:
> >
> > http://www.mathworks.com/matlabcentral/fileexchange/28938-trapz-and-quad-basic-examples
>
> I have tried your code and it appears that not much different between trapz and quad. I have another question:
> In the example, I mentioned a,d,w and what if a and d is not constant and the value is based on t so a and d is function of t:
>
> I tried to modify your code:
>
> for n=1:length(t) % Loop over the vector of upper limits, t
>
> intgr1_quad(n) = quad(@(t)funct1(t,a(n),d(n),w,ps),1e-6,t(n));
> % @(t) tells quad that t is the independent variable. The function
> % being integrated is "funct1", which is a function of t,a,d and w,
> % which we defined above. "0" means the lower limit of integration is
> % always zero and t(n) means the upper limit of integration is the
> % current (nth) value from the vector of upper limits, t.
>
> intgr2_quad(n) = quad(@(t)funct2(t,a(n),d(n),w,ps),1e-6,t(n)); % Same for func2.
>
> end
>
> It is giving me answer but I don't know if it is correct or not.
> How to check my answer?
>
> Thanks

Tags for this Thread

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us