Solved – Mathematically, how are the critical t values calculated

Critical values are quantiles of the distribution. There are numerous ways to approximate quantiles (inverse cdf) of the t-distribution.

Since you have R already, one way is readily found by looking at how R does it; this is referenced in the help (?qt)), under "Source":

For central qt, a C translation of

Hill, G. W. (1970) Algorithm 396: Student's t-quantiles. Communications of the ACM, 13(10), 619–620.

altered to take account of

Hill, G. W. (1981) Remark on Algorithm 396, ACM Transactions on Mathematical Software, 7, 250–1.

To my understanding Hill uses a modified version of a Cornish-Fisher expansion.

Other approximations are in Abramowitz and Stegun, and if I recall correctly, the tables by Pearson and Hartley as well as a number of more recent ones (both those references are even older than Hill).

If you already have an inverse beta cdf (via inverse incomplete beta function, which is available in some libraries), you can do that way — the relationship between the t cdf and the regularized incomplete beta function is given at the Wikipedia page for the t-distribution.