On $\int_0^1\arctan\,_6F_5\left(\frac17,\frac27,\frac37,\frac47,\frac57,\frac67;\,\frac26,\frac36,\frac46,\frac56,\frac76;\frac{n}{6^6}\,x\right)\,dx$
I think you should start by the following general form
$${}_{k}F_{k-1}\left(\frac{1}{k+1} ,\cdots ,\frac{k}{k+1};\frac{2}{k} \cdots ,\frac{k-1}{k},\frac{k+1}{k};\left( \frac{m(1-m^k)}{f_k}\right)^k \right) = \frac{1}{1-m^k}$$
Where
$$f_k \equiv \frac{k}{(1+k)^{1+1/k}}$$
Put $k=4$
$${}_{4}F_{3}\left(\frac{1}{5} ,\cdots ,\frac{4}{5};\frac{2}{4} \cdots ,\frac{3}{4},\frac{5}{4};\left( \frac{m(1-m^4)}{f_4}\right)^4 \right) = \frac{1}{1-m^4}$$
The argument simplifies to
$$ \left(\frac{m(1-m^4)}{\frac{4}{5^{1+1/4}}}\right)^4 = \frac{5^5 \,m^4(1-m^4)^4}{4^4}$$
Hence we have
$${}_{4}F_{3}\left(\frac{1}{5} ,\cdots ,\frac{4}{5};\frac{2}{4} \cdots ,\frac{3}{4},\frac{5}{4};\frac{5^5 (1-m)m^4}{4^4} \right) = \frac{1}{m}$$
Now suppose we want to find
$$\int \arctan\,_4F_3\left(\frac15,\frac25,\frac35,\frac45;\,\frac24,\frac34,\frac54;\,\frac{n}{4^4}\,x\right)\,\mathrm dx$$
Use $nx = 5^5m^4(1-m)$
$$\frac{5^5}{n}\int (4m^3(1-m)-m^ 4)\arctan\,_4F_3\left(\frac15,\frac25,\frac35,\frac45;\,\frac24,\frac34,\frac54;\,\frac{5^5 (1-m)m^4}{4^4}\right)\,\mathrm dm$$
This simplifies to
$$\frac{5^5}{n}\int (4m^3(1-m)-m^ 4)\arctan\,\frac{1}{m}\,\mathrm dm$$
The anti-derivative of this is elementary.
Note that when we use substittution we will need the roots of
$$m^5-m^4+\frac{n}{5^5}= 0$$
As conjectured by the OP.
Using the value $n=4$ we can verify Reshetnikov result
\begin{align}\int_0^1\arctan{_4F_3}\left(\frac15,\frac25,\frac35,\frac45;\frac24,\frac34,\frac54;\frac{1}{64}\,x\right)\,dx &= \frac{5^5}{4}\int^{\alpha}_1 (4x^3(1-x)-x^ 4)\arctan \left(\frac{1}{x} \right)\,\mathrm dx\\&=0.7857194\dots \end{align}
Where $\alpha$ is the largest positive root of
$$m^5-m^4+\frac{4}{5^5}= 0$$
The result does generalize. We will derive a closed formula for the integral $$I_k(n) = \int_0^1\arctan\left[_{k}F_{k-1}\left(\frac{1}{k+1},\cdots, \frac{k}{k+1},\frac{2}{k},\cdots,\frac{k+1}{k}; \frac{nx}{k^k}\right)\right]\,{\rm d}x$$ valid for $|n| \leq k^k$.
For this particular combination of arguments the hypergeometrical function simplifies greatly, see property $(25)$ here$^{(1)}$. Using this we can perform a substitution $\frac{nx}{(k+1)^{k+1}} = y^k(1-y^k)^k$ and the integrand simplifies to $\frac{(k+1)^{k+1}}{n}\frac{1}{1-y^k}\frac{d}{dy}\left(y^k(1-y^k)^k\right)$. Now performing integration by parts and the substitutions $y^k\to y \to 1-y$ we arrive at
$$I_k(n) = \frac{\pi}{2} - \arctan(z_0) + \frac{(k+1)^{k+1}}{n}\int_{z_0}^1\frac{w^{k+1}-w^k}{1+w^2}{\rm d}w$$
where $z_0$ is the largest positive root of $z_0^{k+1} - z_0^k + \frac{n}{(k+1)^{k+1}} = 0$. The integral above can be evaluated in closed form using
$$\matrix{\int \frac{w^{2k}}{1+w^2}{\rm d}w &=& (-1)^k\left[\arctan(w) + F_k(w)\right] + C\\\int \frac{w^{2k+1}}{1+w^2}{\rm d}w &=& (-1)^k\left[\frac{1}{2}\log(1+w^2) + G_k(w)\right] + C}$$ where $F_k(x) = \sum_{i=1}^{k}\frac{(-1)^i x^{2i-1}}{2i-1}$ and $G_k(x) = \sum_{i=1}^{k}\frac{(-1)^i x^{2i}}{2i}$. These relations are easily proven by differentation and summing a geometrical series. This leads to
$$I_{2k}(n) = \frac{\pi}{2} - \arctan(z_0) + \frac{(-1)^k(2k+1)^{2k+1}}{n}\left[\frac{1}{2}\log\left(\frac{2}{1+z_0^2}\right) - \frac{\pi}{4} + \arctan(z_0)\right] + \frac{(-1)^k(2k+1)^{2k+1}}{n}\left[G_k(1) - G_k(z_0) - F_k(1) + F_k(z_0)\right]$$ $$I_{2k+1}(n) = \frac{\pi}{2} - \arctan(z_0) + \frac{(-1)^{k}(2k+2)^{2k+2}}{n}\left[-\frac{1}{2}\log\left(\frac{2}{1+z_0^2}\right) - \frac{\pi}{4} + \arctan(z_0)\right] + \frac{(-1)^{k}(2k+2)^{2k+2}}{n}\left[- G_k(1) + G_k(z_0) -F_{k+1}(1) + F_{k+1}(z_0)\right]$$
For the particular case $k = 2$ ($p=3$ in the question) you asked for this reduces to
$$I_2(1) = \frac{27}{2}\log \left(\frac{z_0^2+1}{2}\right) - \frac{27}{2}(z_0-1)^2 - 28\arctan(z_0) + \frac{29 \pi}{4} \simeq 0.795296 \\ \text{where } z_0 = \frac{1}{3} \left(1+2 \cos \left(\frac{\pi }{9}\right)\right)$$
Finlly since the expression has many terms we should check that it gives the correct result. Here is some code used to check this in Mathematica:
(* Working precision *)
ndigits = 100;
Block[{$MinPrecision = ndigits},
(* Pick k and n *)
k = 4; n = 1;
(* Solve for z0 *)
z0 = Max[N[x /. Solve[{x^(k + 1) - x^k + n/(k + 1)^(k + 1) == 0, x > 0}, x], ndigits]];
(* Define the integral *)
ak = Table[i/(k + 1), {i, 1, k}];
bk = Table[If[i == k, (i + 1)/k, i/k], {i, 2, k}];
exact = NIntegrate[ArcTan[HypergeometricPFQ[ak, bk, (x n)/k^k]], {x, 0, 1}, WorkingPrecision -> ndigits];
(* Our solution *)
result = If[Mod[k, 2] == 0,
(* For even k *)
\[Pi]/2 - ArcTan[z0] + ((-1)^(k/2) (k + 1)^(k + 1))/ n (1/2 Log[2/(1 + z0^2)] - \[Pi]/4 + ArcTan[z0] + Sum[(-1)^i ((1 - z0^(2 i))/(2 i) - (1 - z0^(2 i - 1))/(2 i - 1)), {i, 1, k/2}])
,
(* For odd k *)
\[Pi]/2 - ArcTan[z0] + ((-1)^((k - 1)/2) (k + 1)^(k + 1))/n (-(1/2) Log[2/(1 + z0^2)] - \[Pi]/4 + ArcTan[z0] - Sum[(-1)^i ((1 - z0^(2 i))/(2 i)If[i == (k - 1)/2 + 1, 0, 1] + (1 - z0^(2 i - 1))/(2 i - 1)), {i, 1, (k - 1)/2 + 1}])
];
(* Difference between integral and our formula *)
Print["Difference: ", N[result - exact]];
];
This demonstrates that we get the same result (to numerical precision) for a wide range of $k$ and $n$.
$^{(1)}$ It's easy to check that this holds to arbitrary precisison (checked to a few hundred digits) at least for the smallest values of $k$. Marty Cohen found this reference discussing a very similar type of property which might be useful if anyone wants to try to derive it (which likely is a very messy affair).