You are not logged in.
(* https://resources.wolframcloud.com/PacletRepository/resources/JanM/Dixon/ *)
<< "JanM`Dixon`";
Hold[{f, u, v, α, s, c, x, y, z, b, γ, w1, w2, λ, λ0, λ1, λ2, κ, κ0, κ1, κ2, ϕ, pr, sqrt},
pr = 50;
α = -RandomReal[WorkingPrecision -> pr];
b = (9*(1 - α + α^2))^(1/3);
γ = CubeRoot[1 + α^3];
sqrt[z_] := RandomChoice[{1, -1}]*Sqrt[z];
w1 = Exp[2*sqrt[-1]*Pi/3];
w2 = Conjugate[w1];
λ0 = Gamma[1/3]^3/(2*Sqrt[3]*Pi);
κ0 = 4*Pi^2/(3*Gamma[1/3]^3);
Assert[N[DixonF[λ0, 0] == κ0, pr]];
x = λ0*Hypergeometric2F1[1/3, 1/3, 2/3, -α^3];
y = κ0*α*Hypergeometric2F1[2/3, 2/3, 4/3, -α^3];
λ = Quiet@Re[x + y];
λ1 = w1*x + w2*y;
λ2 = w2*x + w1*y;
x = κ0*Hypergeometric2F1[-1/3, 2/3, 1/3, -α^3];
y = λ0*α^2*Hypergeometric2F1[1/3, 4/3, 5/3, -α^3]/2;
κ = Quiet@Re[x + y];
κ1 = α + w2*x + w1*y;
κ2 = α + w1*x + w2*y;
If[α < 0 && !VectorQ[{λ, κ}, NumericQ], With[{b = 1 + α^(-3)},
{λ, κ} = 2*Pi*3^(-3/2)*{
-Hypergeometric2F1[1/3, 2/3, 1, b]/α,
Hypergeometric2F1[2/3, 4/3, 1, b]/α^2
};
]];
κ += α;
{u, v} = RandomReal[{0, λ}, 2, WorkingPrecision -> pr];
x = DixonSM[0, α];
y = DixonCM[0, α];
Assert[x == 0];
Assert[y == 1];
x = DixonSM[λ, α];
y = DixonCM[λ, α];
Assert[x == 1];
Assert[y == 0];
x = DixonSM[λ/2, α];
y = DixonCM[λ/2, α];
Assert[x == y];
Assert[2*x^3 == 1 + 3*α*x^2];
x = DixonSM[λ/3, α];
y = DixonCM[λ/3, α];
Assert[x^2*y + x - y^2 == 0];
Assert[y^2*x + y - x^2 == b*x*y];
Assert[x^3 + 1 == b*x];
Assert[y^3 + 1 == b*y^2];
x = DixonSM[λ/4, α];
y = DixonCM[λ/4, α];
Assert[x + x^3 == y^3*(1 - x)];
Assert[DixonF[λ, α] == κ];
If[α != -1,
Assert[DixonF[λ1, α] == κ1];
Assert[DixonF[λ2, α] == κ2];
];
s = DixonSM[u, α];
c = DixonCM[u, α];
Assert[s^3 + c^3 == 1 + 3*α*s*c];
Assert[DixonSM[λ - u, α] == c];
Assert[DixonCM[λ - u, α] == s];
Assert[DixonSM[u + λ, α] == DixonCM[-u, α] == 1/c];
Assert[DixonCM[u + λ, α] == DixonSM[-u, α] == -s/c];
Assert[DixonSM[u - λ, α] == DixonCM[-u - λ, α] == -c/s];
Assert[DixonCM[u - λ, α] == DixonSM[-u - λ, α] == 1/s];
Assert[DixonSM[2*u, α] == s*(1 + c^3)/(c*(1 + s^3))];
Assert[DixonCM[2*u, α] == (c^3 - s^3)/(c*(1 + s^3))];
Assert[DixonSM[3*u, α] == s*c*(1 + s^3 + c^3 + s^6 + c^6 - s^3*c^3)/(c^3 - s^6 + 3*s^3*c^3 + c^6*s^3)];
Assert[DixonCM[3*u, α] == -(s^3 - c^6 + 3*s^3*c^3 + s^6*c^3)/(c^3 - s^6 + 3*s^3*c^3 + c^6*s^3)];
With[{v = λ1 - λ2}, If[NumericQ[v],
Assert[DixonSM[u + v, α] == w2*s];
Assert[DixonCM[u + v, α] == w1*c];
]];
Assert[Derivative[1, 0][DixonSM][u, α] == c^2 - α*s];
Assert[Derivative[1, 0][DixonCM][u, α] == -s^2 + α*c];
f[u_, α_] := DixonSM[u, α];
x = f[u, α];
y = Derivative[1, 0][f][u, α];
Assert[y^3 - 3*α*x*y^2 - x^6 + (2 + 4*α^3)*x^3 - 1 == 0];
If[α < 0 && x > 0, Assert[y*Hypergeometric2F1[1/3, 2/3, 1/2, (1 - x^3)^2/(4*(α*x)^3)] == -α*x]];
f[u_, α_] := DixonCM[u, α];
x = f[u, α];
y = Derivative[1, 0][f][u, α];
Assert[y^3 + 3*α*x*y^2 + x^6 - (2 + 4*α^3)*x^3 + 1 == 0];
If[α < 0 && x > 0, Assert[y*Hypergeometric2F1[1/3, 2/3, 1/2, (1 - x^3)^2/(4*(α*x)^3)] == α*x]];
f[u_, α_] := DixonSM[u, α] + DixonCM[u, α];
x = f[u, α];
y = Derivative[1, 0][f][u, α];
Assert[3*y^2 + (x + α)*(x^3 - 3*α*x^2 - 4) == 0];
f[u_, α_] := DixonSM[u, α]*DixonCM[u, α];
x = f[u, α];
y = Derivative[1, 0][f][u, α];
Assert[y^2 + 4*x^3 == (1 + 3*α*x)^2];
f[u_, α_] := If[PossibleZeroQ[1 + α], -1/3, 1/(DixonSM[u, α] + DixonCM[u, α] + α)];
x = f[u, α];
y = Derivative[1, 0][f][u, α];
Assert[3*y^2 == 4*(1 + α^3)*x^3 - (1 - 3*α*x)^2];
If[α < 0,
f[u_, α_] := With[{s = DixonSM[u, α], c = DixonCM[u, α]}, ArcTan[(1 - c - 2*α*s)/(Sqrt[3]*(1 + c))]];
x = f[u, α];
y = Derivative[1, 0][f][u, α];
Assert[y*Hypergeometric2F1[1/3, 2/3, 1/2, (1 + α^(-3))*Sin[x]^2] == -Sqrt[3]*α/2];
];
Clear[f];
If[Context[CarlsonRF] === "System`",
f[u_, v_] :=
Block[{x, y, z, b = Abs[α]^(3/2), ϕ},
ϕ = 2*I*Pi*{1, -1, 0};
ϕ = Which[
α > 0, (ArcCsch[b] + ϕ)/3,
α < 0, (ArcSech[b] + ϕ)/3
];
b = v*α;
{x, y, z} = Which[
b > 0, 1 + Csch[ϕ]^2/4,
b < 0, 1 - Sech[ϕ]^2/4,
True, Table[1, 3]
];
Assert[y == Conjugate[x] && z >= 3/4];
Assert[u == Re[CarlsonRF[x, y, z]]];
Assert[u == Sqrt[2]*CarlsonRF[Re[x] + Abs[x], z + Abs[x] + Abs[z - x], z + Abs[x] - Abs[z - x]]];
]
];
If[α != -1, With[{u = HypergeometricPFQ[{1/2, 1/2, 1}, {5/6, 7/6}, -α^3], s = sqrt[3*α]},
Assert[DixonSM[s*u, α] == s];
Assert[DixonCM[s*u, α] == 1];
Assert[DixonF[s*u, α] == (3*α)/2 + s*(3/5)*α^2*HypergeometricPFQ[{1/2, 1, 3/2}, {7/6, 11/6}, -α^3]];
Assert[DixonF[s*u + λ, α] == κ - (3*α)/2 + s*(3/5)*α^2*HypergeometricPFQ[{1/2, 1, 3/2}, {7/6, 11/6}, -α^3]];
Assert[PossibleZeroQ[α] || DixonF[s*u - λ, α] == -κ + (3*α)/2 + s^(-1)*HypergeometricPFQ[{-1/2, 1/2, 1}, {1/6, 5/6}, -α^3]];
If[α > -1, With[{v = u*Sqrt[1 + α^3]},
With[{b = SetPrecision[(3*α)^3/4, Infinity]},
Assert[u == (1/2)*NIntegrate[(x^3 + b*(x - 1)^2)^(-1/2), {x, 1, Infinity}, WorkingPrecision -> pr]]; (* hard case at x = 3 *)
Assert[u == NIntegrate[(1 + b*x^2*(1 - x^2)^2)^(-1/2), {x, 0, 1}, WorkingPrecision -> pr]]; (* hard case at x = 3^(-1/2) *)
];
If[!PossibleZeroQ[α], With[{b = SetPrecision[1 + α^(-3), Infinity], ϕ = SetPrecision[ArcTan[α^(3/2)], Infinity]},
Assert[v == NIntegrate[(4*x^3 - (1 - 3*x)^2/b)^(-1/2), {x, 1, Infinity}, WorkingPrecision -> pr]];
Assert[v == NIntegrate[(1 - x^2*(3 - x^2)^2/(4*b))^(-1/2), {x, 0, 1}, WorkingPrecision -> pr]];
Assert[u == (2/3)*α^(-3/2)*N[NIntegrate[Hypergeometric2F1[1/3, 2/3, 1/2, b*Sin[x]^2], {x, 0, ϕ},
WorkingPrecision -> pr, PrecisionGoal -> pr/2
], pr/2]];
]];
f[v, -1];
]];
]];
If[α > -1, With[{u = HypergeometricPFQ[{1/2, 1/2, 1}, {5/6, 7/6}, α^3/(1 + α^3)]/Sqrt[1 + α^3], s = sqrt[α]},
Assert[DixonSM[s*u, α] == s*E^((1/3)*ArcSinh[s^3])];
Assert[DixonCM[s*u, α] == E^((2/3)*ArcSinh[s^3])];
Assert[Derivative[1, 0][DixonSM][s*u, α] == Sqrt[1 + α^3]*E^((1/3)*ArcSinh[s^3])];
Assert[Derivative[1, 0][DixonCM][s*u, α] == 0];
If[!PossibleZeroQ[α], With[{b = SetPrecision[α^(-3), Infinity], ϕ = SetPrecision[ArcSinh[α^(3/2)], Infinity]},
Assert[u == NIntegrate[(4*x^3 + (1 - 3*x)^2/b)^(-1/2), {x, 1, Infinity}, WorkingPrecision -> pr]];
Assert[u == NIntegrate[(1 + x^2*(3 - x^2)^2/(4*b))^(-1/2), {x, 0, 1}, WorkingPrecision -> pr]];
Assert[u == (2/3)*α^(-3/2)*N[NIntegrate[Hypergeometric2F1[1/3, 2/3, 1/2, b*Sinh[x]^2], {x, 0, ϕ},
WorkingPrecision -> pr, PrecisionGoal -> pr/2
], pr/2]];
]];
f[u, 1];
]];
If[α <= 0, With[{b = SetPrecision[α^3, Infinity]},
Assert[λ == NIntegrate[(1 - b*x^3)^(-1/3)*(1 + x^3)^(-2/3), {x, 0, Infinity}, WorkingPrecision -> pr]];
Assert[λ == NIntegrate[x*(x^3 - b)^(-1/3)*(1 + x^3)^(-2/3), {x, 0, Infinity}, WorkingPrecision -> pr]];
]];
f[u_, α_] := Through[{DixonSM, DixonCM}[u, α]]; (* ComapApply *)
If[!PossibleZeroQ[γ],
Assert[f[w1*u, w2*α] == {w1*s, c}];
Assert[f[u*(1 + α)/(w2 - w1), (2 - α)/(1 + α)] == -{s + c - 1, s + w2*c - w1}/(s + w1*c - w2)];
Assert[f[3*u*γ/(w2 - w1), -α/γ] == -{3*γ*s*c, s^3 + w2*c^3 - w1}/(s^3 + w1*c^3 - w2)];
Assert[f[u*γ/(w2 - w1), -α/γ]^3 == -{α*s - c + 1, α*s - w2*c + w1}/(α*s - w1*c + w2)];
Assert[f[b*u, (α - 2)/b] == {b*s*c, -(s - c^2 + c*s^2)}/(c - s^2 + s*c^2)];
Assert[f[b*u/3, (α - 2)/b]^3 == {α*s - c + 1, α*c - s + 1}/(s + c + α)];
Assert[f[(γ - α)*u, (2*γ + α)/(γ - α)] == {s^3 + (γ - α)*s*c, c^3 + (γ - α)*s*c}/(1 - (γ - α)*s*c)];
];
Block[{u, v},
{u, v} = RandomReal[{0, λ0}, 2, WorkingPrecision -> pr];
With[{s1 = DixonSM[u], c1 = DixonCM[u], s2 = DixonSM[v], c2 = DixonCM[v]},
Assert[DixonSM[u + w1*v] == (s1*c2*(c2 + s1*c1*s2^2) + w1*s2*c1*(c1 + s2*c2*s1^2))/(1 - s1^3*s2^3)];
Assert[DixonCM[u + w1*v] == (c1*(c2 + s1*c1*s2^2) + w1*s1*s2*(s2*c1^2 - s1*c2^2))/(1 - s1^3*s2^3)];
];
];
] /.
HoldPattern[_[vars:{___Symbol}, expr_]] :>
Block[{$AssertFunction = Automatic},
With[{temp = Unique[Unevaluated[vars], Temporary]},
WithCleanup[
Unevaluated[expr] /. Dispatch[Thread@Hold[vars, temp] /. Hold[lhs_, rhs_] :> (HoldPattern[lhs] -> rhs)],
Remove @@ Unevaluated[temp]
] /; ListQ[temp]
]
];Theta functions:
Module[{u, q, a, b},
u = RandomReal[WorkingPrecision -> 50];
q = RandomReal[WorkingPrecision -> 50];
a = ((EllipticTheta[4, q]/EllipticTheta[4, Pi/3, q])^2 + 2*EllipticTheta[4, Pi/3, q]/EllipticTheta[4, q])/3;
b = EllipticTheta[2, q]*EllipticTheta[3, q]*EllipticTheta[4, q]/(2*Sqrt[3]*EllipticTheta[1, Pi/3, q]);
Assert[EllipticTheta[1, u, q^3]^3 == b*(a*EllipticTheta[1, u, q] - EllipticTheta[1, u + Pi/3, q] + EllipticTheta[1, u + 2*Pi/3, q])];
Assert[EllipticTheta[2, u, q^3]^3 == b*(a*EllipticTheta[2, u, q] - EllipticTheta[2, u + Pi/3, q] + EllipticTheta[2, u + 2*Pi/3, q])];
Assert[EllipticTheta[3, u, q^3]^3 == b*(a*EllipticTheta[3, u, q] + EllipticTheta[3, u + Pi/3, q] + EllipticTheta[3, u + 2*Pi/3, q])];
Assert[EllipticTheta[4, u, q^3]^3 == b*(a*EllipticTheta[4, u, q] + EllipticTheta[4, u + Pi/3, q] + EllipticTheta[4, u + 2*Pi/3, q])];
Assert[EllipticTheta[1, u, q^3] == EllipticTheta[1, u/3, q]*EllipticTheta[1, (u + Pi)/3, q]*EllipticTheta[1, (u + 2*Pi)/3, q]/(3*b)];
Assert[EllipticTheta[2, u, q^3] == -EllipticTheta[2, u/3, q]*EllipticTheta[2, (u + Pi)/3, q]*EllipticTheta[2, (u + 2*Pi)/3, q]/(3*b)];
Assert[EllipticTheta[3, u, q^3] == EllipticTheta[3, u/3, q]*EllipticTheta[3, (u + Pi)/3, q]*EllipticTheta[3, (u + 2*Pi)/3, q]/(3*b)];
Assert[EllipticTheta[4, u, q^3] == EllipticTheta[4, u/3, q]*EllipticTheta[4, (u + Pi)/3, q]*EllipticTheta[4, (u + 2*Pi)/3, q]/(3*b)];
Assert[EllipticTheta[1, u, q^9] == -(EllipticTheta[1, u/3, q] - EllipticTheta[1, (u + Pi)/3, q] + EllipticTheta[1, (u + 2*Pi)/3, q])/3];
Assert[EllipticTheta[2, u, q^9] == (EllipticTheta[2, u/3, q] - EllipticTheta[2, (u + Pi)/3, q] + EllipticTheta[2, (u + 2*Pi)/3, q])/3];
Assert[EllipticTheta[3, u, q^9] == (EllipticTheta[3, u/3, q] + EllipticTheta[3, (u + Pi)/3, q] + EllipticTheta[3, (u + 2*Pi)/3, q])/3];
Assert[EllipticTheta[4, u, q^9] == (EllipticTheta[4, u/3, q] + EllipticTheta[4, (u + Pi)/3, q] + EllipticTheta[4, (u + 2*Pi)/3, q])/3];
With[{s = EllipticTheta[1, u, q], m = EllipticTheta[1, u + Pi/3, q], c = EllipticTheta[1, u + 2*Pi/3, q]},
Assert[EllipticTheta[1, u+Pi/9, q]^3 == ( EllipticTheta[1, Pi/9, q]^3*c*m^2 + EllipticTheta[1, 2*Pi/9, q]^3*s*c^2 + EllipticTheta[1, 4*Pi/9, q]^3*m*s^2)/EllipticTheta[1, Pi/3, q]^3]
Assert[EllipticTheta[1, u-Pi/9, q]^3 == (-EllipticTheta[1, Pi/9, q]^3*m*c^2 + EllipticTheta[1, 2*Pi/9, q]^3*s*m^2 - EllipticTheta[1, 4*Pi/9, q]^3*c*s^2)/EllipticTheta[1, Pi/3, q]^3];
];
] & // Block[{$AssertFunction = Automatic}, #[]] &;Last edited by lanxiyu (2025-04-08 12:51:10)
Offline