Reminder: The OEIS is hiring a new managing editor, and the application deadline is January 26.
%I #31 Sep 20 2023 07:46:51
%S 6,8,7,8,4,1,8,1,0,3,2,8,3,8,9,2,6,3,2,7,1,3,4,4,0,4,4,0,9,8,8,3,3,4,
%T 8,6,1,1,5,8,3,9,7,9,4,8,7,6,6,8,9,5,4,1,1,7,4,7,5,8,6,6,9,4,4,1,0,7,
%U 8,5,2,8,1,7,2,1,2,4,7,5,3,8,9,1,0,8,7,9,1,2,6,5,7,8,9,7,8,5,3,6
%N Put a positive charge at 0 and a negative charge at 1, then keep adding alternating charges at points of zero potential; this is the decimal expansion of the limit.
%C The potential is inversely proportional to the distance: a positive charge at p has potential 1/abs(x-p). Starting from a positive charge at 0 and a negative charge at 1, the potential is zero at the midpoint 1/2, so another positive charge is added there. Then the potential is zero at 0.7886... and a negative charge is added there, and so on. The first few zero potential points are: 0.5, 0.7886..., 0.6325..., 0.7179..., 0.6714..., ... The limit is the constant of this sequence.
%C After initial points p_0 = 0 and p_1 = 1, each new zero point fits between previous two points and is the solution to Sum_{n>=0} 1/(x-p_n) = 0.
%H Rok Cestnik, <a href="/A365125/a365125.png">Visualization</a>
%e 0.68784181032838926327134404409883348611583979...
%t p={0,1};
%t For[r=1,r<39,++r,p=Append[p,x/.NSolve[{Sum[1/(x-p[[k]]),{k,1,Length[p]}]==0,(x-p[[-1]])*(x-p[[-2]])<0},x,WorkingPrecision -> 30][[1,1]]];];
%t A365125 = RealDigits[Floor[10^10*p[[-1]]]][[1]]
%o (Python)
%o from decimal import Decimal, getcontext, ROUND_UP, ROUND_DOWN
%o getcontext().prec = 100
%o def nm(f, df, x):
%o for i in range(10):
%o x -= f(x)/df(x)
%o return x
%o def flip_rounding():
%o if getcontext().rounding == ROUND_UP: getcontext().rounding = ROUND_DOWN
%o else: getcontext().rounding = ROUND_UP
%o def get_zero(vs, rounding):
%o getcontext().rounding = rounding
%o def p(x,v):
%o flip_rounding(); t = x-v
%o flip_rounding(); return 1/t
%o def dp(x,v):
%o flip_rounding(); t = x-v; t = t**2
%o flip_rounding(); return -1/t
%o def f(x): return sum(p(x,vs[n]) for n in range(len(vs)))
%o def df(x): return sum(dp(x,vs[n]) for n in range(len(vs)))
%o sign = -1 if rounding == ROUND_DOWN else 1
%o return nm(f, df, (vs[-1]+vs[-2])/2+sign*abs(vs[-1]-vs[-2])/3)
%o v_lo = [Decimal(0), Decimal(1)]
%o v_up = [Decimal(0), Decimal(1)]
%o for r in range(150):
%o v_lo.append(get_zero(v_lo, ROUND_DOWN))
%o v_up.append(get_zero(v_up, ROUND_UP))
%o lower_bounds = [v_lo[i] for i in range(0, len(v_lo), 2)]
%o upper_bounds = [v_up[i] for i in range(1, len(v_up), 2)]
%o right = True
%o A365125 = [int(l) for l, u in zip(str(lower_bounds[-1])[2:], str(upper_bounds[-1])[2:]) if right and (right := (l == u))]
%K nonn,cons
%O 0,1
%A _Rok Cestnik_, Aug 22 2023