http://code.jsoftware.com/wiki/Vocabulary/dcapdot <http://code.jsoftware.com/wiki/Vocabulary/dcapdot>
Towards the bottom it says that u D.n with n > 1 is less accurate than u D.1 D.1 . I am assuming that u D.1 D.1 isn’t very accurate either though, as your example is correct (2nd der. of +/@:*: is 2: ). ——— Indeed, multinewt does work (better than other options in fact) for +/@:*: . Even though vsnewt seems to be the best so far, as you said it won’t work if the function returns a vector with the same shape as its argument (it doubles its shape at each iteration). Making it work would involve running some axis together (taking diagonals) to eliminate partial derivatives which result from differentiating with respect to two different elements of the argument. I can’t figure out which ones though… Louis > On 14 Feb 2016, at 01:13, Jose Mario Quintana <[email protected]> > wrote: > > Consider the following, > > (+/@:*:D.1 D.1) 0 0 > 1 0 > 0 1 > (+/@:*:D.1 D.1) 0 111 > 0 0 > 0 1.9861 > > (+/@:*:D.2 ) 0 0 > 2 0 > 0 2 > (+/@:*:D.2 ) 0 111 > 0 _304.804 > 1.00391 1.99305 > > As far as I know, the second derivative of the sum of the squares is, in > principle, two times the identity matrix regardless of the arguments. Am I > missing something? > > > > On Sat, Feb 13, 2016 at 4:09 PM, Marshall Lochbaum <[email protected]> > wrote: > >> If u is (+/@:*:), then (u D. 1) returns a vector, not a matrix. From >> your first email it looked like you wanted to find a zero of >> (+/@:*: D. 1), that is, minimize (+/@:*:). multinewt works for that >> case: >> >> +/@:*: D. 1 multinewt (<5) 3 4 >> 3 4 >> _0.142552 _0.0529989 >> 0.00158427 0.00128811 >> _1.42249e_5 2.08734e_5 >> _2.57043e_7 _1.37868e_7 >> >> There's a more complicated formula if the function's range has fewer >> dimensions than its domain. For a scalar function of a vector like >> (+/@:*:) it is >> >> v_{n+1} = v_n - f(v_n) * f'(v_n) / |f'(v_n)|^2 >> >> or >> >> vsnewt =: 2 : '(- u * [: (% +/@:*:) u D.1)^:n' >> >> +/@:*: vsnewt (<5) 3 4 >> 3 4 >> 1.5 2 >> 0.75 1 >> 0.375 0.5 >> 0.1875 0.25 >> >> This converges slowly because the root has order two. >> >> Marshall >> >> On Sat, Feb 13, 2016 at 09:32:15PM +0100, Louis de Forcrand wrote: >>> I’m sorry about that Raul, that should read “y”. >>> >>> ---- >>> >>> Thanks Roger. You were right, after a quick search, I’ve found that >>> indeed I would need to perform a matrix division. Several documents >>> give this as a multi-variable Newton method: >>> >>> v_n+1 = v_n - J_f (v_n) ^-1 f(v_n) >>> >>> In J, this would be: >>> >>> multinewt=: 2 : ‘(- u %. u D.1)^:n’ >>> >>> However, if I take a paraboloid z =: +/@:*: x y , and look for its zeros >> with >>> >>> +/@:*: multinewt (<5) 3 4 >>> >>> I get this: >>> >>> 3 4 >>> _0.5 0.5 >>> _0.5 0.5 >>> _0.5 0.5 >>> _0.5 0.5 >>> >>> which is obviously wrong. Would I need some sort of transposition? This >>> could just be due to precision errors. Or because I wrote an erroneous >>> conjunction? >>> >>> Cheers, >>> Louis >>> >>>> On 13 Feb 2016, at 20:52, Raul Miller <[email protected]> wrote: >>>> >>>> Something is wrong here. >>>> >>>> You say "where A is the original argument, or the result of the >>>> previous iteration." >>>> >>>> ... but I can find no use of A in the expression you are describing. >>>> >>>> Can you look over what you wrote and decide if it is what you meant to >> say? >>>> >>>> Thanks, >>>> >>>> -- >>>> Raul >>>> >>>> On Sat, Feb 13, 2016 at 10:50 AM, Louis de Forcrand <[email protected]> >> wrote: >>>>> 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 >> ---------------------------------------------------------------------- >> 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
