NB. J utilities to creat Rayleigh Quotient iteration images NB. Cliff Reiter, December 21, 2003 NB. Requires the use of the Lapack and Image3 addons NB. run script; then sample usage is: NB. NB. m2 200 mk_rqi_pic _15 9j12 NB. create an image NB. m7 200 mk_rqi_pic _8 8j8 NB. (m7 200 mk_rqi _8 8j8) write_image 'm7.png' NB. save an image NB. NB. updated for J601 Nov 8, 2007 require '~addons\math\lapack\lapack.ijs' require '~addons\math\lapack\dgeev.ijs' require '~addons\media\image3\view_m.ijs' NB. matrix product mp=: +/ . * NB. Euclidean length elen=: (+/&.:*:)@:| NB. unit vector unit=: % elen NB. A rqik u0 rqik=: 4 : 0"_ 0 A=. x k=. 0 x=. x0_rqi u=. y while. k<18 do. try. x=. unit x %. A - u *I_rqi catch. break. end. u=. (+x) mp A mp x k=. k+1 end. ({./:|EV_rqi-u),k ) NB. Makes an image mk_rqi=: 1 : 0 : EV_rqi=: (/: |.&.|:@:*.)>1{dgeev_jlapack_ x I_rqi=: =i.#x x0_rqi=: (#x)#1 z=. m zl_clur0 y 'b k'=. 0 1|:x rqik z p18p;(18~:k)*1+(3|k)+3*b ) NB. show the image mk_rqi_pic=: 1 : 0 : view_image x m mk_rqi y ) NB. get a mesh in the complex plane zl_clur0=: 4 : 0 w=. -~/9 o.y h=. -~/11 o.y (j.h-(+:h)*(i.%<:)x) +/ ({.y)+w*(i.%<:)x ) NB. a palette $p18p=: 0,,/0 2|:<.1 0.72 0.45 */ Hue 0 2 4 1 3 5{(i.%])6 NB. the matrices from [1] m1=: 1|.=i.5 m2=: 1,9 8 4,:5 6 2 m3=: 1,9j1 8 4,:5 6 2 m4=: ".;.(_2) 0 : 0 _6 7 1 _5 7 10 3 _4 1 3 8 _7 _5 _4 _7 8 ) m5=: ".;.(_2) 0 : 0 2 _9j_6 _3j_1 _3j1 _9j6 _4 _9j2 _6j6 _3j1 _9j_2 2 0j_3 _3j_1 _6j_6 0j3 2 ) m6=: ".;.(_2) 0 : 0 1 2 3 4 2 5 6 7 3 6 8 9 4 7 9 10 ) m7=: ".;.(_2) 0 : 0 _14 _10 0 0 21 12 _4 0 _15 _9 3 1 10 7 _1 _1 ) m8=: ".;.(_2) 0 : 0 12 _14 0 4 7 _6 _2 4 _10 12 _2 _2 _18 22 _3 _4 ) NB. [1] C. Reiter, Visualizing the Dynamics of the Rayleigh Quotient Iteration, Computers & Graphics, 16 3 (1992) 341-344