[gfortran] Compute x**(-n) as (1/x)**n

[gfortran] Compute x**(-n) as (1/x)**n

[Toon Moene]

[-- Attachment #1: Type: text/plain, Size: 1056 bytes --]

The following change, which I bootstrapped on powerpc-unknown-linux-gnu 
(C and f95) without regressions, leads the Fortran front-end to generate 
expressions x**(-n) with n constant as if they were written (1/x)**n.

This will enable more opportunities for CSE in case multiple powers of 
the same base are needed in an expression.  E.g., the following:

function force(d)
force = d**(-2) + d**(-6) + d**(-8) + d**(-12)
end

will lead to the following floating point operations (when using 
-ffast-math, whereas the current method needs four divisions and uses no 
fused multiply-add.

         fdivs 1,1,13
         fmuls 0,1,1
         fmuls 1,1,0
         fmuls 13,0,0
         fmuls 1,1,1
         fadds 0,0,1
         fmadds 13,13,13,0
         fmadds 1,1,1,13

Paul, is this OK ?

-- 
Toon Moene - mailto:toon@moene.indiv.nluug.nl - phoneto: +31 346 214290
Saturnushof 14, 3738 XG  Maartensdijk, The Netherlands
Maintainer, GNU Fortran 77: http://gcc.gnu.org/onlinedocs/g77_news.html
GNU Fortran 95: http://gcc.gnu.org/fortran/ (under construction)

[-- Attachment #2: rexpnegint.patch --]
[-- Type: text/plain, Size: 917 bytes --]

2004-06-08  Toon Moene  <toon@moene.indiv.nluug.nl>

	* trans-expr.c (gfc_conv_cst_int_power): Compute
	x**(-n) by converting it to (1/x)**n instead of
	1/x**n.

*** trans-expr.c.orig	Sat Jun  5 13:45:07 2004
--- trans-expr.c	Sat Jun  5 13:48:38 2004
*************** gfc_conv_cst_int_power (gfc_se * se, tre
*** 526,538 ****
  
    memset (vartmp, 0, sizeof (vartmp));
    vartmp[1] = lhs;
- 
-   se->expr = gfc_conv_powi (se, n, vartmp);
    if (sgn == -1)
      {
        tmp = gfc_build_const (type, integer_one_node);
!       se->expr = build (RDIV_EXPR, type, tmp, se->expr);
      }
    return 1;
  }
  
--- 526,539 ----
  
    memset (vartmp, 0, sizeof (vartmp));
    vartmp[1] = lhs;
    if (sgn == -1)
      {
        tmp = gfc_build_const (type, integer_one_node);
!       vartmp[1] = build (RDIV_EXPR, type, tmp, vartmp[1]);
      }
+ 
+   se->expr = gfc_conv_powi (se, n, vartmp);
+ 
    return 1;
  }
  

Re: [gfortran] Compute x**(-n) as (1/x)**n

[Tobias Schlüter]

Toon Moene wrote:
> The following change, which I bootstrapped on powerpc-unknown-linux-gnu 
> (C and f95) without regressions, leads the Fortran front-end to generate 
> expressions x**(-n) with n constant as if they were written (1/x)**n.
> 
Does this deal correctly with integer x? I remember there were some 
special cases somewhere.

- Tobi

Re: [gfortran] Compute x**(-n) as (1/x)**n

[Toon Moene]

Tobias Schlüter wrote:

> Toon Moene wrote:
> 
>> The following change, which I bootstrapped on 
>> powerpc-unknown-linux-gnu (C and f95) without regressions, leads the 
>> Fortran front-end to generate expressions x**(-n) with n constant as 
>> if they were written (1/x)**n.
>>
> Does this deal correctly with integer x? I remember there were some 
> special cases somewhere.

It doesn't have to :-)  By the time control reaches the point where I 
made my change, INTEGER x has already been dealt with.

Hope this helps,

-- 
Toon Moene - mailto:toon@moene.indiv.nluug.nl - phoneto: +31 346 214290
Saturnushof 14, 3738 XG  Maartensdijk, The Netherlands
Maintainer, GNU Fortran 77: http://gcc.gnu.org/onlinedocs/g77_news.html
GNU Fortran 95: http://gcc.gnu.org/fortran/ (under construction)

Re: [gfortran] Compute x**(-n) as (1/x)**n

[Paul Brook]

On Tuesday 08 June 2004 18:33, Toon Moene wrote:
> The following change, which I bootstrapped on powerpc-unknown-linux-gnu
> (C and f95) without regressions, leads the Fortran front-end to generate
> expressions x**(-n) with n constant as if they were written (1/x)**n.
>
> This will enable more opportunities for CSE in case multiple powers of
> the same base are needed in an expression.  E.g., the following:
>
> function force(d)
> force = d**(-2) + d**(-6) + d**(-8) + d**(-12)
> end

Well, it also has the reverse effect on the following code:
force = d**(-2) + d**2

I guess this is probably a less common case, so your patch is ok.
The real solution is to teach the optimizers how to do this transformation 
themselves. I'm not offering to do that though ;) 

Paul