default(realprecision, 80);
bestappr(lfun(lfuncreate(x^4 - x^3 - 3*x^2 + x + 1), -1), 30) \\ 2/15, for D = 725The generator enumerates the totally real quartic fields with $D\leq 30000$ by PARI's nflist [7] over the five transitive groups of degree $4$, and refuses the run unless it finds $204$ fields.
Each value is recognised from PARI lfun computations at several working precisions, with the denominator bounded by Serre's theorem [3]: in degree four the possible denominators divide $30$, $1020$ and $8190$ at $s=-1,-3,-5$. The generator refuses a value unless the same rational is obtained at all working precisions and the functional equation predicts PARI's computations of $\zeta_K(2)$, $\zeta_K(4)$ and $\zeta_K(6)$ to relative error below $10^{-60}$. For abelian quartic fields the value is also checked in exact rational arithmetic from the Dirichlet-character factorisation (2).