Skip to content

Commit

Permalink
Improve HypergeometricPFQ
Browse files Browse the repository at this point in the history
  • Loading branch information
axkr committed Jul 29, 2023
1 parent 271b0fe commit 36534ca
Show file tree
Hide file tree
Showing 7 changed files with 105 additions and 42 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -471,52 +471,20 @@ private static IExpr binaryD(final IExpr functionOfX, IExpr x, final IAST ast,
} else if (function.isTimes()) {
return function.map(F.PlusAlloc(16), new BinaryBindIth1st(function, F.D(S.Null, x)));
} else if (function.isPower()) {
// f ^ g
final IExpr f = function.base();
final IExpr g = function.exponent();
if (g.isFree(x)) {
// g*D(f,y)*f^(g-1)
return F.Times(g, F.D(f, x), F.Power(f, g.dec()));
}
if (f.isFree(x)) {
if (f.isE()) {
return F.Times(F.D(g, x), F.Exp(g));
}
// D(g,y)*Log(f)*f^g
return F.Times(F.D(g, x), F.Log(f), F.Power(f, g));
}

// D[f_^g_,y_]:= f^g*(((g*D[f,y])/f)+Log[f]*D[g,y])
final IASTAppendable resultList = F.TimesAlloc(2);
resultList.append(F.Power(f, g));
resultList
.append(F.Plus(F.Times(g, F.D(f, x), F.Power(f, F.CN1)), F.Times(F.Log(f), F.D(g, x))));
return resultList;
return power(function, x);
} else if (function.isAST(S.Surd, 3)) {
// Surd[f,g]
final IExpr f = function.base();

if (function.exponent().isInteger()) {
final IInteger g = (IInteger) function.exponent();
if (g.isMinusOne()) {
return F.Times(F.CN1, F.D(f, x), F.Power(f, F.CN2));
}
final IRational gInverse = g.inverse();
if (g.isNegative()) {
if (g.isEven()) {
return F.Times(gInverse, F.D(f, x), F.Power(F.Surd(f, g.negate()), g.dec()));
}
return F.Times(gInverse, F.D(f, x), F.Power(f, F.CN1),
F.Power(F.Surd(f, g.negate()), F.CN1));
}
return F.Times(gInverse, F.D(f, x), F.Power(F.Surd(f, g), g.dec().negate()));
}
return surd(function, x);
} else if ((header == S.Log) && (function.isAST2())) {
if (function.isFreeAt(1, x)) {
// D[Log[i_FreeQ(x), x_], z_]:= (x*Log[a])^(-1)*D[x,z];
return F.Times(F.Power(F.Times(function.arg2(), F.Log(function.arg1())), F.CN1),
F.D(function.arg2(), x));
}
} else if (function.isAST(S.HypergeometricPFQ, 4)//
&& function.first().isList()//
&& function.second().isList()) {
return hypergeometricPFQ(function, x);
// } else if (header == F.LaplaceTransform && (listArg1.size()
// == 4)) {
// if (listArg1.arg3().equals(x) && listArg1.arg1().isFree(x,
Expand Down Expand Up @@ -566,6 +534,76 @@ private static IExpr binaryD(final IExpr functionOfX, IExpr x, final IAST ast,
return F.NIL;
}

private static IExpr power(final IAST function, IExpr x) {
// f ^ g
final IExpr f = function.base();
final IExpr g = function.exponent();
if (g.isFree(x)) {
// g*D(f,y)*f^(g-1)
return F.Times(g, F.D(f, x), F.Power(f, g.dec()));
}
if (f.isFree(x)) {
if (f.isE()) {
return F.Times(F.D(g, x), F.Exp(g));
}
// D(g,y)*Log(f)*f^g
return F.Times(F.D(g, x), F.Log(f), F.Power(f, g));
}

// D[f_^g_,y_]:= f^g*(((g*D[f,y])/f)+Log[f]*D[g,y])
final IASTAppendable resultList = F.TimesAlloc(2);
resultList.append(F.Power(f, g));
resultList
.append(F.Plus(F.Times(g, F.D(f, x), F.Power(f, F.CN1)), F.Times(F.Log(f), F.D(g, x))));
return resultList;
}

private static IExpr surd(final IAST function, IExpr x) {
final IExpr f = function.base();
if (function.exponent().isInteger()) {
final IInteger g = (IInteger) function.exponent();
if (g.isMinusOne()) {
return F.Times(F.CN1, F.D(f, x), F.Power(f, F.CN2));
}
final IRational gInverse = g.inverse();
if (g.isNegative()) {
if (g.isEven()) {
return F.Times(gInverse, F.D(f, x), F.Power(F.Surd(f, g.negate()), g.dec()));
}
return F.Times(gInverse, F.D(f, x), F.Power(f, F.CN1),
F.Power(F.Surd(f, g.negate()), F.CN1));
}
return F.Times(gInverse, F.D(f, x), F.Power(F.Surd(f, g), g.dec().negate()));
}
return F.NIL;
}

private static IExpr hypergeometricPFQ(final IAST function, IExpr x) {
IAST list1 = (IAST) function.first();
IAST list2 = (IAST) function.second();
if (list1.isFree(x)&&list2.isFree(x)) {
IExpr arg3 = function.arg3();
if (list1.isEmpty() && list2.isEmpty()) {
return F.Times(F.Exp(arg3), F.D(arg3, x));
}
IExpr timesNumerator = list1.argSize() == 0 ? F.C1 : list1.apply(S.Times, 1);
IExpr timesDenominator = list2.argSize() == 0 ? F.C1 : list2.apply(S.Times, 1);
IASTAppendable newList1 = F.ListAlloc(list1.argSize());

for (int i = 1; i < list1.size(); i++) {
newList1.append(F.Plus(F.C1, list1.get(i)));
}
IASTAppendable newList2 = F.ListAlloc(list2.argSize());
for (int i = 1; i < list2.size(); i++) {
newList2.append(F.Plus(F.C1, list2.get(i)));
}
return F.Times(timesNumerator, F.Power(timesDenominator, F.CN1), //
F.HypergeometricPFQ(newList1, newList2, arg3), //
F.D(arg3, x));
}
return F.NIL;
}

private static IExpr dPiecewise(int[] dim, final IAST piecewiseFunction, final IAST ast,
EvalEngine engine) {

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -164,6 +164,9 @@ public interface FunctionExpandRules {
// HypergeometricPFQ({a_},{b_},z_):=Hypergeometric1F1(a,b,z)
SetDelayed(HypergeometricPFQ(list(a_),list(b_),z_),
Hypergeometric1F1(a,b,z)),
// HypergeometricPFQ({a_,b_},{c_},z_):=Hypergeometric2F1(a,b,c,z)
SetDelayed(HypergeometricPFQ(list(a_,b_),list(c_),z_),
Hypergeometric2F1(a,b,c,z)),
// Hypergeometric2F1(2,b_,c_,-1/2):=1/3*(3-b)/;5/2-b/2==Expand(c)
SetDelayed(Hypergeometric2F1(C2,b_,c_,CN1D2),
Condition(Times(C1D3,Subtract(C3,b)),Equal(Plus(QQ(5L,2L),Times(CN1D2,b)),Expand(c)))),
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -15,10 +15,19 @@ public interface HypergeometricPFQRules {
* <li>index 0 - number of equal rules in <code>RULES</code></li>
* </ul>
*/
final public static int[] SIZES = { 0, 2 };
final public static int[] SIZES = { 0, 5 };

final public static IAST RULES = List(
IInit(HypergeometricPFQ, SIZES),
// HypergeometricPFQ({},{},z_):=E^z
ISetDelayed(HypergeometricPFQ(List(),List(),z_),
Exp(z)),
// HypergeometricPFQ({a_},{},z_):=(1-z)^(-a)
ISetDelayed(HypergeometricPFQ(list(a_),List(),z_),
Power(Subtract(C1,z),Negate(a))),
// HypergeometricPFQ({},{b_},z_):=z^(1/2-b/2)*BesselI(-1+b,2*Sqrt(z))*Gamma(b)
ISetDelayed(HypergeometricPFQ(List(),list(b_),z_),
Times(Power(z,Plus(C1D2,Times(CN1D2,b))),BesselI(Plus(CN1,b),Times(C2,Sqrt(z))),Gamma(b))),
// HypergeometricPFQ({a_,b_},{c_,b_},z_):=HypergeometricPFQ({a},{c},z)
ISetDelayed(HypergeometricPFQ(list(a_,b_),list(c_,b_),z_),
HypergeometricPFQ(list(a),list(c),z)),
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -4658,6 +4658,15 @@ public void testCyclotomic() {
}

public void testD() {
check("D(HypergeometricPFQ({}, {}, x), x)", //
"E^x");
check("D(HypergeometricPFQ({a}, {}, x), x)", //
"a/(1-x)^(1+a)");
check("D(HypergeometricPFQ({a,b}, {u,v,w}, x), x)", //
"(a*b*HypergeometricPFQ({1+a,1+b},{1+u,1+v,1+w},x))/(u*v*w)");
check("D(HypergeometricPFQ({a, b,c,d,e,f}, {u,v,w,y,z}, x), x)", //
"(a*b*c*d*e*f*HypergeometricPFQ({1+a,1+b,1+c,1+d,1+e,1+f},{1+u,1+v,1+w,1+y,1+z},x))/(u*v*w*y*z)");

check("D(HarmonicNumber(Sin(Cos(x)), y), x)", //
"y*Cos(Cos(x))*Sin(x)*(HarmonicNumber(Sin(Cos(x)),1+y)-Zeta(1+y))");
check("D(Sin(f(x)),{x,3})", //
Expand Down
2 changes: 1 addition & 1 deletion symja_android_library/rules/DRules.m
Original file line number Diff line number Diff line change
Expand Up @@ -47,7 +47,7 @@
/; FreeQ({f},x),
D(Zeta(f_,g_),x_?NotListQ):=(-f)*Zeta(1+f, g)*D(g,x)
/; FreeQ({f},x),

D(Hypergeometric2F1(a_, b_, c_, f_), x_?NotListQ) := (a*b*Hypergeometric2F1(1 + a, 1 + b, 1 + c, f)*D(f,x))/c
/; FreeQ({a,b,c},x),
D(Hypergeometric2F1(a_, b_, c_, x_), {x_,n_}) := Hypergeometric2F1(a + n, b + n, c + n, x)*(Pochhammer(a, n)*Pochhammer(b, n))/Pochhammer(c, n)
Expand Down
1 change: 1 addition & 0 deletions symja_android_library/rules/FunctionExpandRules.m
Original file line number Diff line number Diff line change
Expand Up @@ -73,6 +73,7 @@

HypergeometricPFQ({1/2}, {1, 1}, z_) := BesselI(0, Sqrt(z))^2,
HypergeometricPFQ({a_}, {b_}, z_) := Hypergeometric1F1(a,b,z),
HypergeometricPFQ({a_,b_}, {c_}, z_) := Hypergeometric2F1(a,b,c,z),
Hypergeometric2F1(2, b_, c_, -1/2) := (3-b)/3
/; (5/2 - 1/2*b)==Expand(c),
Hypergeometric2F1(a_, a_ + 1/2, c_, z_) := (2^(-1+2*a)*(1+Sqrt(1-z))^(1-2*a))/Sqrt(1-z)
Expand Down
7 changes: 5 additions & 2 deletions symja_android_library/rules/HypergeometricPFQRules.m
Original file line number Diff line number Diff line change
@@ -1,5 +1,8 @@
{
HypergeometricPFQ({a_,b_},{c _,b_},z_) := HypergeometricPFQ({a},{c},z),
{
HypergeometricPFQ[{},{},z_] := E^(z),
HypergeometricPFQ[{a_},{},z_] := (1-z)^(-a),
HypergeometricPFQ[{},{b_},z_] := z^(1/2-b/2)*BesselI(-1+b,2*Sqrt(z))*Gamma(b),
HypergeometricPFQ({a_,b_},{c_,b_},z_) := HypergeometricPFQ({a},{c},z),
HypergeometricPFQ({1/2, b_}, {3/2, c_}, z_) := (b/(2*b - 1))*(Sqrt(Pi/z)*Erfi(Sqrt(z)) - (Gamma(b) - Gamma(b,-z))/(-z)^b)
/; PossibleZeroQ(b+1-c)
}

0 comments on commit 36534ca

Please sign in to comment.