You are right on target - you need a matrix of derivatives which is
called a Jacobian matrix
http://mathworld.wolfram.com/Jacobian.html
In the case of an example: that of a power system flow of real and
reactive power as related to voltage phase and magnitude The Newton
Raphson approach involves the Jacobian matrix.
This approach has been used since the late 1950's for systems with up
to000 to 10000 variables (something that one numerical analysis text
dismissed, erroneously, as not being a realistic problem ) and
provides fast convergence (partially due to the sin/cos differential cycle.
I am not a mathematician but an engineer using this as a tool. I have
written this in Turbo Basic as well as in APL (STSC version11-no
complex number capability). The APL version took about 11 lines of
actual solution with added lines of checking limits and dressing up
output. The other version involves about 3 pages of code. I have
partially developed a J version but overall I may not have the skill to
make it more compact than the APL coding.
https://en.wikipedia.org/wiki/Power-flow_study
Don Kelly
On 2/13/2016 8:00 AM, Roger Hui wrote:
If the atoms of the function are independent, then it's easiest to do f"0
Newton (or f Newton"0) and be done with it. If they are not independent,
that is f is what's called a vector function in conventional mathematical
notation, then you need to use %. (matrix divide) instead of % (plain old
divide), among other things.
My knowledge of vector calculus is minimal. Perhaps the mathematicians in
this Forum can help you if you really do have a vector function.
On Sat, Feb 13, 2016 at 7: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
