crypto/bn/asm/ppc64-mont.pl

   1 #!/usr/bin/env perl
   2
   3 # ====================================================================
   4 # Written by Andy Polyakov <appro@fy.chalmers.se> for the OpenSSL
   5 # project. The module is, however, dual licensed under OpenSSL and
   6 # CRYPTOGAMS licenses depending on where you obtain it. For further
   7 # details see http://www.openssl.org/~appro/cryptogams/.
   8 # ====================================================================
   9
  10 # December 2007
  11
  12 # The reason for undertaken effort is basically following. Even though
  13 # Power 6 CPU operates at incredible 4.7GHz clock frequency, its PKI
  14 # performance was observed to be less than impressive, essentially as
  15 # fast as 1.8GHz PPC970, or 2.6 times(!) slower than one would hope.
  16 # Well, it's not surprising that IBM had to make some sacrifices to
  17 # boost the clock frequency that much, but no overall improvement?
  18 # Having observed how much difference did switching to FPU make on
  19 # UltraSPARC, playing same stunt on Power 6 appeared appropriate...
  20 # Unfortunately the resulting performance improvement is not as
  21 # impressive, ~30%, and in absolute terms is still very far from what
  22 # one would expect from 4.7GHz CPU. There is a chance that I'm doing
  23 # something wrong, but in the lack of assembler level micro-profiling
  24 # data or at least decent platform guide I can't tell... Or better
  25 # results might be achieved with VMX... Anyway, this module provides
  26 # *worse* performance on other PowerPC implementations, ~40-15% slower
  27 # on PPC970 depending on key length and ~40% slower on Power 5 for all
  28 # key lengths. As it's obviously inappropriate as "best all-round"
  29 # alternative, it has to be complemented with run-time CPU family
  30 # detection. Oh! It should also be noted that unlike other PowerPC
  31 # implementation IALU ppc-mont.pl module performs *suboptimaly* on
  32 # >=1024-bit key lengths on Power 6. It should also be noted that
  33 # *everything* said so far applies to 64-bit builds! As far as 32-bit
  34 # application executed on 64-bit CPU goes, this module is likely to
  35 # become preferred choice, because it's easy to adapt it for such
  36 # case and *is* faster than 32-bit ppc-mont.pl on *all* processors.
  37
  38 # February 2008
  39
  40 # Micro-profiling assisted optimization results in ~15% improvement
  41 # over original ppc64-mont.pl version, or overall ~50% improvement
  42 # over ppc.pl module on Power 6. If compared to ppc-mont.pl on same
  43 # Power 6 CPU, this module is 5-150% faster depending on key length,
  44 # [hereafter] more for longer keys. But if compared to ppc-mont.pl
  45 # on 1.8GHz PPC970, it's only 5-55% faster. Still far from impressive
  46 # in absolute terms, but it's apparently the way Power 6 is...
  47
  48 $flavour = shift;
  49
  50 if ($flavour =~ /32/) {
  51         $SIZE_T=4;
  52         $RZONE= 224;
  53         $FRAME= $SIZE_T*12+8*12;
  54         $fname= "bn_mul_mont_ppc64";
  55
  56         $STUX=  "stwux";        # store indexed and update
  57         $PUSH=  "stw";
  58         $POP=   "lwz";
  59         die "not implemented yet";
  60 } elsif ($flavour =~ /64/) {
  61         $SIZE_T=8;
  62         $RZONE= 288;
  63         $FRAME= $SIZE_T*12+8*12;
  64         $fname= "bn_mul_mont";
  65
  66         # same as above, but 64-bit mnemonics...
  67         $STUX=  "stdux";        # store indexed and update
  68         $PUSH=  "std";
  69         $POP=   "ld";
  70 } else { die "nonsense $flavour"; }
  71
  72 $0 =~ m/(.*[\/\\])[^\/\\]+$/; $dir=$1;
  73 ( $xlate="${dir}ppc-xlate.pl" and -f $xlate ) or
  74 ( $xlate="${dir}../../perlasm/ppc-xlate.pl" and -f $xlate) or
  75 die "can't locate ppc-xlate.pl";
  76
  77 open STDOUT,"| $^X $xlate $flavour ".shift || die "can't call $xlate: $!";
  78
  79 $FRAME=($FRAME+63)&~63;
  80 $TRANSFER=16*8;
  81
  82 $carry="r0";
  83 $sp="r1";
  84 $toc="r2";
  85 $rp="r3";       $ovf="r3";
  86 $ap="r4";
  87 $bp="r5";
  88 $np="r6";
  89 $n0="r7";
  90 $num="r8";
  91 $rp="r9";       # $rp is reassigned
  92 $tp="r10";
  93 $j="r11";
  94 $i="r12";
  95 # non-volatile registers
  96 $nap_d="r14";   # interleaved ap and np in double format
  97 $a0="r15";      # ap[0]
  98 $t0="r16";      # temporary registers
  99 $t1="r17";
 100 $t2="r18";
 101 $t3="r19";
 102 $t4="r20";
 103 $t5="r21";
 104 $t6="r22";
 105 $t7="r23";
 106
 107 # PPC offers enough register bank capacity to unroll inner loops twice
 108 #
 109 #     ..A3A2A1A0
 110 #           dcba
 111 #    -----------
 112 #            A0a
 113 #           A0b
 114 #          A0c
 115 #         A0d
 116 #          A1a
 117 #         A1b
 118 #        A1c
 119 #       A1d
 120 #        A2a
 121 #       A2b
 122 #      A2c
 123 #     A2d
 124 #      A3a
 125 #     A3b
 126 #    A3c
 127 #   A3d
 128 #    ..a
 129 #   ..b
 130 #
 131 $ba="f0";       $bb="f1";       $bc="f2";       $bd="f3";
 132 $na="f4";       $nb="f5";       $nc="f6";       $nd="f7";
 133 $dota="f8";     $dotb="f9";
 134 $A0="f10";      $A1="f11";      $A2="f12";      $A3="f13";
 135 $N0="f14";      $N1="f15";      $N2="f16";      $N3="f17";
 136 $T0a="f18";     $T0b="f19";
 137 $T1a="f20";     $T1b="f21";
 138 $T2a="f22";     $T2b="f23";
 139 $T3a="f24";     $T3b="f25";
 140 \f
 141 # sp----------->+-------------------------------+
 142 #               | saved sp                      |
 143 #               +-------------------------------+
 144 #               |                               |
 145 #               +-------------------------------+
 146 #               | 10 saved gpr, r14-r23         |
 147 #               .                               .
 148 #               .                               .
 149 #   +12*size_t  +-------------------------------+
 150 #               | 12 saved fpr, f14-f25         |
 151 #               .                               .
 152 #               .                               .
 153 #   +12*8       +-------------------------------+
 154 #               | padding to 64 byte boundary   |
 155 #               .                               .
 156 #   +X          +-------------------------------+
 157 #               | 16 gpr<->fpr transfer zone    |
 158 #               .                               .
 159 #               .                               .
 160 #   +16*8       +-------------------------------+
 161 #               | __int64 tmp[-1]               |
 162 #               +-------------------------------+
 163 #               | __int64 tmp[num]              |
 164 #               .                               .
 165 #               .                               .
 166 #               .                               .
 167 #   +(num+1)*8  +-------------------------------+
 168 #               | padding to 64 byte boundary   |
 169 #               .                               .
 170 #   +X          +-------------------------------+
 171 #               | double nap_d[4*num]           |
 172 #               .                               .
 173 #               .                               .
 174 #               .                               .
 175 #               +-------------------------------+
 176 \f
 177 $code=<<___;
 178 .machine "any"
 179 .text
 180
 181 .globl  .$fname
 182 .align  5
 183 .$fname:
 184         cmpwi   $num,4
 185         mr      $rp,r3          ; $rp is reassigned
 186         li      r3,0            ; possible "not handled" return code
 187         bltlr-
 188         andi.   r0,$num,1       ; $num has to be even
 189         bnelr-
 190
 191         slwi    $num,$num,3     ; num*=8
 192         li      $i,-4096
 193         slwi    $tp,$num,2      ; place for {an}p_{lh}[num], i.e. 4*num
 194         add     $tp,$tp,$num    ; place for tp[num+1]
 195         addi    $tp,$tp,`$FRAME+$TRANSFER+8+64+$RZONE`
 196         subf    $tp,$tp,$sp     ; $sp-$tp
 197         and     $tp,$tp,$i      ; minimize TLB usage
 198         subf    $tp,$sp,$tp     ; $tp-$sp
 199         $STUX   $sp,$sp,$tp     ; alloca
 200
 201         $PUSH   r14,`2*$SIZE_T`($sp)
 202         $PUSH   r15,`3*$SIZE_T`($sp)
 203         $PUSH   r16,`4*$SIZE_T`($sp)
 204         $PUSH   r17,`5*$SIZE_T`($sp)
 205         $PUSH   r18,`6*$SIZE_T`($sp)
 206         $PUSH   r19,`7*$SIZE_T`($sp)
 207         $PUSH   r20,`8*$SIZE_T`($sp)
 208         $PUSH   r21,`9*$SIZE_T`($sp)
 209         $PUSH   r22,`10*$SIZE_T`($sp)
 210         $PUSH   r23,`11*$SIZE_T`($sp)
 211         stfd    f14,`12*$SIZE_T+0`($sp)
 212         stfd    f15,`12*$SIZE_T+8`($sp)
 213         stfd    f16,`12*$SIZE_T+16`($sp)
 214         stfd    f17,`12*$SIZE_T+24`($sp)
 215         stfd    f18,`12*$SIZE_T+32`($sp)
 216         stfd    f19,`12*$SIZE_T+40`($sp)
 217         stfd    f20,`12*$SIZE_T+48`($sp)
 218         stfd    f21,`12*$SIZE_T+56`($sp)
 219         stfd    f22,`12*$SIZE_T+64`($sp)
 220         stfd    f23,`12*$SIZE_T+72`($sp)
 221         stfd    f24,`12*$SIZE_T+80`($sp)
 222         stfd    f25,`12*$SIZE_T+88`($sp)
 223
 224         ld      $a0,0($ap)      ; pull ap[0] value
 225         ld      $n0,0($n0)      ; pull n0[0] value
 226         ld      $t3,0($bp)      ; bp[0]
 227
 228         addi    $tp,$sp,`$FRAME+$TRANSFER+8+64`
 229         li      $i,-64
 230         add     $nap_d,$tp,$num
 231         and     $nap_d,$nap_d,$i        ; align to 64 bytes
 232 \f
 233         mulld   $t7,$a0,$t3     ; ap[0]*bp[0]
 234         ; nap_d is off by 1, because it's used with stfdu/lfdu
 235         addi    $nap_d,$nap_d,-8
 236         srwi    $j,$num,`3+1`   ; counter register, num/2
 237         mulld   $t7,$t7,$n0     ; tp[0]*n0
 238         addi    $j,$j,-1
 239         addi    $tp,$sp,`$FRAME+$TRANSFER-8`
 240         li      $carry,0
 241         mtctr   $j
 242
 243         ; transfer bp[0] to FPU as 4x16-bit values
 244         extrdi  $t0,$t3,16,48
 245         extrdi  $t1,$t3,16,32
 246         extrdi  $t2,$t3,16,16
 247         extrdi  $t3,$t3,16,0
 248         std     $t0,`$FRAME+0`($sp)
 249         std     $t1,`$FRAME+8`($sp)
 250         std     $t2,`$FRAME+16`($sp)
 251         std     $t3,`$FRAME+24`($sp)
 252         ; transfer (ap[0]*bp[0])*n0 to FPU as 4x16-bit values
 253         extrdi  $t4,$t7,16,48
 254         extrdi  $t5,$t7,16,32
 255         extrdi  $t6,$t7,16,16
 256         extrdi  $t7,$t7,16,0
 257         std     $t4,`$FRAME+32`($sp)
 258         std     $t5,`$FRAME+40`($sp)
 259         std     $t6,`$FRAME+48`($sp)
 260         std     $t7,`$FRAME+56`($sp)
 261         lwz     $t0,4($ap)              ; load a[j] as 32-bit word pair
 262         lwz     $t1,0($ap)
 263         lwz     $t2,12($ap)             ; load a[j+1] as 32-bit word pair
 264         lwz     $t3,8($ap)
 265         lwz     $t4,4($np)              ; load n[j] as 32-bit word pair
 266         lwz     $t5,0($np)
 267         lwz     $t6,12($np)             ; load n[j+1] as 32-bit word pair
 268         lwz     $t7,8($np)
 269         lfd     $ba,`$FRAME+0`($sp)
 270         lfd     $bb,`$FRAME+8`($sp)
 271         lfd     $bc,`$FRAME+16`($sp)
 272         lfd     $bd,`$FRAME+24`($sp)
 273         lfd     $na,`$FRAME+32`($sp)
 274         lfd     $nb,`$FRAME+40`($sp)
 275         lfd     $nc,`$FRAME+48`($sp)
 276         lfd     $nd,`$FRAME+56`($sp)
 277         std     $t0,`$FRAME+64`($sp)
 278         std     $t1,`$FRAME+72`($sp)
 279         std     $t2,`$FRAME+80`($sp)
 280         std     $t3,`$FRAME+88`($sp)
 281         std     $t4,`$FRAME+96`($sp)
 282         std     $t5,`$FRAME+104`($sp)
 283         std     $t6,`$FRAME+112`($sp)
 284         std     $t7,`$FRAME+120`($sp)
 285         fcfid   $ba,$ba
 286         fcfid   $bb,$bb
 287         fcfid   $bc,$bc
 288         fcfid   $bd,$bd
 289         fcfid   $na,$na
 290         fcfid   $nb,$nb
 291         fcfid   $nc,$nc
 292         fcfid   $nd,$nd
 293
 294         lfd     $A0,`$FRAME+64`($sp)
 295         lfd     $A1,`$FRAME+72`($sp)
 296         lfd     $A2,`$FRAME+80`($sp)
 297         lfd     $A3,`$FRAME+88`($sp)
 298         lfd     $N0,`$FRAME+96`($sp)
 299         lfd     $N1,`$FRAME+104`($sp)
 300         lfd     $N2,`$FRAME+112`($sp)
 301         lfd     $N3,`$FRAME+120`($sp)
 302         fcfid   $A0,$A0
 303         fcfid   $A1,$A1
 304         fcfid   $A2,$A2
 305         fcfid   $A3,$A3
 306         fcfid   $N0,$N0
 307         fcfid   $N1,$N1
 308         fcfid   $N2,$N2
 309         fcfid   $N3,$N3
 310         addi    $ap,$ap,16
 311         addi    $np,$np,16
 312
 313         fmul    $T1a,$A1,$ba
 314         fmul    $T1b,$A1,$bb
 315         stfd    $A0,8($nap_d)           ; save a[j] in double format
 316         stfd    $A1,16($nap_d)
 317         fmul    $T2a,$A2,$ba
 318         fmul    $T2b,$A2,$bb
 319         stfd    $A2,24($nap_d)          ; save a[j+1] in double format
 320         stfd    $A3,32($nap_d)
 321         fmul    $T3a,$A3,$ba
 322         fmul    $T3b,$A3,$bb
 323         stfd    $N0,40($nap_d)          ; save n[j] in double format
 324         stfd    $N1,48($nap_d)
 325         fmul    $T0a,$A0,$ba
 326         fmul    $T0b,$A0,$bb
 327         stfd    $N2,56($nap_d)          ; save n[j+1] in double format
 328         stfdu   $N3,64($nap_d)
 329
 330         fmadd   $T1a,$A0,$bc,$T1a
 331         fmadd   $T1b,$A0,$bd,$T1b
 332         fmadd   $T2a,$A1,$bc,$T2a
 333         fmadd   $T2b,$A1,$bd,$T2b
 334         fmadd   $T3a,$A2,$bc,$T3a
 335         fmadd   $T3b,$A2,$bd,$T3b
 336         fmul    $dota,$A3,$bc
 337         fmul    $dotb,$A3,$bd
 338
 339         fmadd   $T1a,$N1,$na,$T1a
 340         fmadd   $T1b,$N1,$nb,$T1b
 341         fmadd   $T2a,$N2,$na,$T2a
 342         fmadd   $T2b,$N2,$nb,$T2b
 343         fmadd   $T3a,$N3,$na,$T3a
 344         fmadd   $T3b,$N3,$nb,$T3b
 345         fmadd   $T0a,$N0,$na,$T0a
 346         fmadd   $T0b,$N0,$nb,$T0b
 347
 348         fmadd   $T1a,$N0,$nc,$T1a
 349         fmadd   $T1b,$N0,$nd,$T1b
 350         fmadd   $T2a,$N1,$nc,$T2a
 351         fmadd   $T2b,$N1,$nd,$T2b
 352         fmadd   $T3a,$N2,$nc,$T3a
 353         fmadd   $T3b,$N2,$nd,$T3b
 354         fmadd   $dota,$N3,$nc,$dota
 355         fmadd   $dotb,$N3,$nd,$dotb
 356
 357         fctid   $T0a,$T0a
 358         fctid   $T0b,$T0b
 359         fctid   $T1a,$T1a
 360         fctid   $T1b,$T1b
 361         fctid   $T2a,$T2a
 362         fctid   $T2b,$T2b
 363         fctid   $T3a,$T3a
 364         fctid   $T3b,$T3b
 365
 366         stfd    $T0a,`$FRAME+0`($sp)
 367         stfd    $T0b,`$FRAME+8`($sp)
 368         stfd    $T1a,`$FRAME+16`($sp)
 369         stfd    $T1b,`$FRAME+24`($sp)
 370         stfd    $T2a,`$FRAME+32`($sp)
 371         stfd    $T2b,`$FRAME+40`($sp)
 372         stfd    $T3a,`$FRAME+48`($sp)
 373         stfd    $T3b,`$FRAME+56`($sp)
 374 \f
 375 .align  5
 376 L1st:
 377         lwz     $t0,4($ap)              ; load a[j] as 32-bit word pair
 378         lwz     $t1,0($ap)
 379         lwz     $t2,12($ap)             ; load a[j+1] as 32-bit word pair
 380         lwz     $t3,8($ap)
 381         lwz     $t4,4($np)              ; load n[j] as 32-bit word pair
 382         lwz     $t5,0($np)
 383         lwz     $t6,12($np)             ; load n[j+1] as 32-bit word pair
 384         lwz     $t7,8($np)
 385         std     $t0,`$FRAME+64`($sp)
 386         std     $t1,`$FRAME+72`($sp)
 387         std     $t2,`$FRAME+80`($sp)
 388         std     $t3,`$FRAME+88`($sp)
 389         std     $t4,`$FRAME+96`($sp)
 390         std     $t5,`$FRAME+104`($sp)
 391         std     $t6,`$FRAME+112`($sp)
 392         std     $t7,`$FRAME+120`($sp)
 393         ld      $t0,`$FRAME+0`($sp)
 394         ld      $t1,`$FRAME+8`($sp)
 395         ld      $t2,`$FRAME+16`($sp)
 396         ld      $t3,`$FRAME+24`($sp)
 397         ld      $t4,`$FRAME+32`($sp)
 398         ld      $t5,`$FRAME+40`($sp)
 399         ld      $t6,`$FRAME+48`($sp)
 400         ld      $t7,`$FRAME+56`($sp)
 401         lfd     $A0,`$FRAME+64`($sp)
 402         lfd     $A1,`$FRAME+72`($sp)
 403         lfd     $A2,`$FRAME+80`($sp)
 404         lfd     $A3,`$FRAME+88`($sp)
 405         lfd     $N0,`$FRAME+96`($sp)
 406         lfd     $N1,`$FRAME+104`($sp)
 407         lfd     $N2,`$FRAME+112`($sp)
 408         lfd     $N3,`$FRAME+120`($sp)
 409         fcfid   $A0,$A0
 410         fcfid   $A1,$A1
 411         fcfid   $A2,$A2
 412         fcfid   $A3,$A3
 413         fcfid   $N0,$N0
 414         fcfid   $N1,$N1
 415         fcfid   $N2,$N2
 416         fcfid   $N3,$N3
 417         addi    $ap,$ap,16
 418         addi    $np,$np,16
 419
 420         fmul    $T1a,$A1,$ba
 421         fmul    $T1b,$A1,$bb
 422         fmul    $T2a,$A2,$ba
 423         fmul    $T2b,$A2,$bb
 424         stfd    $A0,8($nap_d)           ; save a[j] in double format
 425         stfd    $A1,16($nap_d)
 426         fmul    $T3a,$A3,$ba
 427         fmul    $T3b,$A3,$bb
 428         fmadd   $T0a,$A0,$ba,$dota
 429         fmadd   $T0b,$A0,$bb,$dotb
 430         stfd    $A2,24($nap_d)          ; save a[j+1] in double format
 431         stfd    $A3,32($nap_d)
 432
 433         fmadd   $T1a,$A0,$bc,$T1a
 434         fmadd   $T1b,$A0,$bd,$T1b
 435         fmadd   $T2a,$A1,$bc,$T2a
 436         fmadd   $T2b,$A1,$bd,$T2b
 437         stfd    $N0,40($nap_d)          ; save n[j] in double format
 438         stfd    $N1,48($nap_d)
 439         fmadd   $T3a,$A2,$bc,$T3a
 440         fmadd   $T3b,$A2,$bd,$T3b
 441          add    $t0,$t0,$carry          ; can not overflow
 442         fmul    $dota,$A3,$bc
 443         fmul    $dotb,$A3,$bd
 444         stfd    $N2,56($nap_d)          ; save n[j+1] in double format
 445         stfdu   $N3,64($nap_d)
 446          srdi   $carry,$t0,16
 447          add    $t1,$t1,$carry
 448          srdi   $carry,$t1,16
 449
 450         fmadd   $T1a,$N1,$na,$T1a
 451         fmadd   $T1b,$N1,$nb,$T1b
 452          insrdi $t0,$t1,16,32
 453         fmadd   $T2a,$N2,$na,$T2a
 454         fmadd   $T2b,$N2,$nb,$T2b
 455          add    $t2,$t2,$carry
 456         fmadd   $T3a,$N3,$na,$T3a
 457         fmadd   $T3b,$N3,$nb,$T3b
 458          srdi   $carry,$t2,16
 459         fmadd   $T0a,$N0,$na,$T0a
 460         fmadd   $T0b,$N0,$nb,$T0b
 461          insrdi $t0,$t2,16,16
 462          add    $t3,$t3,$carry
 463          srdi   $carry,$t3,16
 464
 465         fmadd   $T1a,$N0,$nc,$T1a
 466         fmadd   $T1b,$N0,$nd,$T1b
 467          insrdi $t0,$t3,16,0            ; 0..63 bits
 468         fmadd   $T2a,$N1,$nc,$T2a
 469         fmadd   $T2b,$N1,$nd,$T2b
 470          add    $t4,$t4,$carry
 471         fmadd   $T3a,$N2,$nc,$T3a
 472         fmadd   $T3b,$N2,$nd,$T3b
 473          srdi   $carry,$t4,16
 474         fmadd   $dota,$N3,$nc,$dota
 475         fmadd   $dotb,$N3,$nd,$dotb
 476          add    $t5,$t5,$carry
 477          srdi   $carry,$t5,16
 478          insrdi $t4,$t5,16,32
 479
 480         fctid   $T0a,$T0a
 481         fctid   $T0b,$T0b
 482          add    $t6,$t6,$carry
 483         fctid   $T1a,$T1a
 484         fctid   $T1b,$T1b
 485          srdi   $carry,$t6,16
 486         fctid   $T2a,$T2a
 487         fctid   $T2b,$T2b
 488          insrdi $t4,$t6,16,16
 489         fctid   $T3a,$T3a
 490         fctid   $T3b,$T3b
 491          add    $t7,$t7,$carry
 492          insrdi $t4,$t7,16,0            ; 64..127 bits
 493          srdi   $carry,$t7,16           ; upper 33 bits
 494
 495         stfd    $T0a,`$FRAME+0`($sp)
 496         stfd    $T0b,`$FRAME+8`($sp)
 497         stfd    $T1a,`$FRAME+16`($sp)
 498         stfd    $T1b,`$FRAME+24`($sp)
 499         stfd    $T2a,`$FRAME+32`($sp)
 500         stfd    $T2b,`$FRAME+40`($sp)
 501         stfd    $T3a,`$FRAME+48`($sp)
 502         stfd    $T3b,`$FRAME+56`($sp)
 503          std    $t0,8($tp)              ; tp[j-1]
 504          stdu   $t4,16($tp)             ; tp[j]
 505         bdnz-   L1st
 506 \f
 507         fctid   $dota,$dota
 508         fctid   $dotb,$dotb
 509
 510         ld      $t0,`$FRAME+0`($sp)
 511         ld      $t1,`$FRAME+8`($sp)
 512         ld      $t2,`$FRAME+16`($sp)
 513         ld      $t3,`$FRAME+24`($sp)
 514         ld      $t4,`$FRAME+32`($sp)
 515         ld      $t5,`$FRAME+40`($sp)
 516         ld      $t6,`$FRAME+48`($sp)
 517         ld      $t7,`$FRAME+56`($sp)
 518         stfd    $dota,`$FRAME+64`($sp)
 519         stfd    $dotb,`$FRAME+72`($sp)
 520
 521         add     $t0,$t0,$carry          ; can not overflow
 522         srdi    $carry,$t0,16
 523         add     $t1,$t1,$carry
 524         srdi    $carry,$t1,16
 525         insrdi  $t0,$t1,16,32
 526         add     $t2,$t2,$carry
 527         srdi    $carry,$t2,16
 528         insrdi  $t0,$t2,16,16
 529         add     $t3,$t3,$carry
 530         srdi    $carry,$t3,16
 531         insrdi  $t0,$t3,16,0            ; 0..63 bits
 532         add     $t4,$t4,$carry
 533         srdi    $carry,$t4,16
 534         add     $t5,$t5,$carry
 535         srdi    $carry,$t5,16
 536         insrdi  $t4,$t5,16,32
 537         add     $t6,$t6,$carry
 538         srdi    $carry,$t6,16
 539         insrdi  $t4,$t6,16,16
 540         add     $t7,$t7,$carry
 541         insrdi  $t4,$t7,16,0            ; 64..127 bits
 542         srdi    $carry,$t7,16           ; upper 33 bits
 543         ld      $t6,`$FRAME+64`($sp)
 544         ld      $t7,`$FRAME+72`($sp)
 545
 546         std     $t0,8($tp)              ; tp[j-1]
 547         stdu    $t4,16($tp)             ; tp[j]
 548
 549         add     $t6,$t6,$carry          ; can not overflow
 550         srdi    $carry,$t6,16
 551         add     $t7,$t7,$carry
 552         insrdi  $t6,$t7,48,0
 553         srdi    $ovf,$t7,48
 554         std     $t6,8($tp)              ; tp[num-1]
 555
 556         slwi    $t7,$num,2
 557         subf    $nap_d,$t7,$nap_d       ; rewind pointer
 558 \f
 559         li      $i,8                    ; i=1
 560 .align  5
 561 Louter:
 562         ldx     $t3,$bp,$i      ; bp[i]
 563         ld      $t6,`$FRAME+$TRANSFER+8`($sp)   ; tp[0]
 564         mulld   $t7,$a0,$t3     ; ap[0]*bp[i]
 565
 566         addi    $tp,$sp,`$FRAME+$TRANSFER`
 567         add     $t7,$t7,$t6     ; ap[0]*bp[i]+tp[0]
 568         li      $carry,0
 569         mulld   $t7,$t7,$n0     ; tp[0]*n0
 570         mtctr   $j
 571
 572         ; transfer bp[i] to FPU as 4x16-bit values
 573         extrdi  $t0,$t3,16,48
 574         extrdi  $t1,$t3,16,32
 575         extrdi  $t2,$t3,16,16
 576         extrdi  $t3,$t3,16,0
 577         std     $t0,`$FRAME+0`($sp)
 578         std     $t1,`$FRAME+8`($sp)
 579         std     $t2,`$FRAME+16`($sp)
 580         std     $t3,`$FRAME+24`($sp)
 581         ; transfer (ap[0]*bp[i]+tp[0])*n0 to FPU as 4x16-bit values
 582         extrdi  $t4,$t7,16,48
 583         extrdi  $t5,$t7,16,32
 584         extrdi  $t6,$t7,16,16
 585         extrdi  $t7,$t7,16,0
 586         std     $t4,`$FRAME+32`($sp)
 587         std     $t5,`$FRAME+40`($sp)
 588         std     $t6,`$FRAME+48`($sp)
 589         std     $t7,`$FRAME+56`($sp)
 590
 591         lfd     $A0,8($nap_d)           ; load a[j] in double format
 592         lfd     $A1,16($nap_d)
 593         lfd     $A2,24($nap_d)          ; load a[j+1] in double format
 594         lfd     $A3,32($nap_d)
 595         lfd     $N0,40($nap_d)          ; load n[j] in double format
 596         lfd     $N1,48($nap_d)
 597         lfd     $N2,56($nap_d)          ; load n[j+1] in double format
 598         lfdu    $N3,64($nap_d)
 599
 600         lfd     $ba,`$FRAME+0`($sp)
 601         lfd     $bb,`$FRAME+8`($sp)
 602         lfd     $bc,`$FRAME+16`($sp)
 603         lfd     $bd,`$FRAME+24`($sp)
 604         lfd     $na,`$FRAME+32`($sp)
 605         lfd     $nb,`$FRAME+40`($sp)
 606         lfd     $nc,`$FRAME+48`($sp)
 607         lfd     $nd,`$FRAME+56`($sp)
 608
 609         fcfid   $ba,$ba
 610         fcfid   $bb,$bb
 611         fcfid   $bc,$bc
 612         fcfid   $bd,$bd
 613         fcfid   $na,$na
 614         fcfid   $nb,$nb
 615         fcfid   $nc,$nc
 616         fcfid   $nd,$nd
 617
 618         fmul    $T1a,$A1,$ba
 619         fmul    $T1b,$A1,$bb
 620         fmul    $T2a,$A2,$ba
 621         fmul    $T2b,$A2,$bb
 622         fmul    $T3a,$A3,$ba
 623         fmul    $T3b,$A3,$bb
 624         fmul    $T0a,$A0,$ba
 625         fmul    $T0b,$A0,$bb
 626
 627         fmadd   $T1a,$A0,$bc,$T1a
 628         fmadd   $T1b,$A0,$bd,$T1b
 629         fmadd   $T2a,$A1,$bc,$T2a
 630         fmadd   $T2b,$A1,$bd,$T2b
 631         fmadd   $T3a,$A2,$bc,$T3a
 632         fmadd   $T3b,$A2,$bd,$T3b
 633         fmul    $dota,$A3,$bc
 634         fmul    $dotb,$A3,$bd
 635
 636         fmadd   $T1a,$N1,$na,$T1a
 637         fmadd   $T1b,$N1,$nb,$T1b
 638          lfd    $A0,8($nap_d)           ; load a[j] in double format
 639          lfd    $A1,16($nap_d)
 640         fmadd   $T2a,$N2,$na,$T2a
 641         fmadd   $T2b,$N2,$nb,$T2b
 642          lfd    $A2,24($nap_d)          ; load a[j+1] in double format
 643          lfd    $A3,32($nap_d)
 644         fmadd   $T3a,$N3,$na,$T3a
 645         fmadd   $T3b,$N3,$nb,$T3b
 646         fmadd   $T0a,$N0,$na,$T0a
 647         fmadd   $T0b,$N0,$nb,$T0b
 648
 649         fmadd   $T1a,$N0,$nc,$T1a
 650         fmadd   $T1b,$N0,$nd,$T1b
 651         fmadd   $T2a,$N1,$nc,$T2a
 652         fmadd   $T2b,$N1,$nd,$T2b
 653         fmadd   $T3a,$N2,$nc,$T3a
 654         fmadd   $T3b,$N2,$nd,$T3b
 655         fmadd   $dota,$N3,$nc,$dota
 656         fmadd   $dotb,$N3,$nd,$dotb
 657
 658         fctid   $T0a,$T0a
 659         fctid   $T0b,$T0b
 660         fctid   $T1a,$T1a
 661         fctid   $T1b,$T1b
 662         fctid   $T2a,$T2a
 663         fctid   $T2b,$T2b
 664         fctid   $T3a,$T3a
 665         fctid   $T3b,$T3b
 666
 667         stfd    $T0a,`$FRAME+0`($sp)
 668         stfd    $T0b,`$FRAME+8`($sp)
 669         stfd    $T1a,`$FRAME+16`($sp)
 670         stfd    $T1b,`$FRAME+24`($sp)
 671         stfd    $T2a,`$FRAME+32`($sp)
 672         stfd    $T2b,`$FRAME+40`($sp)
 673         stfd    $T3a,`$FRAME+48`($sp)
 674         stfd    $T3b,`$FRAME+56`($sp)
 675 \f
 676 .align  5
 677 Linner:
 678         fmul    $T1a,$A1,$ba
 679         fmul    $T1b,$A1,$bb
 680         fmul    $T2a,$A2,$ba
 681         fmul    $T2b,$A2,$bb
 682         lfd     $N0,40($nap_d)          ; load n[j] in double format
 683         lfd     $N1,48($nap_d)
 684         fmul    $T3a,$A3,$ba
 685         fmul    $T3b,$A3,$bb
 686         fmadd   $T0a,$A0,$ba,$dota
 687         fmadd   $T0b,$A0,$bb,$dotb
 688         lfd     $N2,56($nap_d)          ; load n[j+1] in double format
 689         lfdu    $N3,64($nap_d)
 690
 691         fmadd   $T1a,$A0,$bc,$T1a
 692         fmadd   $T1b,$A0,$bd,$T1b
 693         fmadd   $T2a,$A1,$bc,$T2a
 694         fmadd   $T2b,$A1,$bd,$T2b
 695          lfd    $A0,8($nap_d)           ; load a[j] in double format
 696          lfd    $A1,16($nap_d)
 697         fmadd   $T3a,$A2,$bc,$T3a
 698         fmadd   $T3b,$A2,$bd,$T3b
 699         fmul    $dota,$A3,$bc
 700         fmul    $dotb,$A3,$bd
 701          lfd    $A2,24($nap_d)          ; load a[j+1] in double format
 702          lfd    $A3,32($nap_d)
 703
 704         fmadd   $T1a,$N1,$na,$T1a
 705         fmadd   $T1b,$N1,$nb,$T1b
 706          ld     $t0,`$FRAME+0`($sp)
 707          ld     $t1,`$FRAME+8`($sp)
 708         fmadd   $T2a,$N2,$na,$T2a
 709         fmadd   $T2b,$N2,$nb,$T2b
 710          ld     $t2,`$FRAME+16`($sp)
 711          ld     $t3,`$FRAME+24`($sp)
 712         fmadd   $T3a,$N3,$na,$T3a
 713         fmadd   $T3b,$N3,$nb,$T3b
 714          add    $t0,$t0,$carry          ; can not overflow
 715          ld     $t4,`$FRAME+32`($sp)
 716          ld     $t5,`$FRAME+40`($sp)
 717         fmadd   $T0a,$N0,$na,$T0a
 718         fmadd   $T0b,$N0,$nb,$T0b
 719          srdi   $carry,$t0,16
 720          add    $t1,$t1,$carry
 721          srdi   $carry,$t1,16
 722          ld     $t6,`$FRAME+48`($sp)
 723          ld     $t7,`$FRAME+56`($sp)
 724
 725         fmadd   $T1a,$N0,$nc,$T1a
 726         fmadd   $T1b,$N0,$nd,$T1b
 727          insrdi $t0,$t1,16,32
 728          ld     $t1,8($tp)              ; tp[j]
 729         fmadd   $T2a,$N1,$nc,$T2a
 730         fmadd   $T2b,$N1,$nd,$T2b
 731          add    $t2,$t2,$carry
 732         fmadd   $T3a,$N2,$nc,$T3a
 733         fmadd   $T3b,$N2,$nd,$T3b
 734          srdi   $carry,$t2,16
 735          insrdi $t0,$t2,16,16
 736         fmadd   $dota,$N3,$nc,$dota
 737         fmadd   $dotb,$N3,$nd,$dotb
 738          add    $t3,$t3,$carry
 739          ldu    $t2,16($tp)             ; tp[j+1]
 740          srdi   $carry,$t3,16
 741          insrdi $t0,$t3,16,0            ; 0..63 bits
 742          add    $t4,$t4,$carry
 743
 744         fctid   $T0a,$T0a
 745         fctid   $T0b,$T0b
 746          srdi   $carry,$t4,16
 747         fctid   $T1a,$T1a
 748         fctid   $T1b,$T1b
 749          add    $t5,$t5,$carry
 750         fctid   $T2a,$T2a
 751         fctid   $T2b,$T2b
 752          srdi   $carry,$t5,16
 753          insrdi $t4,$t5,16,32
 754         fctid   $T3a,$T3a
 755         fctid   $T3b,$T3b
 756          add    $t6,$t6,$carry
 757          srdi   $carry,$t6,16
 758          insrdi $t4,$t6,16,16
 759
 760         stfd    $T0a,`$FRAME+0`($sp)
 761         stfd    $T0b,`$FRAME+8`($sp)
 762          add    $t7,$t7,$carry
 763          addc   $t3,$t0,$t1
 764         stfd    $T1a,`$FRAME+16`($sp)
 765         stfd    $T1b,`$FRAME+24`($sp)
 766          insrdi $t4,$t7,16,0            ; 64..127 bits
 767          srdi   $carry,$t7,16           ; upper 33 bits
 768         stfd    $T2a,`$FRAME+32`($sp)
 769         stfd    $T2b,`$FRAME+40`($sp)
 770          adde   $t5,$t4,$t2
 771         stfd    $T3a,`$FRAME+48`($sp)
 772         stfd    $T3b,`$FRAME+56`($sp)
 773          addze  $carry,$carry
 774          std    $t3,-16($tp)            ; tp[j-1]
 775          std    $t5,-8($tp)             ; tp[j]
 776         bdnz-   Linner
 777 \f
 778         fctid   $dota,$dota
 779         fctid   $dotb,$dotb
 780         ld      $t0,`$FRAME+0`($sp)
 781         ld      $t1,`$FRAME+8`($sp)
 782         ld      $t2,`$FRAME+16`($sp)
 783         ld      $t3,`$FRAME+24`($sp)
 784         ld      $t4,`$FRAME+32`($sp)
 785         ld      $t5,`$FRAME+40`($sp)
 786         ld      $t6,`$FRAME+48`($sp)
 787         ld      $t7,`$FRAME+56`($sp)
 788         stfd    $dota,`$FRAME+64`($sp)
 789         stfd    $dotb,`$FRAME+72`($sp)
 790
 791         add     $t0,$t0,$carry          ; can not overflow
 792         srdi    $carry,$t0,16
 793         add     $t1,$t1,$carry
 794         srdi    $carry,$t1,16
 795         insrdi  $t0,$t1,16,32
 796         add     $t2,$t2,$carry
 797         ld      $t1,8($tp)              ; tp[j]
 798         srdi    $carry,$t2,16
 799         insrdi  $t0,$t2,16,16
 800         add     $t3,$t3,$carry
 801         ldu     $t2,16($tp)             ; tp[j+1]
 802         srdi    $carry,$t3,16
 803         insrdi  $t0,$t3,16,0            ; 0..63 bits
 804         add     $t4,$t4,$carry
 805         srdi    $carry,$t4,16
 806         add     $t5,$t5,$carry
 807         srdi    $carry,$t5,16
 808         insrdi  $t4,$t5,16,32
 809         add     $t6,$t6,$carry
 810         srdi    $carry,$t6,16
 811         insrdi  $t4,$t6,16,16
 812         add     $t7,$t7,$carry
 813         insrdi  $t4,$t7,16,0            ; 64..127 bits
 814         srdi    $carry,$t7,16           ; upper 33 bits
 815         ld      $t6,`$FRAME+64`($sp)
 816         ld      $t7,`$FRAME+72`($sp)
 817
 818         addc    $t3,$t0,$t1
 819         adde    $t5,$t4,$t2
 820         addze   $carry,$carry
 821
 822         std     $t3,-16($tp)            ; tp[j-1]
 823         std     $t5,-8($tp)             ; tp[j]
 824
 825         add     $carry,$carry,$ovf      ; comsume upmost overflow
 826         add     $t6,$t6,$carry          ; can not overflow
 827         srdi    $carry,$t6,16
 828         add     $t7,$t7,$carry
 829         insrdi  $t6,$t7,48,0
 830         srdi    $ovf,$t7,48
 831         std     $t6,0($tp)              ; tp[num-1]
 832
 833         slwi    $t7,$num,2
 834         addi    $i,$i,8
 835         subf    $nap_d,$t7,$nap_d       ; rewind pointer
 836         cmpw    $i,$num
 837         blt-    Louter
 838 \f
 839         subf    $np,$num,$np    ; rewind np
 840         addi    $j,$j,1         ; restore counter
 841         subfc   $i,$i,$i        ; j=0 and "clear" XER[CA]
 842         addi    $tp,$sp,`$FRAME+$TRANSFER+8`
 843         addi    $t4,$sp,`$FRAME+$TRANSFER+16`
 844         addi    $t5,$np,8
 845         addi    $t6,$rp,8
 846         mtctr   $j
 847
 848 .align  4
 849 Lsub:   ldx     $t0,$tp,$i
 850         ldx     $t1,$np,$i
 851         ldx     $t2,$t4,$i
 852         ldx     $t3,$t5,$i
 853         subfe   $t0,$t1,$t0     ; tp[j]-np[j]
 854         subfe   $t2,$t3,$t2     ; tp[j+1]-np[j+1]
 855         stdx    $t0,$rp,$i
 856         stdx    $t2,$t6,$i
 857         addi    $i,$i,16
 858         bdnz-   Lsub
 859
 860         li      $i,0
 861         subfe   $ovf,$i,$ovf    ; handle upmost overflow bit
 862         and     $ap,$tp,$ovf
 863         andc    $np,$rp,$ovf
 864         or      $ap,$ap,$np     ; ap=borrow?tp:rp
 865         addi    $t7,$ap,8
 866         mtctr   $j
 867
 868 .align  4
 869 Lcopy:                          ; copy or in-place refresh
 870         ldx     $t0,$ap,$i
 871         ldx     $t1,$t7,$i
 872         std     $i,8($nap_d)    ; zap nap_d
 873         std     $i,16($nap_d)
 874         std     $i,24($nap_d)
 875         std     $i,32($nap_d)
 876         std     $i,40($nap_d)
 877         std     $i,48($nap_d)
 878         std     $i,56($nap_d)
 879         stdu    $i,64($nap_d)
 880         stdx    $t0,$rp,$i
 881         stdx    $t1,$t6,$i
 882         stdx    $i,$tp,$i       ; zap tp at once
 883         stdx    $i,$t4,$i
 884         addi    $i,$i,16
 885         bdnz-   Lcopy
 886 \f
 887         $POP    r14,`2*$SIZE_T`($sp)
 888         $POP    r15,`3*$SIZE_T`($sp)
 889         $POP    r16,`4*$SIZE_T`($sp)
 890         $POP    r17,`5*$SIZE_T`($sp)
 891         $POP    r18,`6*$SIZE_T`($sp)
 892         $POP    r19,`7*$SIZE_T`($sp)
 893         $POP    r20,`8*$SIZE_T`($sp)
 894         $POP    r21,`9*$SIZE_T`($sp)
 895         $POP    r22,`10*$SIZE_T`($sp)
 896         $POP    r23,`11*$SIZE_T`($sp)
 897         lfd     f14,`12*$SIZE_T+0`($sp)
 898         lfd     f15,`12*$SIZE_T+8`($sp)
 899         lfd     f16,`12*$SIZE_T+16`($sp)
 900         lfd     f17,`12*$SIZE_T+24`($sp)
 901         lfd     f18,`12*$SIZE_T+32`($sp)
 902         lfd     f19,`12*$SIZE_T+40`($sp)
 903         lfd     f20,`12*$SIZE_T+48`($sp)
 904         lfd     f21,`12*$SIZE_T+56`($sp)
 905         lfd     f22,`12*$SIZE_T+64`($sp)
 906         lfd     f23,`12*$SIZE_T+72`($sp)
 907         lfd     f24,`12*$SIZE_T+80`($sp)
 908         lfd     f25,`12*$SIZE_T+88`($sp)
 909         $POP    $sp,0($sp)
 910         li      r3,1    ; signal "handled"
 911         blr
 912         .long   0
 913 .asciz  "Montgomery Multiplication for PPC64, CRYPTOGAMS by <appro\@fy.chalmers.se>"
 914 ___
 915
 916 $code =~ s/\`([^\`]*)\`/eval $1/gem;
 917 print $code;
 918 close STDOUT;