# Why does squaring generate a different result?

**URL:** <https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578>\
**Category:** Help\
**Created:** [April 15, 2023, 6:38pm UTC](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578 "2023-04-15T18:38:36Z")\
**Posts on this page:** 11\
**Page:** 2

<div class="post-metadata">

**Author:** ![garynewport](https://avatars.discourse-cdn.com/v4/letter/g/898d66/32.png) [@garynewport](https://fortran-lang.discourse.group/u/garynewport)\
**Post date:** [April 17, 2023, 8:46pm UTC](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578/21 "2023-04-17T20:46:00Z")

</div>

Oh no, I get that. It was specific to the real\*real and real\*\*2 (and other powers).

I have found several times that ‘simplifying’ the maths results in some ‘interesting’ results.

---

<div class="post-metadata">

**Author:** ![msz59](https://avatars.discourse-cdn.com/v4/letter/m/3d9bf3/32.png) [@msz59](https://fortran-lang.discourse.group/u/msz59)\
**Post date:** [April 17, 2023, 9:05pm UTC](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578/22 "2023-04-17T21:05:10Z")

</div>

> [@sblionel](#):
>
> I’d expect most compilers to perform real\*\*2 by doing real\*real inline.

Isn’t that required by the standard? If I remember well, _real\*\*integer_ used to be treated as an exception to the general rule of converting integer to real in any other _real (op) intreger_ expression, just because it should be done as repeated multiplication (plus one reciprocal in case of negative power), instead of taking logarithm, multiplication and exponent.

---

<div class="post-metadata">

**Author:** ![sblionel](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/sblionel/32/853_2.png) [@sblionel](https://fortran-lang.discourse.group/u/sblionel)\
**Post date:** [April 17, 2023, 9:20pm UTC](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578/23 "2023-04-17T21:20:58Z")

</div>

> [@msz59](#):
>
> Isn’t that required by the standard?

No. All the standard says is that the interpretation of x1\*\*x2 is “Raise x1 to the power x2” and, “When a value of type real or complex is raised to an integer power, the integer operand need not be converted.”

It would be standard-conforming to convert the integer power to a real and then use logarithms to compute the result. I don’t know of any compiler that does this.

---

<div class="post-metadata">

**Author:** ![RonShepard](https://avatars.discourse-cdn.com/v4/letter/r/a3d4f5/32.png) [@RonShepard](https://fortran-lang.discourse.group/u/RonShepard)\
**Post date:** [April 17, 2023, 10:29pm UTC](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578/24 "2023-04-17T22:29:43Z")

</div>

> [@msz59](#):
>
> > [@sblionel](#):
> >
> > I’d expect most compilers to perform real\*\*2 by doing real\*real inline.
> 
> Isn’t that required by the standard?

I think what @sblionel meant is that the optimization performed is to inline the multiplication, as if it were written as `x*x` in the fortran code in the first place, rather than to call a general procedure that computes the value. Some compilers do this inline optimization, i.e, avoid the procedure call overhead, for integer powers up to 8 or so. For example, `x**8` might be computed as `t=x*x; t=t*t; t=t*t;` with only three multiplications and no procedure call.

That general procedure most likely has two separate algorithms within it. For small integer powers, the result is computed recursively by looking at the bits of the integer exponent. When a bit is 1 then two multiplications are done (one for the recursive squaring of the intermediate, the other to accumulate the result), and when the bit is 0, then only one multiplication is done (for the recursive squaring). That is what that inline `x**8` code does above. When the integer power is large enough (i.e. few leading zeros), then it is cheaper to compute the result as `exp(k*log(x))` by scaling `x` to be near 1 (so the `log(x)` expansion converges quickly) and scaling `k` accordingly.

---

<div class="post-metadata">

**Author:** ![mecej4](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/mecej4/32/855_2.png) [@mecej4](https://fortran-lang.discourse.group/u/mecej4)\
**Post date:** [April 17, 2023, 10:41pm UTC](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578/25 "2023-04-17T22:41:19Z")

</div>

> [@garynewport](#):
>
> There is a distinct possibility that there is no difference and that an error elsewhere was just being brought out by the use of squaring.

There is also a possibility, however slight, that not squaring made a real and significant error go into hiding. In other words, the changes made to the code could be the equivalent of converting “bad news” into “no news” by jailing the messenger.

---

<div class="post-metadata">

**Author:** ![sblionel](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/sblionel/32/853_2.png) [@sblionel](https://fortran-lang.discourse.group/u/sblionel)\
**Post date:** [April 18, 2023, 12:03am UTC](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578/27 "2023-04-18T00:03:59Z")

</div>

> [@kargl](#):
>
> Yes, but what @sblionel failed to mention but obviously knows, is that `a = b * c**2` is likely inlined as `a = b * (c * c)`. The `**` operator has higher precedence than `*` and must be evaluated first

Oh, but I did mention it!

> [@sblionel](#):
>
> However, you must also consider operator precedence (see [Doctor Fortran in “Order! Order!” - Doctor Fortran (stevelionel.com)](https://stevelionel.com/drfortran/2021/04/03/doctor-fortran-in-order-order/)) Since exponentiation has higher precedence than multiplication, it will always be done first

---

<div class="post-metadata">

**Author:** ![RonShepard](https://avatars.discourse-cdn.com/v4/letter/r/a3d4f5/32.png) [@RonShepard](https://fortran-lang.discourse.group/u/RonShepard)\
**Post date:** [April 18, 2023, 12:12am UTC](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578/28 "2023-04-18T00:12:53Z")

</div>

> [@sblionel](#):
>
> Oh, but I did mention it!

I think what @kargl meant is that when the ` **` expression is inlined, it is inlined** with** the parentheses in order to maintain the same evaluation order. The original expression with multiplication was without the parentheses, thus the original expression and the inlined version might not necessarily agree since, as you pointed out, the compiler is free to evaluate it in either order in one case, but a specific order in the other.

---

<div class="post-metadata">

**Author:** ![mecej4](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/mecej4/32/855_2.png) [@mecej4](https://fortran-lang.discourse.group/u/mecej4)\
**Post date:** [April 18, 2023, 12:39am UTC](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578/29 "2023-04-18T00:39:39Z")

</div>

There is a type of expression that some compilers allow as an extension to the standard, and compile without warning unless asked to, that may have been used in the simulation code that @garynewport is using. I have not come across an explanation of why this extension exists nor have read an explanation of when it is useful. Here is an example:

```auto
program ffcub
   implicit none
   real x,y,z
   x = 12.0
   y = 5.0
   z = x **3+-y** 3
   print *,x,y,z
end

```

---

<div class="post-metadata">

**Author:** ![RonShepard](https://avatars.discourse-cdn.com/v4/letter/r/a3d4f5/32.png) [@RonShepard](https://fortran-lang.discourse.group/u/RonShepard)\
**Post date:** [April 18, 2023, 5:15am UTC](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578/33 "2023-04-18T05:15:20Z")

</div>

> [@sblionel](#):
>
> (Some compilers will not honor the parentheses by default, however.)

Could you say which compilers?

---

<div class="post-metadata">

**Author:** ![garynewport](https://avatars.discourse-cdn.com/v4/letter/g/898d66/32.png) [@garynewport](https://fortran-lang.discourse.group/u/garynewport)\
**Post date:** [April 18, 2023, 6:59pm UTC](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578/34 "2023-04-18T18:59:16Z")

</div>

Thanks everyone for your comments.  
As I stated in my long reply, I had another instance of this where it didn’t make small variations in the data but caused the program to crash. I, therefore, decided to test this out by making those changes - and everything worked fine!

In fact, not simply fine but the code is running on average faster.

The good news is that I am almost finished with this model now. After 2 years of adjustments, updates, rebuilds and some interesting stressful times, the model is stable and fast. I have just run through checking for unused variables, checked off all the warnings and run maximum optimisation. Doing a test on my Mac, then will test under Windows before accepting that I have a working model for the basic stages. I will then begin testing for long-term models and, hopefully, if this works, I will consider the system finished.

Now onto the interface! Woohoo!

---

<div class="post-metadata">

**Author:** ![sblionel](https://yyz2.discourse-cdn.com/free1/user_avatar/fortran-lang.discourse.group/sblionel/32/853_2.png) [@sblionel](https://fortran-lang.discourse.group/u/sblionel)\
**Post date:** [April 19, 2023, 12:39am UTC](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578/35 "2023-04-19T00:39:27Z")

</div>

> [@RonShepard](#):
>
> Could you say which compilers?

ifort, for one, will associate across parens unless you add `-assume protect-parens`. The Fortran team argued with the back-end team a long time about this. I don’t know if the LLVM-based ifx is any different.

[Previous page](https://fortran-lang.discourse.group/t/why-does-squaring-generate-a-different-result/5578.md?page=1)
