Hi - I was just looking at this but this line gives me an error:
red=: points&(p. chi2) D.1 newt 0.5 1000 ] 1 1 1 1 1 1 1 1 1 1 |domain error | red=: points&(p.chi2)D.1 newt 0.5 1000]1 1 1 1 1 1 1 1 1 1 Am I missing something here? Thanks, Devon On Mon, Feb 15, 2016 at 11:27 AM, Louis de Forcrand <[email protected]> wrote: > Thanks everyone! > > The lab was pretty much exactly what I was looking for, and indeed > matrix division was the way to go. For satisfactory results, ~1000 > iterations are needed in some cases, but they are spectacular: > > > > The code is probably minimal, but it’s surprisingly versatile; the code > used to generate the curves is: > > newt=: 2 : '(- (0{::n) * u %. u D.1)^:(1{::n)’ NB. fetch is used here to > be able to box the 2nd parameter > NB. for > a list of iterations > chi2=: 1 : '[: +/@:*: {:@[ - ] u {.@[‘ > points=: (i: 5) ,: ? 10 $ 10 > green=: points&(p. chi2) D.1 newt 0.5 1000 ] 1 1 1 1 1 > red=: points&(p. chi2) D.1 newt 0.5 1000 ] 1 1 1 1 1 1 1 1 1 1 > > The optimised parameters were then stored in green or red, and I simply > plotted the two functions on top of the points with the plot package. > > Theoretically (ignoring possible precision problems) any function which > takes parameters on its left > and a value to its right can be used as a placeholder to minimise the > square of the differences with > the points to fit; with only two operators! > I really am glad I decided to look into APL and then J; you won’t see me > doing this any time soon in C. > > Thanks again for your generous help, > Louis > > > > On 14 Feb 2016, at 15:14, Ben Gorte - CITG <[email protected]> > wrote: > > > > Hello Louis, > > > > did you have a look at the Best Fit lab in the math section of Help -> > Studio -> Labs ? > > I think towards the end it is getting pretty close to what you are > looking for. > > > > Ben > > > > ________________________________________ > > From: Programming [[email protected]] on behalf > of Louis de Forcrand [[email protected]] > > Sent: Saturday, February 13, 2016 16:50 > > To: [email protected] > > Subject: [Jprogramming] Non-scalar Newton-Raphson > > > > I’ve been trying to write a conjunction that will find the zeros to a > function using the > > Newton-Raphson method. The simplest way to do this is probably: > > > > English: > > x_n+1 = x_n - f ( x_n ) / f ‘ ( x_n ) > > > > J: > > eznewt=: 2 : ‘ ( - u % u D.1 ) ^: n ‘ > > > > This works fine for scalar -> scalar functions, but won’t work if the > rank of the > > result is AND the rank of the arguments are above 1. > > Probably the most evident situation where this would be a problem is if > one > > were searching for a minimum instead of a zero, in which case the > algorithm > > would be applied to the derivative of the function: +/@:*: D.1 newt 30 ] > 3 4 > > As I understand it, this would give a result shape at each iteration > (with a function u) of: > > > > ( $y ) , $ u eznewt 1 y > > > > where A is the original argument, or the result of the previous > iteration. What we would > > want is a result with shape $y . > > > > First, let’s get some rules clear: > > - the syntax should be: (verb) newt (# of iterations) (argument) > > - the argument can be of any rank > > - the shape of the argument matches the shape of the result at any > iteration > > - this implies that $ u y matches i. 0 or $ y > > > > The hard part is getting u % u D.1 to have the shape of the argument. > > If y is a vector, and u y is a scalar, then u D.1 y will be a vector, and > > u D.1 D.1 y will be a $y by $y matrix. All that is needed then is to take > > its diagonal with (<0 1)&|: . > > > > But what if y is a matrix? Since $ u D.1 y -: ( $ y ) , $ u y , I though > maybe running > > y's axis together with |: might work > > > > newt=: 2 : '(- u diag_and_divide u D.1)^:n’ > > diag_and_divide=: [ % ] |:~ -@#@$@[ <@}. i.@#@$@] > > > > but something about it doesn’t. > > My head is $pinning now, and I figured I’d send this and then take a > break. > > > > Thanks in advance! > > Louis > > ---------------------------------------------------------------------- > > For information about J forums see http://www.jsoftware.com/forums.htm > > ---------------------------------------------------------------------- > > For information about J forums see http://www.jsoftware.com/forums.htm > > ---------------------------------------------------------------------- > For information about J forums see http://www.jsoftware.com/forums.htm > -- Devon McCormick, CFA Quantitative Consultant ---------------------------------------------------------------------- For information about J forums see http://www.jsoftware.com/forums.htm
