sysdeps/x86_64/fpu/multiarch/svml_d_exp8_core_avx512.S

   1 /* Function exp vectorized with AVX-512. KNL and SKX versions.
   2    Copyright (C) 2014-2023 Free Software Foundation, Inc.
   3    This file is part of the GNU C Library.
   4
   5    The GNU C Library is free software; you can redistribute it and/or
   6    modify it under the terms of the GNU Lesser General Public
   7    License as published by the Free Software Foundation; either
   8    version 2.1 of the License, or (at your option) any later version.
   9
  10    The GNU C Library is distributed in the hope that it will be useful,
  11    but WITHOUT ANY WARRANTY; without even the implied warranty of
  12    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13    Lesser General Public License for more details.
  14
  15    You should have received a copy of the GNU Lesser General Public
  16    License along with the GNU C Library; if not, see
  17    <https://www.gnu.org/licenses/>.  */
  18
  19 #include <sysdep.h>
  20 #include "svml_d_exp_data.h"
  21 #include "svml_d_wrapper_impl.h"
  22
  23         .section .text.evex512, "ax", @progbits
  24 ENTRY (_ZGVeN8v_exp_knl)
  25 /*
  26    ALGORITHM DESCRIPTION:
  27
  28      Argument representation:
  29      N = rint(X*2^k/ln2) = 2^k*M+j
  30      X = N*ln2/2^k + r = M*ln2 + ln2*(j/2^k) + r
  31      then -ln2/2^(k+1) < r < ln2/2^(k+1)
  32      Alternatively:
  33      N = trunc(X*2^k/ln2)
  34      then 0 < r < ln2/2^k
  35
  36      Result calculation:
  37      exp(X) = exp(M*ln2 + ln2*(j/2^k) + r)
  38      = 2^M * 2^(j/2^k) * exp(r)
  39      2^M is calculated by bit manipulation
  40      2^(j/2^k) is stored in table
  41      exp(r) is approximated by polynomial
  42
  43      The table lookup is skipped if k = 0.  */
  44
  45         pushq     %rbp
  46         cfi_adjust_cfa_offset (8)
  47         cfi_rel_offset (%rbp, 0)
  48         movq      %rsp, %rbp
  49         cfi_def_cfa_register (%rbp)
  50         andq      $-64, %rsp
  51         subq      $1280, %rsp
  52         movq      __svml_dexp_data@GOTPCREL(%rip), %rax
  53
  54 /* dR = X - dN*dbLn2hi, dbLn2hi is 52-8-k hi bits of ln2/2^k */
  55         vmovaps   %zmm0, %zmm8
  56
  57 /* iAbsX = (int)(lX>>32), lX = *(longlong*)&X */
  58         vpsrlq    $32, %zmm0, %zmm1
  59
  60 /* iAbsX = iAbsX&iAbsMask */
  61         movl      $255, %edx
  62         vpmovqd   %zmm1, %ymm2
  63         kmovw     %edx, %k2
  64
  65 /* iRangeMask = (iAbsX>iDomainRange) */
  66         movl      $-1, %ecx
  67
  68 /* table lookup for dT[j] = 2^(j/2^k) */
  69         vpxord    %zmm11, %zmm11, %zmm11
  70         vmovups __dbInvLn2(%rax), %zmm5
  71         vmovups __dbLn2hi(%rax), %zmm7
  72         kxnorw    %k3, %k3, %k3
  73
  74 /* dM = X*dbInvLn2+dbShifter, dbInvLn2 = 2^k/Ln2 */
  75         vfmadd213pd __dbShifter(%rax), %zmm0, %zmm5
  76         vmovups __dPC2(%rax), %zmm12
  77
  78 /* dN = dM-dbShifter, dN = rint(X*2^k/Ln2) */
  79         vsubpd __dbShifter(%rax), %zmm5, %zmm9
  80         vmovups __lIndexMask(%rax), %zmm4
  81         vfnmadd231pd %zmm9, %zmm7, %zmm8
  82         vpandd __iAbsMask(%rax), %zmm2, %zmm2{%k2}
  83
  84 /* lIndex = (*(longlong*)&dM)&lIndexMask, lIndex is the lower K bits of lM */
  85         vpandq    %zmm4, %zmm5, %zmm10
  86         vgatherqpd (%rax,%zmm10,8), %zmm11{%k3}
  87         vpcmpgtd __iDomainRange(%rax), %zmm2, %k1{%k2}
  88
  89 /* lM = (*(longlong*)&dM)&(~lIndexMask) */
  90         vpandnq   %zmm5, %zmm4, %zmm6
  91         vpbroadcastd %ecx, %zmm3{%k1}{z}
  92
  93 /* lM = lM<<(52-K), 2^M */
  94         vpsllq    $42, %zmm6, %zmm14
  95
  96 /* dR = dR - dN*dbLn2lo, dbLn2lo is 40..94 bits of lo part of ln2/2^k */
  97         vfnmadd132pd __dbLn2lo(%rax), %zmm8, %zmm9
  98
  99 /* Mask = iRangeMask?1:0, set mask for overflow/underflow */
 100         vptestmd  %zmm3, %zmm3, %k0{%k2}
 101
 102 /* exp(r) = b0+r*(b0+r*(b1+r*b2)) */
 103         vfmadd213pd __dPC1(%rax), %zmm9, %zmm12
 104         kmovw     %k0, %ecx
 105         movzbl    %cl, %ecx
 106         vfmadd213pd __dPC0(%rax), %zmm9, %zmm12
 107         vfmadd213pd __dPC0(%rax), %zmm9, %zmm12
 108
 109 /* 2^(j/2^k) * exp(r) */
 110         vmulpd    %zmm12, %zmm11, %zmm13
 111
 112 /* multiply by 2^M through integer add */
 113         vpaddq    %zmm14, %zmm13, %zmm1
 114         testl     %ecx, %ecx
 115         jne       .LBL_1_3
 116
 117 .LBL_1_2:
 118         cfi_remember_state
 119         vmovaps   %zmm1, %zmm0
 120         movq      %rbp, %rsp
 121         cfi_def_cfa_register (%rsp)
 122         popq      %rbp
 123         cfi_adjust_cfa_offset (-8)
 124         cfi_restore (%rbp)
 125         ret
 126
 127 .LBL_1_3:
 128         cfi_restore_state
 129         vmovups   %zmm0, 1152(%rsp)
 130         vmovups   %zmm1, 1216(%rsp)
 131         je        .LBL_1_2
 132
 133         xorb      %dl, %dl
 134         kmovw     %k4, 1048(%rsp)
 135         xorl      %eax, %eax
 136         kmovw     %k5, 1040(%rsp)
 137         kmovw     %k6, 1032(%rsp)
 138         kmovw     %k7, 1024(%rsp)
 139         vmovups   %zmm16, 960(%rsp)
 140         vmovups   %zmm17, 896(%rsp)
 141         vmovups   %zmm18, 832(%rsp)
 142         vmovups   %zmm19, 768(%rsp)
 143         vmovups   %zmm20, 704(%rsp)
 144         vmovups   %zmm21, 640(%rsp)
 145         vmovups   %zmm22, 576(%rsp)
 146         vmovups   %zmm23, 512(%rsp)
 147         vmovups   %zmm24, 448(%rsp)
 148         vmovups   %zmm25, 384(%rsp)
 149         vmovups   %zmm26, 320(%rsp)
 150         vmovups   %zmm27, 256(%rsp)
 151         vmovups   %zmm28, 192(%rsp)
 152         vmovups   %zmm29, 128(%rsp)
 153         vmovups   %zmm30, 64(%rsp)
 154         vmovups   %zmm31, (%rsp)
 155         movq      %rsi, 1064(%rsp)
 156         movq      %rdi, 1056(%rsp)
 157         movq      %r12, 1096(%rsp)
 158         cfi_offset_rel_rsp (12, 1096)
 159         movb      %dl, %r12b
 160         movq      %r13, 1088(%rsp)
 161         cfi_offset_rel_rsp (13, 1088)
 162         movl      %ecx, %r13d
 163         movq      %r14, 1080(%rsp)
 164         cfi_offset_rel_rsp (14, 1080)
 165         movl      %eax, %r14d
 166         movq      %r15, 1072(%rsp)
 167         cfi_offset_rel_rsp (15, 1072)
 168         cfi_remember_state
 169
 170 .LBL_1_6:
 171         btl       %r14d, %r13d
 172         jc        .LBL_1_12
 173
 174 .LBL_1_7:
 175         lea       1(%r14), %esi
 176         btl       %esi, %r13d
 177         jc        .LBL_1_10
 178
 179 .LBL_1_8:
 180         addb      $1, %r12b
 181         addl      $2, %r14d
 182         cmpb      $16, %r12b
 183         jb        .LBL_1_6
 184
 185         kmovw     1048(%rsp), %k4
 186         movq      1064(%rsp), %rsi
 187         kmovw     1040(%rsp), %k5
 188         movq      1056(%rsp), %rdi
 189         kmovw     1032(%rsp), %k6
 190         movq      1096(%rsp), %r12
 191         cfi_restore (%r12)
 192         movq      1088(%rsp), %r13
 193         cfi_restore (%r13)
 194         kmovw     1024(%rsp), %k7
 195         vmovups   960(%rsp), %zmm16
 196         vmovups   896(%rsp), %zmm17
 197         vmovups   832(%rsp), %zmm18
 198         vmovups   768(%rsp), %zmm19
 199         vmovups   704(%rsp), %zmm20
 200         vmovups   640(%rsp), %zmm21
 201         vmovups   576(%rsp), %zmm22
 202         vmovups   512(%rsp), %zmm23
 203         vmovups   448(%rsp), %zmm24
 204         vmovups   384(%rsp), %zmm25
 205         vmovups   320(%rsp), %zmm26
 206         vmovups   256(%rsp), %zmm27
 207         vmovups   192(%rsp), %zmm28
 208         vmovups   128(%rsp), %zmm29
 209         vmovups   64(%rsp), %zmm30
 210         vmovups   (%rsp), %zmm31
 211         movq      1080(%rsp), %r14
 212         cfi_restore (%r14)
 213         movq      1072(%rsp), %r15
 214         cfi_restore (%r15)
 215         vmovups   1216(%rsp), %zmm1
 216         jmp       .LBL_1_2
 217
 218 .LBL_1_10:
 219         cfi_restore_state
 220         movzbl    %r12b, %r15d
 221         shlq      $4, %r15
 222         vmovsd    1160(%rsp,%r15), %xmm0
 223         call      JUMPTARGET(exp)
 224         vmovsd    %xmm0, 1224(%rsp,%r15)
 225         jmp       .LBL_1_8
 226
 227 .LBL_1_12:
 228         movzbl    %r12b, %r15d
 229         shlq      $4, %r15
 230         vmovsd    1152(%rsp,%r15), %xmm0
 231         call      JUMPTARGET(exp)
 232         vmovsd    %xmm0, 1216(%rsp,%r15)
 233         jmp       .LBL_1_7
 234 END (_ZGVeN8v_exp_knl)
 235
 236 ENTRY (_ZGVeN8v_exp_skx)
 237 /*
 238    ALGORITHM DESCRIPTION:
 239
 240      Argument representation:
 241      N = rint(X*2^k/ln2) = 2^k*M+j
 242      X = N*ln2/2^k + r = M*ln2 + ln2*(j/2^k) + r
 243      then -ln2/2^(k+1) < r < ln2/2^(k+1)
 244      Alternatively:
 245      N = trunc(X*2^k/ln2)
 246      then 0 < r < ln2/2^k
 247
 248      Result calculation:
 249      exp(X) = exp(M*ln2 + ln2*(j/2^k) + r)
 250      = 2^M * 2^(j/2^k) * exp(r)
 251      2^M is calculated by bit manipulation
 252      2^(j/2^k) is stored in table
 253      exp(r) is approximated by polynomial
 254
 255      The table lookup is skipped if k = 0.  */
 256
 257         pushq     %rbp
 258         cfi_adjust_cfa_offset (8)
 259         cfi_rel_offset (%rbp, 0)
 260         movq      %rsp, %rbp
 261         cfi_def_cfa_register (%rbp)
 262         andq      $-64, %rsp
 263         subq      $1280, %rsp
 264         movq      __svml_dexp_data@GOTPCREL(%rip), %rax
 265
 266 /* table lookup for dT[j] = 2^(j/2^k) */
 267         kxnorw    %k1, %k1, %k1
 268
 269 /* iAbsX = (int)(lX>>32), lX = *(longlong*)&X */
 270         vpsrlq    $32, %zmm0, %zmm1
 271         vmovups __dbInvLn2(%rax), %zmm7
 272         vmovups __dbShifter(%rax), %zmm5
 273         vmovups __lIndexMask(%rax), %zmm6
 274         vmovups __dbLn2hi(%rax), %zmm9
 275         vmovups __dPC0(%rax), %zmm12
 276
 277 /* dM = X*dbInvLn2+dbShifter, dbInvLn2 = 2^k/Ln2 */
 278         vfmadd213pd %zmm5, %zmm0, %zmm7
 279         vpmovqd   %zmm1, %ymm2
 280
 281 /* dN = dM-dbShifter, dN = rint(X*2^k/Ln2) */
 282         vsubpd    %zmm5, %zmm7, %zmm11
 283
 284 /* iAbsX = iAbsX&iAbsMask */
 285         vpand __iAbsMask(%rax), %ymm2, %ymm3
 286
 287 /* dR = X - dN*dbLn2hi, dbLn2hi is 52-8-k hi bits of ln2/2^k */
 288         vmovaps   %zmm0, %zmm10
 289         vfnmadd231pd %zmm11, %zmm9, %zmm10
 290         vmovups __dPC2(%rax), %zmm9
 291
 292 /* dR = dR - dN*dbLn2lo, dbLn2lo is 40..94 bits of lo part of ln2/2^k */
 293         vfnmadd132pd __dbLn2lo(%rax), %zmm10, %zmm11
 294
 295 /* exp(r) = b0+r*(b0+r*(b1+r*b2)) */
 296         vfmadd213pd __dPC1(%rax), %zmm11, %zmm9
 297         vfmadd213pd %zmm12, %zmm11, %zmm9
 298         vfmadd213pd %zmm12, %zmm11, %zmm9
 299
 300 /* iRangeMask = (iAbsX>iDomainRange) */
 301         vpcmpgtd __iDomainRange(%rax), %ymm3, %ymm4
 302
 303 /* Mask = iRangeMask?1:0, set mask for overflow/underflow */
 304         vmovmskps %ymm4, %ecx
 305
 306 /* lIndex = (*(longlong*)&dM)&lIndexMask, lIndex is the lower K bits of lM */
 307         vpandq    %zmm6, %zmm7, %zmm13
 308         vpmovqd   %zmm13, %ymm14
 309         vpxord    %zmm15, %zmm15, %zmm15
 310         vgatherdpd (%rax,%ymm14,8), %zmm15{%k1}
 311
 312 /* 2^(j/2^k) * exp(r) */
 313         vmulpd    %zmm9, %zmm15, %zmm10
 314
 315 /* lM = (*(longlong*)&dM)&(~lIndexMask) */
 316         vpandnq   %zmm7, %zmm6, %zmm8
 317
 318 /* lM = lM<<(52-K), 2^M */
 319         vpsllq    $42, %zmm8, %zmm1
 320
 321 /* multiply by 2^M through integer add */
 322         vpaddq    %zmm1, %zmm10, %zmm1
 323         testl     %ecx, %ecx
 324         jne       .LBL_2_3
 325
 326 .LBL_2_2:
 327         cfi_remember_state
 328         vmovaps   %zmm1, %zmm0
 329         movq      %rbp, %rsp
 330         cfi_def_cfa_register (%rsp)
 331         popq      %rbp
 332         cfi_adjust_cfa_offset (-8)
 333         cfi_restore (%rbp)
 334         ret
 335
 336 .LBL_2_3:
 337         cfi_restore_state
 338         vmovups   %zmm0, 1152(%rsp)
 339         vmovups   %zmm1, 1216(%rsp)
 340         je        .LBL_2_2
 341
 342         xorb      %dl, %dl
 343         xorl      %eax, %eax
 344         kmovw     %k4, 1048(%rsp)
 345         kmovw     %k5, 1040(%rsp)
 346         kmovw     %k6, 1032(%rsp)
 347         kmovw     %k7, 1024(%rsp)
 348         vmovups   %zmm16, 960(%rsp)
 349         vmovups   %zmm17, 896(%rsp)
 350         vmovups   %zmm18, 832(%rsp)
 351         vmovups   %zmm19, 768(%rsp)
 352         vmovups   %zmm20, 704(%rsp)
 353         vmovups   %zmm21, 640(%rsp)
 354         vmovups   %zmm22, 576(%rsp)
 355         vmovups   %zmm23, 512(%rsp)
 356         vmovups   %zmm24, 448(%rsp)
 357         vmovups   %zmm25, 384(%rsp)
 358         vmovups   %zmm26, 320(%rsp)
 359         vmovups   %zmm27, 256(%rsp)
 360         vmovups   %zmm28, 192(%rsp)
 361         vmovups   %zmm29, 128(%rsp)
 362         vmovups   %zmm30, 64(%rsp)
 363         vmovups   %zmm31, (%rsp)
 364         movq      %rsi, 1064(%rsp)
 365         movq      %rdi, 1056(%rsp)
 366         movq      %r12, 1096(%rsp)
 367         cfi_offset_rel_rsp (12, 1096)
 368         movb      %dl, %r12b
 369         movq      %r13, 1088(%rsp)
 370         cfi_offset_rel_rsp (13, 1088)
 371         movl      %ecx, %r13d
 372         movq      %r14, 1080(%rsp)
 373         cfi_offset_rel_rsp (14, 1080)
 374         movl      %eax, %r14d
 375         movq      %r15, 1072(%rsp)
 376         cfi_offset_rel_rsp (15, 1072)
 377         cfi_remember_state
 378
 379 .LBL_2_6:
 380         btl       %r14d, %r13d
 381         jc        .LBL_2_12
 382
 383 .LBL_2_7:
 384         lea       1(%r14), %esi
 385         btl       %esi, %r13d
 386         jc        .LBL_2_10
 387
 388 .LBL_2_8:
 389         incb      %r12b
 390         addl      $2, %r14d
 391         cmpb      $16, %r12b
 392         jb        .LBL_2_6
 393
 394         kmovw     1048(%rsp), %k4
 395         kmovw     1040(%rsp), %k5
 396         kmovw     1032(%rsp), %k6
 397         kmovw     1024(%rsp), %k7
 398         vmovups   960(%rsp), %zmm16
 399         vmovups   896(%rsp), %zmm17
 400         vmovups   832(%rsp), %zmm18
 401         vmovups   768(%rsp), %zmm19
 402         vmovups   704(%rsp), %zmm20
 403         vmovups   640(%rsp), %zmm21
 404         vmovups   576(%rsp), %zmm22
 405         vmovups   512(%rsp), %zmm23
 406         vmovups   448(%rsp), %zmm24
 407         vmovups   384(%rsp), %zmm25
 408         vmovups   320(%rsp), %zmm26
 409         vmovups   256(%rsp), %zmm27
 410         vmovups   192(%rsp), %zmm28
 411         vmovups   128(%rsp), %zmm29
 412         vmovups   64(%rsp), %zmm30
 413         vmovups   (%rsp), %zmm31
 414         vmovups   1216(%rsp), %zmm1
 415         movq      1064(%rsp), %rsi
 416         movq      1056(%rsp), %rdi
 417         movq      1096(%rsp), %r12
 418         cfi_restore (%r12)
 419         movq      1088(%rsp), %r13
 420         cfi_restore (%r13)
 421         movq      1080(%rsp), %r14
 422         cfi_restore (%r14)
 423         movq      1072(%rsp), %r15
 424         cfi_restore (%r15)
 425         jmp       .LBL_2_2
 426
 427 .LBL_2_10:
 428         cfi_restore_state
 429         movzbl    %r12b, %r15d
 430         shlq      $4, %r15
 431         vmovsd    1160(%rsp,%r15), %xmm0
 432         vzeroupper
 433         vmovsd    1160(%rsp,%r15), %xmm0
 434         call      JUMPTARGET(exp)
 435         vmovsd    %xmm0, 1224(%rsp,%r15)
 436         jmp       .LBL_2_8
 437
 438 .LBL_2_12:
 439         movzbl    %r12b, %r15d
 440         shlq      $4, %r15
 441         vmovsd    1152(%rsp,%r15), %xmm0
 442         vzeroupper
 443         vmovsd    1152(%rsp,%r15), %xmm0
 444         call      JUMPTARGET(exp)
 445         vmovsd    %xmm0, 1216(%rsp,%r15)
 446         jmp       .LBL_2_7
 447
 448 END (_ZGVeN8v_exp_skx)