I would like an efficient algorithm for square root of a positive integer. Is there a reference that compares various square root algorithms? Thanks.
-
3Not sure this really qualifies as a research level question, you might be better off posting on MathStackexchage. But a very fast and easy algorithm to compute square root of $A$ to $N$ decimal places for some reasonable $N$ is Newton's algorithm to find a root of $X^2 - A$. So pick an initial $x_0$, say $x_0=A/2$, and then iteratively compute $x_{i+1} = \dfrac{x_i}{2}-\dfrac{A}{2x_i}$. – Jul 08 '13 at 19:02
-
The documentation and source code to GMP could be used as a reference. Also, do note that the size of the integer is important; different methods are better for differently sized integers. Any other special properties the integer might have (such as being a perfect square, or near a power of 2) are also very relevant. – Jul 08 '13 at 19:42
-
3Do you want the integer square root? That is, the integer $i$ such that $i^2 \le n \lt (i+1)^2$? Or do you want an approximate rational number $r$ with $r^2\approx n$? And do you want an algorithm that is suitable for pencil-and-paper, for an electronic digital computer, or for something else? – MJD Jul 08 '13 at 19:42
-
See https://en.wikipedia.org/wiki/Integer_square_root. – lhf Jan 04 '16 at 12:27
-
How large are these integers? – gnasher729 Jul 11 '23 at 10:41
4 Answers
Method 1
Newton's Method converges quadratically. $a \ge 0, x_0 \ne \sqrt{a}$
$$x_{k+1} = \dfrac{1}{2}\left(x_{k} +\dfrac{a}{x_{k}}\right)$$
Method 2
$n \ge 0, x_n \gt \sqrt{a}$ for all $n \gt 0$.
$$x^2_{k+1}-a = \left[\dfrac{x^2_{k} - a}{2x_n}\right]^2$$
Method 3
Third order method, $n \ge 0$
$$x_{n+1} = \dfrac{x_n(x_n^2 + 3a)}{2x^2_n+a}$$
Method 4
See Math World Bhaskara-Brouncker algorithm
Additional Methods (also see references)
- 56,093
-
-
1@Daryl: Thank you for catching that, I accidentally deleted it before posting when I was cleaning up! Thanks! Regards – Amzoti Jul 08 '13 at 22:38
-
-
@amWhy: Thank you. Surprising how many interesting algorithms there are for that! ;-) – Amzoti Jul 09 '13 at 01:19
-
None of this would work if you are limited to integer arithmetic, which OP's supposition that $n$ is an integer often presupposes as well... – gt6989b Aug 03 '23 at 12:38
"efficient" rather depends on what you constraints are. For instance, if you have enough memory to store some floats, but little cpu time, then you can store a look up table for all integers until some power of 10. So say you had to find $\sqrt(1632397825)$. You could write: $$\sqrt{1732397825} =\sqrt{ 17*10^8 + 32397825}\sim \sqrt{17}\cdot 10^4 + \epsilon$$
To calculate $\epsilon$, use the fact that: $$\sqrt{a^2+b}\sim a + \frac{b}{2a} - \frac{b^2}{8a^3} + ...$$ So in our example, $$\sqrt{1732397825} \sim 41622.06575$$ Quite close to the true value of: $$\sqrt{1732397825}=41622.08338$$
- 33,018
Overview of Starting Taking a Square Root
Generally we can consider a square number as $S^2$ that is comprised of multiple digits. For an irrational square root $S$ there is no termination in the sequence of digits that represent the root and therefore we can only consider that root to some degress of precision. But if $S$ is a real number that can be represented by integer digits, we can divide the classes of roots according to if $S^2 \geq 1$ which results in $S \geq 1$ or otherwise $0< S^2 < 1$ which results in $0<S<1$. To find the square root, we start with real numbers of two cases and we reduce with a simple division that root to another generalized root of the form:
$$\text{Case 1: } {(S')}^2=\frac{S^2}{((c_n)^2 \times 10^n)}=\left( { \sqrt{a_{0} . + a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+...} } \right) ^2$$
$$\text{Case 2: } {(S')}^2=\frac{S^2}{((c_n)^2 \times 10^n)}=\left( { \sqrt{ a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+...} } \right) ^2$$
From that generalized format, we can use a very large lookup table (preferred because of large memory available as well as minimal time required to load pre-computed sequences) to find the square root, just storing the digits to the right of the decimal place of course for $S'$. Or we could also use any of the alternative techniques listed in the other answers to arrive at the root $S'$. And after that , we make a single integer multiplication, and we are finished finding our root to the desired level of accuracy.
$$\text{Case 1: } S=(S') \times (c_n)*10^{n/2} $$
$$\text{Case 1: } S=\left( (c_n) \times 10^{\frac{n}{2}} \right) \times \left( { \sqrt{a_{0} . + a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+...} } \right) $$
$$\text{Case 2: } S=(S') \times (c_n)*10^{n/2} $$ $$\text{Case 2: } S= \left( (c_n) \times 10^{\frac{n}{2}} \right) \times \left( { \sqrt{0. + a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+...} } \right)$$
Lemma Trivial Case for 1: If S² = 1 then S = 1
By inspection, we know that if $S^2$ is rational, and not imaginary, and also that $S$ is rational (and not imaginary) then $S=1$. Obviously we have here that $\sqrt{1} = 1 $.
Lemma Case 1: If S² > 1 then S > 1
This can be proven by writing out the formula for $(1+\epsilon)^2=1+(2+\epsilon) \times \epsilon$. In the case that $S^2 \ge 1$ then $(1+\epsilon)^2 \ge 1$. That can only happen if $\epsilon \ge 0$ so that both $(2+\epsilon)>0$ and also as just stated $\epsilon \ge 0$.
Lemma Case 2: If 0 < S² < 1 then 0 < S < 1
In the case that $S^2 \le 1$ then $(1+\epsilon)^2 \le 1$. That can only happen if $-2 < \epsilon \le 0$ so that both $(2+\epsilon)>0$ and also as stated just now $\epsilon \le 0$. For the case $1^2=1$, by inspection the real root is $\sqrt{1^2}=1$. Otherwise $(2+\epsilon) \ge (2+(-2))$ and $(2+\epsilon)\ge 0$ and $\epsilon < 0$ so that the product $(2+\epsilon)\times \epsilon<0$ and $S^2 \le 1$ as well as $S \le 1$.
Analysis
The two non-trivial $\sqrt{1}=1$ root types are then:
$$\text{Case 1: } S^2 > 1 \implies S^2 = \left( \sqrt{ a_{0} . + \sum_{i=i_{min}}^{i_{max}} a_{-i} \times 10^{-i} } \right)^2$$
$$\text{Case 2: } S^2 < 1 \implies S^2 = \left( \sqrt{\sum_{i={-({i_{min}})}}^{-(i_{max})} a_{-i} \times 10^{-i} } \right)^2$$
where $a_i$ are the limited number of digits of the real number $S^2$. Obviously, a condition of infinite digits cannot be computed with a limited computer.
Square Roots of S² Greater or Less than 1
Case 1: S² > 1 and Therefore S > 1
Let $S \in \text{{Integers>0}}$. Then if we want to take the square root of $S^2$ then we can proceed to label the digits $a_i \in \text{{Integers 0..9}}$ and $b_j \in \text{{Integers 0..9}}$ as follows (where $n$ represents the maximum power of 10 for the number $S^2$ and $n/2=\text{Floor(}n/2\text{)})$:
$$S^2 = \sum_{i=0}^n a_i \times 10^i$$
$$S = \sum_{j=0}^{n/2} b_j \times 10^j$$
Then, also we have:
$$S^2= (\sum_{j=0}^{n/2} b_j \times 10^j)^2=(b_{n/2})^2 \times 10^n + \text{ Terms with lower powers of 10}$$
We can then seek the maximum digit $b_{n/2}$ such that $((c_n)^2=(b_{n/2})^2) \le (a_n)$ for $n \in \text{\{Set of Even Integers\}}$.
Or if $n \in \text{\{Set of Odd Integers\}}$, we can then seek the maximum digit $b_{n/2}$ such that $(b_{n/2})^2 \le (a_n)*10+a_{n-1}$ with here $(c_n)^2=(b_n)^2$.
Now Finding the Square Root of 1 Plus a Small Remainder (Case 1)
We go ahead and divide $S^2$ as follows:
$$\text{Case 1: } {(S')}^2=\frac{S^2}{((c_n)^2 \times 10^n)} $$ $$ {(S')}^2 =\left( \sqrt{ a_{0} . + a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+...} \right)^2$$
$$\text{Case 2: } {(S')}^2=\frac{S^2}{((c_n)^2 \times 10^n)} $$ $$ {(S')}^2 =(\sqrt{a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+...})^2$$
Now we can reduce the problem of taking a square root of an integer sequence to finding the square root of
$$\text{Case 1: } \sqrt{a_{0} . + a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+...}$$ or with $$\text{Case 2: } \sqrt{a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+...}$$
With modern computing allowing for the storage of a certain number of digits in a memory look-up table; it is a very fast and efficient means of finding a square root. The initial digit division and final integer-digit multiplication take up a hard-to-measure small amount time with modern computers.
Considering Case 1, there is a question of how many digits $-k_{min}$ to the right of $\sqrt{a_{0}}.$ in the sequences $$\sqrt{a_{0}+\sum_{k=-1}^{k=k_{min}} a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+...+ a_{k_{min}} \times 10^{k_{min}}}$$ to calculate the root in the lookup table. Modern computers have several GB of storage available usually to store a massive lookup table with ease.
Both $\text{Case 1 }$ and $\text{Case 2 }$ can be stored as lookup tables.
Even with smaller memory available, the other question is how can other digits beyond the digits in the lookup table be figured out with a limited lookup table? That question detracts from the main method, so first the main method of taking the root is presented. And what to do to increase resolution beyond the lookup-table method (or other answers to this question already) can be discussed further after.
Case 2: 0 < S² < 1 and Therefore 0 < S < 1
Finding the Square Root of 1 a Fractional Remainder
One prefered option is to store Case 2 roots directly as lookup tables, for those that only have a fractional part to the right of the decimal point. Roots of numbers starting with a long sequence of $0's$ or a long sequence of $9's$ can be solved for by many methods, the simplest being a taylor expansion.
Also, another option is to start from the top as:
$$\text{Case 2: } S^2 < 1 \implies S^2 = \left( \sqrt{\sum_{i={-({i_{min}})}}^{-(i_{max})} a_{-i} \times 10^{-i} } \right)^2$$
For $i_{max} \in \{\text{Even Integers}\}$ let $S^{'} = S \times 10^{\frac{i_{max}}{2}} \in \text{{Integers>0}}$.
For $i_{max} \in \{\text{Odd Integers}\}$ Let $S^{'} = S \times 10^{\frac{i_{max}+1}{2}} \in \text{{Integers>0}}$.
Then if we want to take the square root of $(S^{'})^2$ then we can proceed to label the digits $a_i \in \text{{Integers 0..9}}$ and $b_j \in \text{{Integers 0..9}}$ as follows (where $n$ represents the maximum power of 10 for the number $(S^{'})^2$ and $n/2=\text{Floor(}n/2\text{)}$:
$$(S^{'})^2 = \sum_{i=0}^n (a^{'})_i \times 10^i$$
$$S^{'} = \sum_{j=0}^{n/2} (b^{'})_j \times 10^j$$
Then, also we have:
$$(S^{'})^2= (\sum_{j=0}^{n/2} (b^{'})_j \times 10^j)^2=((b^{'})_{n/2})^2 \times 10^n + \text{ Terms with lower powers of 10}$$
We can then seek the maximum digit $(b^{'})_{n/2}$ such that $(((c^{'})_n)^2=((b^{'})_{n/2})^2) \le ((a^{'})_n)$ for $n \in \text{\{Set of Even Integers\}}$. Or if $n \in \text{\{Set of Odd Integers\}}$, we can then seek the maximum digit $(b^{'})_{n/2}$ such that $((b^{'})_{n/2})^2 \le ((a^{'})_n)*10+(a^{'})_{n-1}$ with here $((c^{'})_n)^2=((b^{'})_n)^2$.
Like with case one, after the lookup is complete, we can perform the relevant multiplication to recover the complete square root.
Finding the Square Root of 1 Plus a Very Small Remainder
We go ahead and devide $S^2$ as follows:
$${(S')}^2= \frac{S^2}{((c_n)^2 \times 10^n)}=\left(\sqrt{ {a_{0}} . + a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+...} \right)^2$$
It may be that $a_0 = 1, a_{-1} = 0, a_{-2}=0, ...., a_{-{(min-1)}}=0$ so that the resolution is below the resolution of the lookup table. To demonstrate this algorithm in that case, we present a hand calculation illustration. This is not to say that other algorithms could be used. Rather it is to show a simple approach to calculating the root that is efficient:
Let us start with a lookup table of five digits (which are functions of $1000 \times a_{0} 1000 \times a_{-1} + 100 \times a_{-2} + 10 \times a{-3} + a{-4}$) as $d_{-1}, d_{-2}, d_{-3}, \text{ and } d_{-4}$ such that:
$$\sqrt{a{0} . + a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+a_{-3} \times 10^{-3}+ a_{-4} \times 10^{-4}}$$ $$=1. + d_{-1} \times 10^{-1} + d_{-2} \times 10^{-2}+d_{-3} \times 10^{-3}+d_{-4} \times 10^{-4}$$
It may be that the root we seek needs additional precision beyond the lookup table. For instance, the square $(S^{'})^2$ might have non-zero terms until $a_{-5} \times 10^{-5}... $ and so forth. So how can the lookup table approach be expanded?
One approach: Taylor Expansion to Extend a Lookup Table's Precision
One approach (not to exclude others) is to use a Taylor Series expansion, with our lookup table being (for convenience in describing the procedure) an even number of digits.
Now we are seeking the root
$$\sqrt{1. + a_{-(i_{min})} \times 10^{-(i_{min})} + a_{-(i_{min+1})} \times 10^{-(i_{min+1})}+a_{-(i_{min+2})} \times 10^{-(i_{min+2})}+ a_{-(i_{min+3})} \times 10^{-(i_{min+3})}}$$
Case A: where $a_{-(i_{min})} \in \{ \text{1,2,3,4,5,6,7,8,9} \}$ corresponds to the first non-zero digit in the series of digits to the right of the decimal point, in the case where $i_{min} \in \{ \text{Set of Even Integers} \}$.
Quick example: $\sqrt{1.000001} \implies i_{min}=6$
Case B: where $a_{-(i_{min})} \in \{ \text{1,2,3,4,5,6,7,8,9} \}$ corresponds to the first non-zero digit in the series of digits to the right of the decimal point, in the case where $i_{min} \in \{ \text{Set of Odd Integers} \}$.
Quick example: $\sqrt{1.0000001} \implies i_{min}=7$
Both Case A and Case B: we will adjust $i_{min}$ as $i_{min\text{ }adj}$ so that the same calculations can be used in Case A and Case B by defining $i_{min\text{ }adj}$ as follows:
$$ \text{Case A: } i_{min} \in \{ \text{Set of Even Integers} \} \implies i_{min\text{ }adj}=i_{min}$$ $$ \text{Case B: } i_{min} \in \{ \text{Set of Odd Integers} \} \implies i_{min\text{ }adj}=i_{min}-1$$
It is now clear that $i_{min\text{ }adj} \in \{ \text{Set of Even Integers} \}$, so that $i_{min\text{ }adj}$ can be evenly dividable by 2.
Of course we previously had a lookup table of several even digits at the least. In this case we had already illustrated with a simple square root lookup table of 4 digits.
Now (starting from scratch with this special case of $1+\text{(A Very Small Remainder)}$ we notice that:
$$S^2 = 1 + \sum_{i=(i_{min\text{ }adj} )}^{i_{max}} a_{-i} \times 10^{-i}$$ $$S = 1 + \sum_{j=\frac{i_{min\text{ }adj}}{2}}^{j_{max}} b_{-j} \times 10^{-j}$$
$$\implies S = 1 + 10^{-\frac{i_{min\text{ }adj}}{2}} \times \sum_{j=\frac{i_{min\text{ }adj}}{2}}^{j_{max}} b_{-j} \times 10^{-j+\frac{i_{min\text{ }adj}}{2}}$$
$$\implies S^2 = \left( 1 + 10^{-\frac{i_{min\text{ }adj}}{2}} \times \sum_{j=\frac{i_{min\text{ }adj}}{2}}^{j_{max}} b_{-j} \times 10^{-j+\frac{i_{min\text{ }adj}}{2}} \right)^2$$
$$ = \left( 1 +2 \times 10^{-\frac{i_{min\text{ }adj}}{2}} \times \sum_{j=\frac{i_{min\text{ }adj}}{2}}^{j_{max}} b_{-j} \times 10^{-j+\frac{i_{min\text{ }adj}}{2}} + 10^{-i_{min\text{ }adj}} \times \left( \sum_{j=\frac{i_{min\text{ }adj}}{2}}^{j_{max}} b_{-j} \times 10^{-j+\frac{i_{min\text{ }adj}}{2}} \right)^2 \right)$$
And thus:
$$S^2 = 1 + \sum_{i=(i_{min\text{ }adj} )}^{i_{max}} a_{-i} \times 10^{-i}$$ $$ = \left( 1 +2 \times 10^{-\frac{i_{min\text{ }adj}}{2}} \times \sum_{j=\frac{i_{min\text{ }adj}}{2}}^{j_{max}} b_{-j} \times 10^{-j+\frac{i_{min\text{ }adj}}{2}} + 10^{-i_{min\text{ }adj}} \times \left( \sum_{j=\frac{i_{min\text{ }adj}}{2}}^{j_{max}} b_{-j} \times 10^{-j+\frac{i_{min\text{ }adj}}{2}} \right)^2 \right)$$
And:
$$S^2-1=\sum_{i=(i_{min\text{ }adj} )}^{i_{max}} a_{-i} \times 10^{-i}$$ $$ = \left( 2 \times 10^{-\frac{i_{min\text{ }adj}}{2}} \times \sum_{j=\frac{i_{min\text{ }adj}}{2}}^{j_{max}} b_{-j} \times 10^{-j+\frac{i_{min\text{ }adj}}{2}} + 10^{-i_{min\text{ }adj}} \times \left( \sum_{j=\frac{i_{min\text{ }adj}}{2}}^{j_{max}} b_{-j} \times 10^{-j+\frac{i_{min\text{ }adj}}{2}} \right)^2 \right)$$
Example 1 (Case 1): Find the approximate square root of $S^2=1.00000008$.
We have $S^2-1=0.00000008=8 \times 10^{-8}$ that we need a solution for. The approximate root is $S-1 \approx 10^{-8} \times 8 \div 2 = 4 \times 10^{-8}$ so we have our approximate solution as $S \approx 1+ 4 \times 10^{-8}=1.00000004$. Indeed the Ubuntu Linux calculator finds that $√1.00000008=1.00000004$. Of course we know from our formula that this is an approximation. However, it is as accurate as the Ubuntu Calculator application takes the root so the approximation is fairly good. The error is within $\approx 8^2 \times 10^{-16}$, the term at the right in the Taylor expansion of the sequences.
Once we have our approximate root, we can of course take it to the left of the equation, so we come back to either something within our lookup table or else like above something that has $1$ with several $0's$ to the right of the decimal point before the digit sequence, making another approximation easier.
Example 2 (Case 2): Find the approximate square root of $S^2=0.99999992$.
We have $S^2-1=-0.00000008=-8 \times 10^{-8}$ that we need a solution for. The approximate root is $S-1 \approx -1 \times 10^{-8} \times 8 \div 2 = -4 \times 10^{-8}$ so we have our approximate solution as $S \approx 1- 4 \times 10^{-8}=0.99999996$. Indeed the Ubuntu Linux calculator finds that $√0.99999992=0.99999996$. Of course we know from our formula that this is an approximation. However, it is as accurate as the Ubuntu Calculator application takes the root so the approximation is fairly good. The error is within $\approx 8^2 \times 10^{-16}$, the term at the right in the Taylor expansion of the sequences.
Final Case 1 Calculations
In the begining of this Answer, we took out factors:
$${(S')}^2=\frac{S^2}{((c_n)^2 \times 10^n)}=\left( { \sqrt{{a_0}=1 . + a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+...} } \right) ^2$$
to get our root into a generic form of (as in the example to be generalized) as:
Let us start with a lookup table of five digits (which are functions of $10000 \times a_{0} + 1000 \times a_{-1} + 100 \times a_{-2} + 10 \times a{-3} + a{-4}$) as $d_{-1}, d_{-2}, d_{-3}, \text{ and } d_{-4}$ such that:
$$\sqrt{a{0}. + a_{-1} \times 10^{-1} + a_{-2} \times 10^{-2}+a_{-3} \times 10^{-3}+ a_{-4} \times 10^{-4}}$$ $$=d_{0}. + d_{-1} \times 10^{-1} + d_{-2} \times 10^{-2}+d_{-3} \times 10^{-3}+d_{-4} \times 10^{-4}$$
Once we conclude our talbe lookups and other calculatons, we again multiply with what we divided out to generalize the technique to get the complete answer.
Basic Calculation Examples
Although these routines are very efficient, concrete examples are of course helpful to make this answer clearer. Time available, I hope to update the answer with this more concrete detail. Right now, though, this answer should be a complete description of how to take the root efficiently, at least as good as the calculator applications on modern computers.
Concrete Calculation Examples √3
First, we observe that $3 \text{ is a Prime Number!} $ So $3 \text{ does not have a rational root!}$ However, we can try to approximate the square root of three by finding a square root sequence that is close, either a little too great or a little to small. As the error in taking the square root becomes very small, the taylor approximations become more and more constrained, a little bit too positive or a little bit too negative, never exactly accurate. Thus, the best we can hope for in taking $√3$ are closer and closer approximations and limiting bounds.
$$ \sqrt{3}= \sqrt{ ((2 \times 2)=4) \times \frac{3}{4}} $$
We can write out approximations as follows (getting the approximate root of $0.8660^2=0.784996$ from the lookup table):
$$ \sqrt{3} =\sqrt{(0.8660^2=0.784996) \times ((2 \times 2)=4) \times \frac{3}{4} \frac{1}{0.784996}} $$
$$ \sqrt{3} = 0.8660 \times 2 \times \sqrt{\frac{3}{4} \frac{1}{0.784996}} $$
When we have sufficient accuracy, we can then move the terms inside of the right-most root to the left as follows (understanding that because 3 is not an integer square, we can only approximate a root with rational numbers or floating point numbers):
$$\sqrt{3} \approx 0.8660 \times 2 \approx 1.732 $$ $$\left( \sqrt{3} \right)^2 \approx {1.732}^2 \approx 2.999824 \text{ This square is within −0.0059% error.}$$
After that we can use taylor series expansions to narrow the root down even futher, since the root estimate is already very close with just these two first steps. The details how to accomplish the taylor expansion are still something that I am actively writing and checking for accurary before updating this answer here with that additional concrete detail. The theoretical detail was given earlier in this answer.
Accuracy for √3 With a Full Precision Calculation
The Full Precision Calculator On-Line can be used to approximate the $\sqrt{3}$ as $$\text{"sqrt(3)" } \approx $$ $$\left( 1.732050807568877293527446341505872366942805253810380628055806979451933016908800037081146186757248575675626141415406703029969945094998952478811655512094373648528093231902305582067974820101084674923265 \right) $$ and $$\left( \text{"sqrt(3)" }\right)^2 = $$ $$\left( 2.9999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999 \right)$$ The approximation that we have of $\sqrt{3} \approx 1.732$ is correct to three digits after the decimal point, using just the initial division, the lookup-table calculation, and the final multiplications, a very efficient process.
Continue With Full Precision Taylor Expansion for Calculation for √3
We could at this point calculate additional terms with software code, but there is no point getting detracted by the limited resolution offered by the C-language (long double) for instance, just now. The limitation in resolution of the Ubuntu Calculator appication also is not being used because of its also limited resolution. Eventually a high-precision data type is needed, but that is left to program implementation.
Let us discuss for a bit the algorithm and how it had proceeded.
We had a value $S^2$ that we do not know what the value of $S$ is.
However, we start with approximations to narrow down on the value of $S$.
For instance, $S_1=2 \text{ and } \left(S_1\right)^2=4$.
Now, since we know exactly the root of $\left(S_1\right)^2=4$, we can then seek from our lookup table the next iteration in calculating the root $\sqrt{3}$.
Our lookup table has that $\sqrt{\frac{3}{4}} \approx 0.8660$ and that $0.8660^2 = 0.749956 $ so $S_2=0.8660$ and exactly $\left(S_2\right)^2=0.749956$ using the full precision calculator. So we have:
$$S \approx S_1 \times S_2$$
We want to go forward with a Taylor series approximation for $S_3$ and $(S_3)^2$. We desire that:
$$\left( \frac{3}{4 \times 0.749956} \right) \times \frac{1}{\left( S_3 \right)^{2}} \approx 1 $$
We have from the Full Precision calculator that:
$$\left( S_3 \right)^{2} \approx \left( \frac{3}{4 \times 0.749956} \right) $$ $$\approx \left( 1.00005867010864637392 \right) $$ $$\text{ already } \approx 1$$
But we want a better reprensentation. So then - already - the $((S_1)^2=4) \times ((S_2)^2=0.784996) \text{ is a fairly good representation of 3. }$ And we are hoping to increase the accuracy even furhter with the Taylor series.
That is so that at the end of determining the Taylor approximate $S_3$ and from there calculate exactly $\left(S_3 \right)^2 $; then we can at the end of the calculation multipy $S_1$ and $S_2$ by $S_3$ so that $S \approx S_1 \times S_2 \times S_3 $ and determine the accuracy of the root.
If we approximate $(S_3)^{2}= 1.00005867010864637392$ by taking the Taylor expansion, we have that: $$1.00005867010864637392-1 = $$ $$=0.00005867010864637392$$ And $$\frac{ \left( 0.00005867010864637392 \right) }{2}=$$ $$=0.00002933505432318696$$
So we have our approximation for $(S_3)$ as:
$$S_3 \approx$$ $$\approx 1.00002933505432318696$$
Now we will find an exact value for $(S_3)^2$ using the Full Precision Calculator; We have:
$$(S_3)^2=$$ $$=1.00002933505432318696^2 $$ $$=1.000058670969191786064329951841491114 $$ $$\text{ to a Limited Number of Digits}$$
And the technique seems to have somewhat succeeded. Note that we are increasing precision of the calculation for $\sqrt{S} = \prod_1^{\left(3 +... \right)} {\left( S_i \right)}$, and now with each iteration of the Taylor series, we are able to obtain lower and upper bounds for the root. This is especially handy if we want to determine if there is an integer square, since if the lower and upper bounds of $S_{min}$ and $S_{max}$ are the same integer and a different small fraction, there then is no integer root and rather the root has either a fractional part or is irrational.
Now we continue verifying our calculations. We had:
$$\left( \frac{3}{ \left( S_1 \right)^2 \times \left( S_1 \right)^2} \right) \times $$ $$\times \frac{1}{\left( S_3 \right)^{2}} \approx 1 $$
$$\left( \frac{3}{4 \times 0.749956} \right) \times $$ $$ \frac{1}{\left( \left( S_3 \right)^{2} = 1.000058670969191786064329951841491114 \right)} $$ $$ = \frac{3}{4 \times 0.749956 \times 1.000058670969191786064329951841491114}$$ $$ = 0.999999999139505$$ $$ \approx 1$$
And
$$\left( \frac{3}{4 \times 0.749956} \right) \times $$ $$ \frac{1}{\left( 1.00002933505432318696 \right)^{2}} \approx 1 $$
We stored away from estimated $S_1, S_2, S_3$ the exactly derived $\left( S_1 \right)^2, \left( S_2 \right)^2 $ and $\left( S_3 \right)^2 $.
Again, using the Full Precision Calculator we should have $\sqrt{3} \approx $ (to three product terms) as follows:
$$\sqrt{3} \approx \left( S_1 \right) \times \left( S_2 \right) \times \left( S_3 \right) \approx $$
$$\approx 2 \times 0.8660 \times $$ $$\times 1.00002933505432318696 = $$ $$=1.73205080831408775981472$$ And $$\left( \approx \sqrt{3} \right)^2 =$$ $$= 3.0000000025814847804386425334529492396887286784$$ With just three multiplicative terms $S_1, S_2, and S_3$ we now have a Relative Error of: $$\text{Relative Error %} = ( 3.0000000025814847804386425334529492396887286784-$$ $$-3 ) /3 \times 100 \%$$ $$= 0.0000000860494926812880844484316413229576226 \%$$ This is a very efficient and accurate algorithm here, and it even lets us know that there is is no perfect integer root!
- 175
-
It is really hard with an ordinary cut and paste into word-processing applications to print all of the answers here. However, using Print to PDF with the latest Firefox browser download at the stackprinter.com site has excellent results to clearly view all of the answers without any formatting difficulties. Thanks stackprinter.com for all your help! And for those thinking it is my site - it's not! I am just really happy with its results to recommend it and hopefully soon donate something if I can! – Stephen Elliott Aug 02 '23 at 16:41
You can use the Binomial Expansion to calculate square roots and other fractional powers. Since most square roots are irrational you get an infinite series, which you truncate to give you whatever level of accuracy/precision you require.
The Binomial Expansion only converges for |x| < 1 , so to get round this restriction take out a factor p where p is the largest perfect square ( or cube and so on... ) below s the number whose fractional index you are trying to evaluate.
For example for s = 5 the largest perfect square below it is 4 so
√5 = √4*√(1+1/4) here x is 1/4 satisfying the condition for convergence.
The first few terms for √5 are:
√5 = √4[1 + 1/2*1/4 - 1/2*(1/4)²/2 + …… ]
√5 = 2(1 + 1/8 - 1/64 + …… )
- 193