Re: st: passing extra information to function evaluator program nl

["Brian P. Poi" <bpoi@stata.com>]

On Tue, 23 Mar 2010, Ben Zipperer wrote:

> I am interested in using Stata to find a real roots to a single
> polynomial several times, where the coefficients vary each time.  I
> basically want to apply Newton's method or something else to vectors
> of coefficients. I originally planned to restrict the solver at
> http://www.stata.com/support/faqs/lang/nl.html to a single polynomial
> and then iterate it, changing the coefficients at each iteration.
> However, I can't figure out how to pass the coefficients to the
> function evaluator program. How do I do this? Do I modify the syntax
> statement of the program? Or should I try a different strategy?
>
> For example, let's say the rows of the following matrix form the
> coefficients for a quadratic polynomial
>
> mat def M = (1,2,3\ 4,5,6\ 7,8,9)
>
> Then the nl program evaluator program is something like
>
> capture program drop nlpolysolve
> program nlpolysolve
>
>        syntax varlist(min=1 max=1) [if], at(name)
>
>        tempname x z
>        scalar `x' = `at'[1, 1]
>        scalar `z' = `at'[1, 2]
>
>        tempvar yh
> 	// the real polynomial
>        gen double `yh' = a_`i' * x^2 + b_`i' * `x' + c_`i' + 1 in 1
>
>        replace `yh' = 1 - `z' in 2
>        replace `varlist' = `yh'
>
> end
>
> And my iterations are something like
>
> forvalues i=1/3 {
> preserve
> clear
> set obs 2
> gen y = 0
> replace y = 1 in 1
>
> local a_`i' = M[`i',1]
> local b_`i' = M[`i',2]
> local c_`i' = M[`i',3]
>
> nl polysolve @ y, parameters(A B) initial(A 0.5 B 1)
> stuffiwanttopass(a_`i' b_`i' c_i')
>
> restore
> }
>
> How do I pass the entries of the matrix to the nlpolysolve program above?

Ben,

Just make your program accept an additional option; any options that 
-nl- doesn't recognize as belonging to it will get passed to your program. 
Try something like this:

capture program drop nlpolysolve
program nlpolysolve

        syntax varlist(min=1 max=1) [if], at(name) mycoefs(numlist)

	local a : word 1 of `mycoefs'
	local b : word 2 of `mycoefs'
	local c : word 3 of `mycoefs'
...

        gen double `yh' = `a' * x^2 + `b' * `x' + `c' + 1 in 1

end


...

local a_`i' = M[`i',1]
local b_`i' = M[`i',2]
local c_`i' = M[`i',3]


nl polysolve @ y, parameters(A B) initial(A 0.5 B 1) ///
                  mycoefs(`a_`i'' `b_`i'' `c_`i'')



   -- Brian Poi
   -- bpoi@stata.com
*
*   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/