(* num-ht3/Prog/a26-ppp.m, 2021.05.18-- *)
(*Orig: num-ht-2019/a26-ppp.m, 2017.06.06 *)

hh[k_,n_,x_,y_]:=NIntegrate[t^k*Exp[-t]*HypergeometricPFQ[{},{n},t*y],{t,0,x}];

testTiming1[]:=Timing[Table[hh[10,1,1,y],{y,{1,2,3,4,5,6,7,8,9,10}}]]
(* 0.1049s on orange3m/Mathematica 11.2 (9 points)
   python3 hkn-10-1.py  takes 0.217s (900 points)
*)

test20191029a[]:=InputForm[Table[{y,N[Log[hh[10,1,1,y]]]},{y,1,1000,10}]]

(* Hmat[I][J] = ( *Hkn_p )(T-II,T-S+1,X,Y[J]|option_list=getopt()); 
   in a26-init.rr *)
matHH[s_,t_,x_,y_]:=Block[{Mat,i,j},
  Mat=IdentityMatrix[s];
  For[i=1,i<=s,i++,
    For[j=1,j<=s,j++,
      Mat[[i,j]] = hh[t-i,t-s+1,x,y[[j]]];
    ]
  ];
  Return[Mat];
]

vandermond[y_]:=Block[{n,Mat,i,j},
  n=Length[y];
  Mat=IdentityMatrix[n];
  For[i=1,i<=n,i++,
    For[j=1,j<=n,j++,
      Mat[[i,j]] = y[[i]]^(n-j);
    ]
  ];
  Return[Mat];
]

a26Error=0;
ppp[s_,t_,x_,y_]:=Block[{mat,d,f},
  f=Exp[-Sum[y[[i]],{i,1,Length[y]}]]/Det[vandermond[y]];
  f=f/((t-s)!)^s;
  mat = matHH[s,t,x,y];
  Print[mat];  a26Error={{"hmat",mat}};
  Print["Det=",Det[mat]];
  Return[N[Det[mat]*f]];
]

(* 2017.06.11, cf. a26-init.rr  init_f3()
   Test. initf3[9,1,3,1/2]
init_f3(9,1,3,1/2);
[[1439933.27389483,0.200476085708227,-0.108662520939692],...]
*)
initf3[k_,n_,t_,y_]:=
  Block[{x,z,h,h2,h2diff,init},
     x=t^2; z=y^(1/2);
     h=hh[k,n,x,y];
     h2=N[HypergeometricPFQ[{},{n},x*y]];
     h2diff=N[D[HypergeometricPFQ[{},{n},xx],xx] /. {xx->(x*y)}];
     init={h,N[Exp[-2*z*t]*h2],
     N[t*(-2*z*Exp[-2*z*t]*h2+Exp[-2*z*t]*2*t*y*h2diff)]};
     Return[init]
]
(* 2017.07.13 *)
hhwithPrec[k_,n_,x_,y_,prec_]:=N[Integrate[t^k*Exp[-t]*HypergeometricPFQ[{},{n},t*y],{t,0,x}],prec];
initf3withPrec[k_,n_,t_,y_,prec_]:=
  Block[{x,z,h,h2,h2diff,init},
     x=t^2; z=y^(1/2);
     h=hhwithPrec[k,n,x,y,prec];
     h2=N[HypergeometricPFQ[{},{n},x*y],prec];
     h2diff=N[D[HypergeometricPFQ[{},{n},xx],xx] /. {xx->(x*y)},prec];
     init={h,N[Exp[-2*z*t]*h2,prec],
     N[t*(-2*z*Exp[-2*z*t]*h2+Exp[-2*z*t]*2*t*y*h2diff),prec]};
     Return[init]
]
logInitf3withPrec[k_,n_,t_,y_,prec_]:=
  Block[{x,z,h,h2,h2diff,init},
     x=t^2; z=y^(1/2);
     h=N[Log[hhwithPrec[k,n,x,y,prec]],prec];
     h2=N[HypergeometricPFQ[{},{n},x*y],prec];
     h2diff=N[D[HypergeometricPFQ[{},{n},xx],xx] /. {xx->(x*y)},prec];
     init={h,N[Exp[-2*z*t]*h2,prec],
     N[t*(-2*z*Exp[-2*z*t]*h2+Exp[-2*z*t]*2*t*y*h2diff),prec]};
     Return[init]
]

matHHwithPrec[s_,t_,x_,y_,prec_]:=Block[{Mat,i,j},
  Mat=IdentityMatrix[s];
  For[i=1,i<=s,i++,
    For[j=1,j<=s,j++,
      Mat[[i,j]] = hhwithPrec[t-i,t-s+1,x,y[[j]],prec];
    ]
  ];
  Return[Mat];
]
pppWithPrec[s_,t_,x_,y_,prec_]:=Block[{mat,d,f},
  f=Exp[-Sum[y[[i]],{i,1,Length[y]}]]/Det[vandermond[y]];
  f=f/((t-s)!)^s;
  mat = matHHwithPrec[s,t,x,y,prec];
  Print[mat];  a26Error={{"hmat",mat},f};
  Print["Det=",Det[mat]];
  Return[N[Det[mat]*N[f,prec],prec]];
]

(* 2023.07.13 *)
initOfGaugedHkn[yval_]:=Block[{hh1,myrule,G,init1,y},
  hh1 = Integrate[t^k*Exp[-t]*HypergeometricPFQ[{},{n},t*y],{t,0,x}] /. {k->10,n->1, x->1};
  (* specify k = 10, n = 1, x = 1 *)
  myrule = {n -> 1, k -> 10, x -> 1};
  G = Exp[-y]*y^(-(k - n + 1));
  init1=Table[N[D[hh1, {y, i}] /. y -> yval, $MachinePrecision], {i, 0, 3}] ;
  Return[N[G /. myrule /. {y -> yval}]*init1]
]

(* Table[initOfGaugedHkn[y], {y, 1, 10}] *)

initOfHkn[yval_] := 
 Block[{hh1, myrule, G, init1, y}, 
  hh1 = Integrate[
     t^k*Exp[-t]*HypergeometricPFQ[{}, {n}, t*y], {t, 0, x}] /. {k -> 
      10, n -> 1, x -> 1};
  (*specify k=10,n=1,x=1*)
  
  init1 = 
   Table[N[D[hh1, {y, i}] /. y -> yval, $MachinePrecision], {i, 0, 3}];
  Return[init1]]

