[VA] Short & Sweet Math Challenge #24: "2019 Spring Special 5-tier"

[Oulan]

Extending the approach of JF, you can use the following HP prime program
--------------------
#pragma mode( separator(., integer(h32) )
#cas
reduc(l,n):= BEGIN
LOCAL m;
m:=MAKELIST(0,z,0,n);
FOR Z FROM 0 TO n DO
m[n+1-Z]:=l[SIZE(l)-Z];
END;
RETURN m;
END;
fff():= BEGIN
LOCAL l,m,n,z,p,q,r;
LOCAL a15,a13,a11,a9,a7,a5,a3;
PURGE(a15,a13,a11,a9,a7,a5,a3);
z:={a15,0,a13,0,a11,0,a9,0,a7,0,a5,0,a3,0,1,0};
n:=15;
p:=poly2symb(z,x);
r:=poly2symb(z,y);
FOR Y FROM 1 TO 2 DO
q:=(p|x=r);
l:=symb2poly(q,y);
l:=reduc(l,n);
p:=poly2symb(l,x);
END;
m:=MAKELIST(0,Z,0,n);
FOR Z FROM 0 to n DO
IF (Z MOD 2) = 1 THEN
m[n+1-Z]:=((-1)^FLOOR(Z/2))/(Z!);
END;
END;
PRINT(l);
PRINT(m);
a3:=eval(solve(l[n-2]=m[n-2],a3)[1]);
// was a3:=solve(l[n-2]=m[n-2],a3);
a5:=eval(solve(l[n-4]=m[n-4],a5)[1]);
a7:=eval(solve(l[n-6]=m[n-6],a7)[1]);
a9:=eval(solve(l[n-8]=m[n-8],a9)[1]);
a11:=eval(solve(l[n-10]=m[n-10],a11)[1]);
a13:=eval(solve(l[n-12]=m[n-12],a13)[1]);
a15:=eval(solve(l[n-14]=m[n-14],a15)[1]);
RETURN [a3,a5,a7,a9,a11,a13,a15];
END;
#end
--------------------
to compute the coefficients of the taylor serie of CIN(x).
Use it in CAS mode : " k:=fff() 'enter' ", then " k 'enter' "

Be careful this program will take some times on a real Prime.

But the converging is very slow, coefficient up to x^15 follows but give only 10^-6 error near 1

-1/18 -7/1080 -643/408240 -13583/29393280 -29957/215550720 -24277937/648499737600 -6382646731/953294614272000

Btw can someone explain the warning displayed when solving for the coefficients ?

"Warning, ^ is ambiguous on non square matrices. Use .^ to apply ^ element by element."

I don't see any matrices solving here ... ok I saw the problem see listing. Sometimes list of list are not displayed with all brackets.

Anyway, there should be a better approach to solve this nice challenge

EDIT new version of program, avoid computing useless power

#pragma mode( separator(., integer(h32) )
#cas
mulT(a,b):= BEGIN
LOCAL n,p,j,k;
n:=SIZE(a);
p:=MAKELIST(0,j,1,n);
FOR j FROM 1 TO n DO
FOR k FROM 1 TO n+1-j DO
p[j+k-1]+=a[j]*b[k]; // test removed Thanks Albert
END;
END;
RETURN simplify(p);
END;

fff2(n):= BEGIN
LOCAL cin,ccin,si;
LOCAL p,q,xn,s;
LOCAL lvar,lexpr;
LOCAL a3,a5,a7,a9,a11,a13,a15,a17;
LOCAL a19,a21,a23,a25,a27,a29,a31,a33;
LOCAL a35,a37,a39,a41,a43,a45,a47,a49;
LOCAL vars;
PURGE(a3,a5,a7,a9,a11,a13,a15,a17);
PURGE(a19,a21,a23,a25,a27,a29,a31,a33);
PURGE(a35,a37,a39,a41,a43,a45,a47,a49);
PURGE(x,y);
vars:={1,a3,a5,a7,a9,a11,a13,a15,a17,a19,a21,a23,a25,a27,a29,a31,a33,a35,a37,a39​,a41,a43,a45,a47,a49};
cin:=MAKELIST(IFTE(p MOD 2,vars[(p+1)/2],0),p,0,n);
ccin:=cin;
FOR q FROM 1 TO 2 DO
s:=MAKELIST(0,p,0,n);
s[1]:=ccin[1];
xn:=cin;
FOR p FROM 2 TO n+1 DO
s:=s+ccin[p]*xn;
IF p<=n THEN xn:=mulT(xn,cin);END;
END;
ccin:=s;
END;
si:=MAKELIST(IFTE(p MOD 2,((-1)^FLOOR(p/2))/(p!),0),p,0,n);
lvar:=MAKELIST(cin[p],p,4,n+1,2);
lexpr:=MAKELIST(ccin[p]=si[p],p,4,n+1,2);
s:=solve(lexpr,lvar);
RETURN s[1];
END;
#end

use with ff2(2*n+1) n from 2 to 24 (fff2(29) start to be long on a real G2), fff2(49) take few seconds on a virtual one.