NB: I wrote this before looking at the J rosetta code which includes the generation of the Ybus matrix -

If your object is to find the equivalent (Thevenin) resistance between any 2 points, I would suggest that it would be better to work with an admittance model. It eliminates the problem with infinite resistances and typically will be sparse. Such methods are used with large power grid analyses. The admittance method chooses one node as reference, typically ground but it could be an other node or vertex. (it doesn't appear in the admittance "Ybus" matrix) . This is easy to set up as the yii terms are simply the sum of admittance's connected to that node and the yij terms are the negative of the sum of admittance's connected between i and j. Inverting this matrix gives what is called the Zbus matrix where the diagonal elements are the Thevenin impedance's (resistive in your case). Off diagonal elements are transfer impedance's. Essentially a two port idea extended to multi ports. A neat thing is that Zbus includes all the information to find the resistance seen between any pair of nodes as once you have it, the resistance between nodes i and j is Rii +Rjj- 2Rij

here is a simple sample:

node      node      r            1/r
  0            ref        0.25        4
  0            1            0.5         2
  0            2            0.2         5
  1            ref         0.25       4
  1             2           0.25       4

]ybus=: > 11 _2 _5;_2 10 _4;_5 _4 9

11 _2 _5

_2 10 _4

_5 _4 9

]zbus=:%.ybus

0.165179 0.0848214 0.129464

0.0848214 0.165179 0.120536

0.129464 0.120536 0.236607

+/ 1 _1 _1 1*,(<0 2;0 2){zbus

0.142857 impedance seen between node 0 and node 2

etc

Don

On 2017-07-27 9:56 PM, Brian Babiak wrote:
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

Reply via email to