Hi all,
I have written a function for a part of a larger program. This
program generates a numeric integral for a function that it calls. I
would like to allow the outter function to pass several optional
arguments (between zero and seven) on to the called function which is
being integrated, which could then pass those arguments on to another
lower-level function...etc. The code for the integration function is
below. In the pasted example, I integrate the function named "trial"
at bottom. As a very simple example, how could I get the overall
function, NDint, to pass the z parameter on to "trial"?
Thanks,
Nathan
//An n-slice integration routine using simpson's approximation
mata
mata clear
mata set matalnum off
scalar NDint(
// to perform integration via simpson'sapproximation
pointer(real scalar function) scalar f,
//points to the functionto be integrated
real scalar lower,
//lower bound of integration
real scalar higher,
//upper bound
real scalar m) // where m = number of slices, must be an EVEN number
{
if ( (m/2) != (floor(m/2)) ) _error("m must be an EVEN number")
real scalar width //I like declarations, your mileage may vary
real colvector Y
real rowvector X
real scalar n //slices plus one
n = m+1
if (lower > higher) _error("limits of integration reversed")
//make sure upper > lower NB: better to use _error() than
return()
width = (higher-lower) / m
//calculates width of each "slice"
Y = J(n,1,.)
//creates an empty 15 x 1 column vector
for (i=1; i<=n; i++) {
x = lower + ((i-1)*width)
Y[i,1] = (*f)(x)
//places function evaluations directly into Y,
avoiding all the copying you were doing previously
}
X = J(1,n,.)
for (i=1;i<=n; i++) {
if (i==1) X[1,i] = 1
if (i==n) X[1,i] = 1
if ((i/2) == (floor(i/2))) X[1,i] = 4
if ( ((i+1)/2) == (floor( ((i+1)/2)) ) & (i!=1) & (i!=n) ) X[1,i] = 2
}
return(X*Y*width/3) //matrix multiplication makes this quicker
and more readable
}
scalar trial(real scalar x) {
//declares the output to be a scalar. otherwise
you could run into problems in the main program
return(z*exp(x))
}
f = &trial()
NDint(f, 0, 4, 60)
end
--
"Imagination is more important than knowledge..."
-- Albert Einstein
*
* For searches and help try:
* http://www.stata.com/help.cgi?search
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/
