> -----Original Message-----
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org]
> On Behalf Of Craig O'Connell
> Sent: Monday, October 11, 2010 8:10 AM
> To: alain.guil...@uclouvain.be
> Cc: r-help@r-project.org; pda...@gmail.com
> Subject: Re: [R] MATLAB vrs. R
>
>
> alain,
>
> Perhaps i'm still entering the code wrong. I tried using your
> result=myquadrature(f,0,2000)
> print(result)
>
> Instead of my:
> val = myquadrature(f,a,b)
> result=myquadrature(val,0,2000)
> print(result)
>
> ...and I am still getting an inf inf inf inf inf...
>
> Did you change any of the previous syntax in addition to changing the
> result statement?
>
> Thank you so much and I think my brain is fried! Happy Holiday.
>
> Craig
>
Craig,
I haven't seen an answer to this yet, so let me jump in. You seem to have some
stuff still leftover from MATLAB. Here is some cleaned up code that produces
the result you expect. I don't think the value of dx was being correctly
computed in your code. I did not change the assignment operator you used (=),
but in R the "preferred" operator is "<-" (without the quotes).
myquadrature <- function(f,a,b){
npts = length(f)
nint = npts-1
if(npts <= 1) error('need at least two points to integrate')
if(b <= a) error('something wrong with the interval, b should be greater than
a') else dx=b/nint
sum(f[-npts]+f[-1])/2*dx
}
#Call my quadrature
x = seq(0,2000,10)
h = 10*(cos(((2*pi)/2000)*(x-mean(x)))+1)
u = (cos(((2*pi)/2000)*(x-mean(x)))+1)
a = x[1]
b = x[length(x)]
plot(x,-h)
a = x[1];
b = x[length(x)]
#call your quadrature function. Hint, the answer should be 30000.
f = u*h
result = myquadrature(f,a,b)
result
Hope this is helpful,
Dan
Daniel Nordlund
Bothell, WA USA
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
