Thanks for the help everyone. LAPACK is very useful. I also appreciate the Rosetta code mesh. It's interesting how many choices were made the same. I appreciate the help on this!
Brian Babiak MD https://drbabiak.com/ <http://drbabiak.com> paypal: *http://www.paypal.me/BrianBabiakMD <https://www.paypal.me/BrianBabiakMD>* cash, check, debit/credit cards, paypal, apple pay, android pay, discover and amex accepted for all payments square: *https://squareup.com/store/brian-babiak-md <https://squareup.com/store/brian-babiak-md>* On Fri, Jul 28, 2017 at 2:51 AM, Lippu Esa <esa.li...@varma.fi> wrote: > I think too. Here is the smallest eigenvalue of a: > > require 'math/lapack/geev' > ev=.geev_jlapack_ a > {.sort|>1{ev > 6.83897e_14 > > Esa > -----Original Message----- > From: Programming [mailto:programming-boun...@forums.jsoftware.com] On > Behalf Of 'Jon Hough' via Programming > Sent: Friday, July 28, 2017 9:17 AM > To: programm...@jsoftware.com > Subject: Re: [Jprogramming] Matrix operations questions > > I think a is singular. You should use LAPACK to decompose it perhaps. > see: http://code.jsoftware.com/wiki/Studio/LAPACK > > -------------------------------------------- > On Fri, 7/28/17, Brian Babiak <bdbab...@gmail.com> wrote: > > Subject: Re: [Jprogramming] Matrix operations questions > To: programm...@jsoftware.com > Date: Friday, July 28, 2017, 2:33 PM > > Very helpful! At least I can use > +/ .* for an iterative solution to Lp=I, solving for p. But > why can’t I always get an inverse? (Even with non > -infinite resistances, I get domain errors.) > > a=.L 2 2;3 1;1 4 > a > 4 _2 0 0 0 > _2 0 0 0 0 0 0 0 0 > 0 0 0 0 0 0 0 0 0 0 0 > _2 8 _2 0 0 _2 _2 0 0 > 0 0 0 0 0 0 0 0 > 0 0 0 0 0 0 0 0 > 0 _2 8 _2 > 0 0 _2 _2 0 0 0 0 0 > 0 0 0 0 0 0 0 0 0 0 0 > 0 > 0 0 _2 8 _2 0 0 _2 > _2 0 0 0 0 0 0 0 0 > 0 0 0 0 0 0 0 0 > 0 0 0 > _2 1004 0 0 0 _2 _1000 0 0 > 0 0 0 0 0 0 0 0 0 0 0 > 0 0 > _2 _2 0 0 0 8 _2 0 > 0 0 _2 0 0 0 0 0 > 0 0 0 0 0 0 0 0 0 > 0 _2 > _2 0 0 _2 12 _2 0 0 _2 _2 > 0 0 0 0 0 0 0 0 0 0 0 > 0 0 > 0 0 _2 _2 0 0 _2 1010 > _2 0 0 _2 _1000 0 0 0 0 > 0 0 0 0 0 0 0 0 > 0 0 0 > _2 _2 0 0 _2 2008 _1000 0 0 _1000 > _2 0 0 0 0 0 0 0 0 0 0 > 0 > 0 0 0 0 _1000 0 0 0 > _1000 4000 0 0 0 _1000 _1000 0 0 > 0 0 0 0 0 0 0 0 > 0 0 0 > 0 0 _2 _2 0 0 0 8 _2 > 0 0 0 _2 0 0 0 0 0 0 0 > 0 0 > 0 0 0 0 0 0 _2 > _2 0 0 _2 1009 _1000 0 0 _2 > _1 0 0 0 0 0 0 0 0 > 0 > 0 0 0 0 0 0 _1000 _1000 0 0 _1000 > 5999 _1000 0 0 _999 _1000 0 0 0 0 0 0 > 0 > 0 0 0 0 0 0 0 0 > _2 _1000 0 0 _1000 2008 _2 0 0 > _2 _2 0 0 0 0 0 0 > 0 0 0 > 0 0 0 0 0 0 _1000 0 0 > 0 _2 1006 0 0 0 _2 _2 0 0 0 0 > 0 > 0 0 0 0 0 0 0 0 > 0 0 _2 _2 0 0 0 7 _1 > 0 0 0 _2 0 0 0 0 > 0 0 > 0 0 0 0 0 0 0 0 0 _1 > _999 0 0 _1 1004 _1 0 0 _1 _1 0 0 > 0 > 0 0 0 0 0 0 0 0 > 0 0 0 0 _1000 _2 0 0 _1 > 1009 _2 0 0 _2 _2 0 0 > 0 0 0 > 0 0 0 0 0 0 0 0 0 > 0 _2 _2 0 0 _2 12 _2 0 0 _2 _2 > 0 > 0 0 0 0 0 0 0 0 > 0 0 0 0 0 0 _2 0 > 0 0 _2 8 0 0 0 _2 _2 > 0 0 > 0 0 0 0 0 0 0 0 0 0 > 0 0 0 _2 _1 0 0 0 5 _2 0 > 0 0 > 0 0 0 0 0 0 0 > 0 0 0 0 0 0 0 0 0 > _1 _2 0 0 _2 7 _2 0 0 > 0 > 0 0 0 0 0 0 0 0 0 0 > 0 0 0 0 0 0 _2 _2 0 0 _2 > 8 _2 0 > 0 0 0 0 0 0 0 > 0 0 0 0 0 0 0 0 0 > 0 0 _2 _2 0 0 _2 8 _2 > 0 > 0 0 0 0 0 0 0 0 0 0 > 0 0 0 0 0 0 0 0 _2 0 > 0 0 _2 4 > -/ . * a > _8.08478e26 > %. a > |domain error > | %.a > > > Brian Babiak > MD > https://drbabiak.com > > > 167 Ridgecrest Rd. > Ithaca > NY 14850 – 9449 > (607) 379 – 1335 > fax (866) 813 – 4158 > > Payment links: > > http://PayPal.Me/BrianBabiakMD http://mkt.com/brian-babiak-md > > > On 7/28/17, > 1:03 AM, "Programming on behalf of 'Jon Hough' > via Programming" <programming-boun...@forums.jsoftware.com > on behalf of programm...@jsoftware.com> > wrote: > > Matrix product > is +/ .* (dyadic) > determinant is -/ .* > (monadic) > > see: > http://code.jsoftware.com/wiki/Vocabulary/dot > > The matrix g is > singular (since the rows are identical) and has no inverse. > Hence %. gives an error. > > i.e. > > (-/ .*) g > > is zero. > > Hope that helps. > > > Jon > > > -------------------------------------------- > On Fri, 7/28/17, Brian Babiak <bdbab...@gmail.com> > wrote: > > > Subject: [Jprogramming] Matrix operations questions > To: programm...@jsoftware.com > Date: Friday, July 28, 2017, 1:56 PM > > I’m working on > a project to compute the > distance > resistance between two vertices in an electronic > grid. I will share that later. But I am > having domain error > issues with > matrix inverse and matrix product operations. > See the following: > > > > > f=.2 2$3 4 1 4 > > > g=.2 2$i.2 > > > h=.2 2$_3 4 2 _3 > > > > 2 > 3$'f';'g';'h';f;g;h > > > ┌───┬───┬─────┐ > > │f │g > │h │ > > > ├───┼───┼─────┤ > > │3 4│0 > 1│_3 4│ > > > │1 4│0 1│ 2 _3│ > > > └───┴───┴─────� > > > > 2 > > 3$'$f';'$g';'$h';($f);($g);($h) > > > ┌───┬───┬───┐ > > > │$f │$g │$h │ > > > ├───┼───┼───┤ > > > │2 2│2 2│2 2│ > > > └───┴───┴───� > > > > > f . g > > |domain error > > | f .g > > > f . h > > > |domain error > > > | f .h > > > %. f > > > 0.5 _0.5 > > > _0.125 0.375 > > > %. g > > > |domain error > > > | %.g > > > f . f > > > |domain error > > > | f .f > > > g . g > > > |domain error > > > | g .g > > > h . h > > > |domain error > > > | h .h > > > h * g > > > 0 4 > > 0 _3 > > > > The answers to all these operations > are > simple. Why are there these > domain errors? > > > > > Here is > the project this issue is > coming up > in, although I believe it is more basic than > that. > > > > vw=: (,@(((_ > > 0)$~#@:])`]`(1$~,~@:}.@[)}))`(,@(((0 > > _)$~#@:])`]`(1$~,~@:}.@[)}))@.(<:@:({.@:[)) > > NB. (moveno > boardsz) nw (x1 y1;x2 y2; > etc) > > NB. Gives vertex > weights for a move > list of a game as > vector > > > > > > rows=: > <"1 > > > cols=: <"1@|: > > diags=: (-.((0;_1)&{))@:(</.) > > a=: 3 3$i.9 > > conndiag=: 3 : > 0 > > r=. r,.|.r > [ r=. >0{z [ i=. 1 [ n=. > #z [ z=. > diags y > > > while. i<n do. > > > col=. >i{z > > i=. >:i > > adds=. >2<\col > > r=. r,adds > > r=. > r,(|."1 adds) > > end. > > ) > > > > connsq=: 3 : 0 > > > r=. 0 2$0 [ i=. 0 [ n=. #z [ z=. > y > > while. > i<n do. > > > col=. >i{z > > > i=. >:i > > > adds=. >2<\col > > r=. r,adds > > r=. r,(|."1 adds) > > end. > > ) > > connrws=: > connsq@:rows > > > conncls=: connsq@:cols > > > > edges=: 3 : 0 NB. takes board size as > argument > > a=. (2$bs)$i.*:bs [ bs=. y > > r=. (connrws > a),(conncls a) > > > r=. r,conndiag a > > > ) > > > > BOARDSZ=: 5 > > EDGES=: edges > BOARDSZ > > > VPLYR=: 1 NB. whether the vertical > > player is first or second > > > > A=: 3 : 0 > > mat=. (2$netsz)$0 [ netsz=. *:BOARDSZ > [ > j=. 0 [ i=. 0 [ arg=. > VPLYR,BOARDSZ > > > if. y=0 do. pcweights=. netsz$1 else. > > pcweights=. arg vw y end. > > while. i<netsz do. > > > while. j<netsz do. > > mat=. > (0{(((1 2$i,j) e. > > EDGES)*((i{pcweights)+j{pcweights))) (<(i,j))}mat > > j=. > >:j > > > end. > > j=. > 0 > > i=. > >:i > > > end. > > mat > > ) > > > > D=: 3 : 0 > > diag=. > (2$netsz)$ 1,netsz#0 [ netsz=. > > *:BOARDSZ [ colsum=. +/y [ y=. A y > > > r=. colsum*diag > > > ) > > > > L=: D-A NB. the Laplacian > > > > NB. %. Laplacian > or . (dot) for > Laplacian times > current gives domain errors > > > > Any help would be appreciated! > > > > Brian Babiak > MD > > https://drbabiak.com > > > > > > > > ---------------------------------------------------------------------- > 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