Copyright © Masayuki Noro 2016. All rights reserved.
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
1 n_wishartd.rr | ||
Index |
• Restrction of matrix 1F1 on diagonal regions | ||
1.2 Numerical comptation of restricted function | ||
1.3 Differential operators with partial fraction coefficients | ||
1.4 Experimental implementation of Runge-Kutta methods |
This manual explains about ‘n_wishartd.rr’, a package for computing a system of differential equations satisfied by the matrix 1F1 on a diagonal region. To use this package one has to load ‘n_wishartd.rr’.
[...] load("n_wishartd.rr");
A prefix n_wishartd.
is necessary to call the functions in this package.
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
• n_wishartd.diagpf | ||
1.1.2 n_wishartd.message |
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
n_wishartd.diagpf
computes a system of PDEs satisfied by the m variate matrix 1F1 on a diagonal region specified by blocks.
A list [E1,E2,...] where each Ei is a differential operator with partial fraction coefficients and it annihilates the restricted 1F1.
A natural number
A list [[s1,e1],[s2,e2],...].
See below.
[2649] Z=n_wishartd.diagpf(5,[[1,3],[4,5]]); [ [[[[-1,[]]],(1)*<<0,0,0,0,3,0>>], [[[-2,[[y1-y4,1]]],[-2,[[y4,1]]]],(1)*<<0,1,0,0,1,0>>], [[[9/2,[[y1-y4,1]]],[-3*c+11/2,[[y4,1]]],[3,[]]],(1)*<<0,0,0,0,2,0>>], ... [[[-6*a,[[y1-y4,1],[y4,1]]],[(4*c-10)*a,[[y4,2]]],[-4*a,[[y4,1]]]], (1)*<<0,0,0,0,0,0>>]], [[[[-1,[]]],(1)*<<0,4,0,0,0,0>>], [[[-6,[[y1-y4,1]]],[-6*c+10,[[y1,1]]],[6,[]]],(1)*<<0,3,0,0,0,0>>], [[[5,[[y1-y4,1]]],[-5,[[y1,1]]]],(1)*<<0,2,0,0,1,0>>], ... [[[21*a,[[y1-y4,2],[y1,1]]],[(36*c-87)*a,[[y1-y4,1],[y1,2]]], [-36*a,[[y1-y4,1],[y1,1]]],[(18*c^2-84*c+96)*a,[[y1,3]]], [-9*a^2+(-36*c+87)*a,[[y1,2]]],[18*a,[[y1,1]]]],(1)*<<0,0,0,0,0,0>>]] ]
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
n_wishartd.message
starts/stops displaying messages during computation.
Start displaying messages if onoff=1. Stop displaying messages if onoff=0.
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
1.2.1 n_wishartd.prob_by_hgm | ||
1.2.2 n_wishartd.prob_by_ps | ||
1.2.3 n_wishartd.ps |
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
n_wishartd.prob_by_hgm
computes the value of the distribution function of the largest eigenvalue of a Wishart matrix.
The number of variables.
The degrees of freedom.
A list of the multiplicities of repeated eigenvalues.
A list of repeated eigenvalues.
n_wishartd.prob_by_hgm
computes the value of the distribution function by using HGM
for a covariance matrix which has repeated eigenvalues si with multiplicity pi (i=1,...,k).
[...] n_wishartd.message(1)$ [...] P=n_wishartd.prob_by_hgm(10,100,[9,1],[1/100,1],100|eps=10^(-6)); ... [x0=,8/25] Step=10000 [0] [8.23700622458446e-17,8.23700622459772e-17] ... Step=1280000 [0][100000][200000][300000]...[900000][1000000][1100000][1200000] [0.516246820120598,0.516246820227214] [log ratio=,4.84611265040128] Step=2560000 [0][100000][200000][300000]...[2200000][2300000][2400000][2500000] [0.516246912003845,0.516246912217004] [log ratio=,4.93705929488356] [diag,18.6292,pfaffian,1.09207,ps,41.0026,rk,213.929] 0.516246912217004 266.4sec + gc : 8.277sec(276.8sec)
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
n_wishartd.prob_by_ps
computes the value of the distribution function of the largest eigenvalue of a Wishart matrix.
The number of variables.
The degrees of freedom.
A list of the multiplicities of repeated eigenvalues.
A list of repeated eigenvalues.
[...] Q=n_wishartd.prob_by_ps(10,100,[9,1],[1/100,1],1/2); ... [I=,109,act,24.9016,actmul,0,gr,19.7852] 2.69026137621748e-165 61.69sec + gc : 2.06sec(64.23sec) [...] R=n_wishartd.prob_by_hgm(10,100,[9,1],[1/100,1],1/2|td=50); [diag,15.957,pfaffian,1.00006,ps,5.92437,rk,1.29208] 2.69026135182769e-165 23.07sec + gc : 1.136sec(24.25sec)
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
n_wishartd.ps
computes a truncated power series solution up to the total degree td.
A list of polynomial
A list of differential operators with partial fraction coefficients.
A list of variables.
[...] Z=n_wishartd.diagpf(10,[[1,5],[6,10]])$ [...] Z0=subst(Z,a,(10+1)/2,c,(10+100+1)/2)$ [...] PS=n_wishartd.ps(Z0,[y1,y6],10)$ [...] PS[0]; 197230789502743383953639/230438384724900975787223158176000*y1^10+ ... +(6738842542131976871672233/1011953706634779427957034268904320*y6^9 ...+3932525/62890602*y6^2+1025/4181*y6+55/111)*y1 +197230789502743383953639/230438384724900975787223158176000*y6^10 +...+1395815/62890602*y6^3+3175/25086*y6^2+55/111*y6+1
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
The coefficients of the PDE satisfied by the matrix 1F1 are written as a sum of 1/yi and 1/(yi-yj) multiplied by constants. Furthermore the result of diagonalization by l’Hopital rule can also be written as a sum of partial fractions.
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
By using partial fractions explained in the previous section, differential operators with partial fraction coefficients are represented. Let f1,...,fk be partial fractions and d1,...,dk distributed monomials such that d1>...>dk) with respected to the current monomial ordering. Then a differential operator f1*d1+...+fk*dk is represented as a list [f1,d1],...[fk,dk]].
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
1.3.3.1 n_wishartd.wsetup | ||
1.3.3.2 n_wishartd.addpf | ||
1.3.3.3 n_wishartd.mulcpf | ||
1.3.3.4 n_wishartd.mulpf | ||
1.3.3.5 n_wishartd.muldpf |
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
n_wishartd.wsetup
A natural number.
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
n_wishartd.addpf
A differential operator with partial fraction coefficients.
Differential operators with partial fraction coefficients.
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
n_wishartd.mulcpf
A differential operator with partial fraction coefficients.
A partial fraction.
Differential operators with partial fraction coefficients.
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
n_wishartd.mulpf
A differential operator with partial fraction coefficients.
Differential operators with partial fraction coefficients.
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
n_wishartd.muldpf
A differential operator with partial fraction coefficients.
A variable.
A differential operator with partial fraction coefficients.
[...] n_wishartd.wsetup(4)$ [...] P=n_wishartd.wishartpf(4,1); [[[[1,[]]],(1)*<<0,2,0,0,0>>],[[[1/2,[[y1-y2,1]]],[1/2,[[y1-y3,1]]], ...,[[[-a,[[y1,1]]]],(1)*<<0,0,0,0,0>>]] [...] Q=n_wishartd.muldpf(y1,P); [[[[1,[]]],(1)*<<0,3,0,0,0>>],[[[1/2,[[y1-y2,1]]],[1/2,[[y1-y3,1]]], ...,[[[a,[[y1,2]]]],(1)*<<0,0,0,0,0>>]]
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
1.4.1 rk_ratmat |
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
In the function n_wishartd.ps_by_hgm
, after computing the Pfaffian matrices for
the sytem of PDEs on a diagonal region, it executes a built-in function
rk_ratmat
which computes an approximate solution of the Pfaffian system
by Runge-Kutta method for a spcified step size.
This function is repeated until the result gets stabilized, by doubling the step size.
rk_ratmat
can be used as a general-purpose Runge-Kutta driver and we explain how to use it.
rk_ratmat
solves a system of linear ODEs with rational function coefficients.
A list of real numbers.
4 or 5.
An array of constant matrices.
A polynomial.
Real numbers.
A natural number.
A real vector.
rk_ratmat
solves an initial value problem
dF/dx = P(x)F, F(x0)=f0 for P(x)=1/den(num[0]+num[1]x+...+num[k-1]x^(k-1)) by a Runge-Kutta method.
rk_ratmat
we multiply the result for the normalized f0 and the 0-th component
of the original f0 to get the desired result.
[...] F=ltov([sin(1/x),cos(1/x),sin(1/x^2),cos(1/x^2)]); [ sin((1)/(x)) cos((1)/(x)) sin((1)/(x^2)) cos((1)/(x^2)) ] [...] F0=map(eval,map(subst,F,x,1/10)); [ -0.54402111088937 -0.839071529076452 -0.506365641109759 0.862318872287684 ] [...] N0=matrix(4,4,[[0,0,0,0],[0,0,0,0],[0,0,0,-2],[0,0,2,0]])$ [...] N1=matrix(4,4,[[0,-1,0,0],[1,0,0,0],[0,0,0,0],[0,0,0,0]])$ [...] N=ltov([N0,N1])$ [...] D=x^3$ [...] R=rk_ratmat(5,N,D,1/10,10,10^4,F0)$ [...] for(T=R,A=1;T!=[];T=cdr(T))A *=car(T)[1]; [...] A; 0.0998334 [...] F1=map(eval,map(subst,F,x,10)); [ 0.0998334166468282 0.995004165278026 0.00999983333416666 0.999950000416665 ]
[ << ] | [ < ] | [ Up ] | [ > ] | [ >> ] | [Top] | [Contents] | [Index] | [ ? ] |
Jump to: | N R |
---|
Jump to: | N R |
---|
[Top] | [Contents] | [Index] | [ ? ] |
[Top] | [Contents] | [Index] | [ ? ] |
[Top] | [Contents] | [Index] | [ ? ] |
This document was generated on September 25, 2018 using texi2html 5.0.
The buttons in the navigation panels have the following meaning:
Button | Name | Go to | From 1.2.3 go to |
---|---|---|---|
[ << ] | FastBack | Beginning of this chapter or previous chapter | 1 |
[ < ] | Back | Previous section in reading order | 1.2.2 |
[ Up ] | Up | Up section | 1.2 |
[ > ] | Forward | Next section in reading order | 1.2.4 |
[ >> ] | FastForward | Next chapter | 2 |
[Top] | Top | Cover (top) of document | |
[Contents] | Contents | Table of contents | |
[Index] | Index | Index | |
[ ? ] | About | About (help) |
where the Example assumes that the current position is at Subsubsection One-Two-Three of a document of the following structure:
This document was generated on September 25, 2018 using texi2html 5.0.