"PSLQ for dilogarithm identities"의 두 판 사이의 차이
둘러보기로 가기
검색하러 가기
introduction==
Implement the PSLQ algorithm first.
encyclopedia==
books==
imported>Pythagoras0 |
imported>Pythagoras0 |
||
| 134번째 줄: | 134번째 줄: | ||
* [http://www.research.att.com/~njas/sequences/index.html The On-Line Encyclopedia of Integer Sequences] | * [http://www.research.att.com/~njas/sequences/index.html The On-Line Encyclopedia of Integer Sequences] | ||
* http://functions.wolfram.com/ | * http://functions.wolfram.com/ | ||
| − | [[분류:개인노트]] | + | [[분류:개인노트]][[분류:개인노트]] |
2012년 10월 28일 (일) 15:58 판
introduction==
Implement the PSLQ algorithm first.
-
PSLQ[inx_List, prec_] :=
Block[
{
x,
n = Length[inx],
\[Gamma] = 2/Sqrt[3],
A, B, H, D, Dinv, t, i, j, k, l, iter,
\[Alpha], \[Beta], \[Lambda], \[Delta], r, R
},
(*Initialize*)
x = N[inx /Sqrt[inx . inx], prec];
s = Sqrt[MapIndexed[Plus @@ Drop[x^2, First[#2] - 1] &, x]];
A = B = IdentityMatrix[n];
H = Table[Which[
i > j, (-xi*xj)/(sj*sj + 1),
i == j, si + 1/si,
i < j, 0
], {i, 1, n}, {j, 1, n - 1}];
(* Reduce H *)
t = HermiteReduce[H];
D = First[t];
Dinv = Inverse[D];
(*Update*)
H = Last[t]; x = x.Dinv; A = D.A; B = B.Dinv;
For[iter = 0, iter < $IterationLimit, ++iter,
(* Step One *)
r = MaxIndex[MapIndexed[\[Gamma]^First[#2] Abs[#1] &, Tr[H, List]]];
If[r < n - 1, \[Alpha] = Hr, r; \[Beta] =
Hr + 1, r; \[Lambda] = Hr + 1, r + 1; \[Delta] =
Sqrt[\[Beta]^2 + \[Lambda]^2]];
R = IdentityMatrix[n]; t = Rr; Rr = Rr + 1;
Rr + 1 = t;
x = x.R; H = R.H; A = R.A; B = B.R;
(* Step Two *)
If[r < n - 1,
H = H.Table[
Which[
i == r && j == r, \[Beta]/\[Delta],
i == r && j == r + 1, -\[Lambda]/\[Delta],
i == r + 1 && j == r, \[Lambda]/\[Delta],
i == r + 1 && j == r + 1, \[Beta]/\[Delta],
i == j && j != r || i == j && j != r + 1, 1,
True, 0],
{i, 1, n - 1}, {j, 1, n - 1}]
];
(* Step Three *)
t = HermiteReduce[H];
D = First[t];
Dinv = Inverse[D];
(*Update*)
H = Last[t]; x = x.Dinv; A = D.A; B = B.Dinv;
(* Step Four *)
If[Min[Abs[Union[x, Tr[H, List]]]] <= 10^(-prec + 5), Break[]]
];(*Main Iteraton*)
Return[Transpose[B][[MaxIndex[-Abs[x]]]]]
]
Then find a relation.
- Clear[a]
f[x_] := x^3 + x - 1
Solve[f[x] == 0, x]
N[%]
a := -(2/(3 (9 + Sqrt[93])))^(1/3) + (1/2 (9 + Sqrt[93]))^(1/3)/3^(2/3)
N[a, 20]
L[x_] := PolyLog[2, x] + 1/2 Log[x] Log[1 - x]
- L[1] := Pi^2/6
- (* choose a expected height of the dilogarithm identity *)
- deg:=6
S := Table[L[a^i],{i,0,deg}]
N[S, 100]
PSLQ[N[S, 100], 1000]
- N[N[S, 50].%, 50]
I found
\(-2L(1)+2L(\alpha)+2L(\alpha^{2})-L(\alpha^{4})=0\) or
\(2L(\alpha)+2L(\alpha^{2})-L(\alpha^{4})=\frac{\pi^2}{3}\)
PSLQ[inx_List, prec_] :=
Block[
{
x,
n = Length[inx],
\[Gamma] = 2/Sqrt[3],
A, B, H, D, Dinv, t, i, j, k, l, iter,
\[Alpha], \[Beta], \[Lambda], \[Delta], r, R
},
(*Initialize*)
x = N[inx /Sqrt[inx . inx], prec];
s = Sqrt[MapIndexed[Plus @@ Drop[x^2, First[#2] - 1] &, x]];
A = B = IdentityMatrix[n];
H = Table[Which[
i > j, (-xi*xj)/(sj*sj + 1),
i == j, si + 1/si,
i < j, 0
], {i, 1, n}, {j, 1, n - 1}];
(* Reduce H *)
t = HermiteReduce[H];
D = First[t];
Dinv = Inverse[D];
(*Update*)
H = Last[t]; x = x.Dinv; A = D.A; B = B.Dinv;
For[iter = 0, iter < $IterationLimit, ++iter,
(* Step One *)
r = MaxIndex[MapIndexed[\[Gamma]^First[#2] Abs[#1] &, Tr[H, List]]];
If[r < n - 1, \[Alpha] = Hr, r; \[Beta] =
Hr + 1, r; \[Lambda] = Hr + 1, r + 1; \[Delta] =
Sqrt[\[Beta]^2 + \[Lambda]^2]];
R = IdentityMatrix[n]; t = Rr; Rr = Rr + 1;
Rr + 1 = t;
x = x.R; H = R.H; A = R.A; B = B.R;
(* Step Two *)
If[r < n - 1,
H = H.Table[
Which[
i == r && j == r, \[Beta]/\[Delta],
i == r && j == r + 1, -\[Lambda]/\[Delta],
i == r + 1 && j == r, \[Lambda]/\[Delta],
i == r + 1 && j == r + 1, \[Beta]/\[Delta],
i == j && j != r || i == j && j != r + 1, 1,
True, 0],
{i, 1, n - 1}, {j, 1, n - 1}]
];
(* Step Three *)
t = HermiteReduce[H];
D = First[t];
Dinv = Inverse[D];
(*Update*)
H = Last[t]; x = x.Dinv; A = D.A; B = B.Dinv;
(* Step Four *)
If[Min[Abs[Union[x, Tr[H, List]]]] <= 10^(-prec + 5), Break[]]
];(*Main Iteraton*)
Return[Transpose[B][[MaxIndex[-Abs[x]]]]]
]
f[x_] := x^3 + x - 1
Solve[f[x] == 0, x]
N[%]
a := -(2/(3 (9 + Sqrt[93])))^(1/3) + (1/2 (9 + Sqrt[93]))^(1/3)/3^(2/3)
N[a, 20]
L[x_] := PolyLog[2, x] + 1/2 Log[x] Log[1 - x]
S := Table[L[a^i],{i,0,deg}]
N[S, 100]
PSLQ[N[S, 100], 1000]