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
