[Mesa-dev] [RFC 1/3] gk110/ir: Add rcp f64 implementation

Boyan Ding boyan.j.ding at gmail.com
Sun Mar 5 15:34:56 UTC 2017


Signed-off-by: Boyan Ding <boyan.j.ding at gmail.com>
---
 src/gallium/drivers/nouveau/codegen/lib/gk110.asm  | 156 ++++++++++++++++++++-
 .../drivers/nouveau/codegen/lib/gk110.asm.h        |  96 ++++++++++++-
 2 files changed, 248 insertions(+), 4 deletions(-)

diff --git a/src/gallium/drivers/nouveau/codegen/lib/gk110.asm b/src/gallium/drivers/nouveau/codegen/lib/gk110.asm
index b9c05a04b9..871571e1c3 100644
--- a/src/gallium/drivers/nouveau/codegen/lib/gk110.asm
+++ b/src/gallium/drivers/nouveau/codegen/lib/gk110.asm
@@ -83,11 +83,165 @@ gk110_div_s32:
    $p0 sub b32 $r1 $r1 $r2
    $p0 add b32 $r0 $r0 0x1
    $p3 cvt s32 $r0 neg s32 $r0
-   sched 0x04 0x2e 0x04 0x28 0x04 0x20 0x2c
+   sched 0x04 0x2e 0x28 0x04 0x28 0x28 0x28
    $p2 cvt s32 $r1 neg s32 $r1
    ret
 
+// RCP F64
+//
+// INPUT:   $r0d
+// OUTPUT:  $r0d
+// CLOBBER: $r2 - $r9, $p0
+//
+// The core of RCP and RSQ implementation is Newton-Raphson step, which is
+// used to find successively better approximation from an imprecise initial
+// value (single precision rcp in RCP and rsqrt64h in RSQ).
+//
+// The formula of Newton-Raphson step used in RCP(x) is:
+//     RCP_{n + 1} = 2 * RCP_{n} - x * RCP_{n} * RCP_{n}
+// The following code uses 2 FMAs for each step, and it will basically
+// look like:
+//     tmp = -src * RCP_{n} + 1
+//     RCP_{n + 1} = RCP_{n} * tmp + RCP_{n}
+//
 gk110_rcp_f64:
+   // Step1: classify input according to exponent and value, and calculate
+   // result for 0/inf/nan, $r2 holds the exponent value. The 0xb14 in ext
+   // instruction means to extract 0xb bits starting from bit 0x14
+   ext u32 $r2 $r1 0xb14
+   add b32 $r3 $r2 0xffffffff
+   joinat #rcp_L3
+   // (exponent-1) > 0x7fd (unsigned) means exponent is either 0x7ff of 0.
+   // There are three cases: nan, inf, and denorm (including 0)
+   set b32 $p0 0x1 gt u32 $r3 0x7fd
+   // $r3: 0 for norms, 0x36 for denorms, -1 for others
+   mov b32 $r3 0x0
+   sched 0x2b 0x04 0x2d 0x2b 0x04 0x2b 0x28
+   (not $p0) bra #rcp_L2
+   // Nan/Inf/denorm goes here
+   mov b32 $r3 0xffffffff
+   // A number is NaN if its abs value is greater than inf
+   set $p0 0x1 gtu f64 abs $r0d 0x7ff0000000000000
+   (not $p0) bra #rcp_L4
+   // NaN -> NaN
+   or b32 $r1 $r1 0x80000
+   bra #rcp_L2
+rcp_L4:
+   and b32 $r4 $r1 0x7ff00000
+   sched 0x28 0x2b 0x04 0x28 0x2b 0x2d 0x2b
+   // Other values with nonzero in exponent field should be inf
+   set b32 $p0 0x1 eq s32 $r4 0x0
+   $p0 bra #rcp_L5
+   // +/-Inf -> +/-0
+   xor b32 $r1 $r1 0x7ff00000
+   mov b32 $r0 0x0
+   bra #rcp_L2
+rcp_L5:
+   set $p0 0x1 gtu f64 abs $r0d 0x0
+   $p0 bra #rcp_L6
+   // +/-0 -> +/-Inf
+   sched 0x28 0x2b 0x20 0x28 0x2f 0x28 0x2b
+   or b32 $r1 $r1 0x7ff00000
+   bra #rcp_L2
+rcp_L6:
+   // non-0 denorms: multiply with 2^54 (the 0x36 in $r3), join with norms
+   mul rn f64 $r0d $r0d 0x4350000000000000
+   mov b32 $r3 0x36
+rcp_L2:
+   join nop
+rcp_L3:
+   // All numbers with -1 in $r3 have their result ready in $r0d, return them
+   // others need further calculation
+   set b32 $p0 0x1 lt s32 $r3 0x0
+   $p0 bra #rcp_end
+   sched 0x28 0x28 0x04 0x28 0x2b 0x04 0x28
+   // Step 2: Before the real calculation goes on, renormalize the values to
+   // range [1, 2) by setting exponent field to 0x3ff (the exponent of 1)
+   // result in $r6d. The exponent will be recovered later.
+   ext u32 $r2 $r1 0xb14
+   and b32 $r7 $r1 0x800fffff
+   add b32 $r7 $r7 0x3ff00000
+   mov b32 $r6 $r0
+   // Step 3: Convert new value to float (no overflow will occur due to step
+   // 2), calculate rcp and do newton-raphson step once
+   cvt rz f32 $r5 f64 $r6d
+   rcp f32 $r4 $r5
+   mov b32 $r0 0xbf800000
+   sched 0x28 0x28 0x2a 0x2b 0x2e 0x28 0x2e
+   fma rn f32 $r5 $r4 $r5 $r0
+   add ftz rn f32 $r5 neg $r5 neg 0x0
+   fma rn f32 $r0 $r4 $r5 $r4
+   // Step 4: convert result $r0 back to double, do newton-raphson steps
+   cvt f64 $r0d f32 $r0
+   cvt f64 $r6d f64 neg $r6d
+   mov b32 $r9 0x3ff00000
+   mov b32 $r8 0x0
+   sched 0x29 0x29 0x29 0x29 0x29 0x29 0x29
+   // 4 Newton-Raphson Steps, tmp in $r4d, result in $r0d
+   fma rn f64 $r4d $r6d $r0d $r8d
+   fma rn f64 $r0d $r0d $r4d $r0d
+   fma rn f64 $r4d $r6d $r0d $r8d
+   fma rn f64 $r0d $r0d $r4d $r0d
+   fma rn f64 $r4d $r6d $r0d $r8d
+   fma rn f64 $r0d $r0d $r4d $r0d
+   fma rn f64 $r4d $r6d $r0d $r8d
+   sched 0x20 0x28 0x28 0x28 0x28 0x28 0x28
+   fma rn f64 $r0d $r0d $r4d $r0d
+   // Step 5: Exponent recovery and final processing
+   // The exponent is recovered by adding what we added to the exponent.
+   // Suppose we want to calculate rcp(x), but we have rcp(cx), then
+   //     rcp(x) = c * rcp(cx)
+   // The delta in exponent comes from two sources:
+   //   1) The renormalization in step 2. The delta is:
+   //      0x3ff - $r2
+   //   2) (For the denorm input) The 2^54 we multiplied at rcp_L6, stored
+   //      in $r3
+   // These 2 sources are calculated in the first two lines below, and then
+   // added to the exponent extracted from the result above.
+   // Note that after processing, the new exponent may >= 0x7ff (inf)
+   // or <= 0 (denorm). Those cases will be handled respectively below
+   subr b32 $r2 $r2 0x3ff
+   add b32 $r4 $r2 $r3
+   ext u32 $r3 $r1 0xb14
+   // New exponent in $r3
+   add b32 $r3 $r3 $r4
+   add b32 $r2 $r3 0xffffffff
+   // (exponent-1) < 0x7fe (unsigned) means the result is in norm range
+   // (same logic as in step 1)
+   set b32 $p0 0x1 lt u32 $r2 0x7fe
+   sched 0x2b 0x28 0x2b 0x28 0x28 0x2b 0x20
+   (not $p0) bra #rcp_L7
+   // Norms: convert exponents back and return
+   shl b32 $r4 $r4 clamp 0x14
+   add b32 $r1 $r4 $r1
+   bra #rcp_end
+rcp_L7:
+   // New exponent >= 0x7ff means that result is inf
+   set b32 $p0 0x1 ge s32 $r3 0x7ff
+   (not $p0) bra #rcp_L8
+   // Infinity
+   and b32 $r1 $r1 0x80000000
+   sched 0x25 0x28 0x2b 0x23 0x25 0x28 0x23
+   mov b32 $r0 0x0
+   add b32 $r1 $r1 0x7ff00000
+   bra #rcp_end
+rcp_L8:
+   // denorms, they only fall within a small range, can't be smaller than
+   // 0x0004000000000000, which means if we set the exponent field to 1,
+   // we can get the final result by mutiplying it with 1/2 or 1/4. Decide
+   // which one of the two is needed with exponent value, if not 0, 1/4 is
+   // used, 1/2 otherwise
+   set b32 $p0 0x1 ne u32 $r3 0x0
+   and b32 $r1 $r1 0x800fffff
+   $p0 mov b32 $r7 0x3fd00000
+   (not $p0) mov b32 $r7 0x3fe00000
+   sched 0x25 0x28 0x2c 0x2e 0x2e 0x00 0x00
+   add b32 $r1 $r1 0x00100000
+   mov b32 $r6 0x0
+   mul rn f64 $r0d $r0d $r6d
+rcp_end:
+   ret
+
 gk110_rsq_f64:
    ret
 
diff --git a/src/gallium/drivers/nouveau/codegen/lib/gk110.asm.h b/src/gallium/drivers/nouveau/codegen/lib/gk110.asm.h
index 8d00e2a224..ce937a71f9 100644
--- a/src/gallium/drivers/nouveau/codegen/lib/gk110.asm.h
+++ b/src/gallium/drivers/nouveau/codegen/lib/gk110.asm.h
@@ -65,11 +65,101 @@ uint64_t gk110_builtin_code[] = {
 	0xe088000001000406,
 	0x4000000000800001,
 	0xe6010000000ce802,
-	0x08b08010a010b810,
+	0x08a0a0a010a0b810,
 	0xe60100000088e806,
 	0x19000000001c003c,
 /* 0x0218: gk110_rcp_f64 */
-/* 0x0218: gk110_rsq_f64 */
+	0xc00000058a1c0409,
+	0x407fffffff9c080d,
+	0x1480000060000000,
+	0xb3401c03fe9c0c1d,
+	0xe4c03c007f9c000e,
+	0x08a0ac10acb410ac,
+	0x120000004c20003c,
+	0x747fffffff9fc00e,
+	0xb4601fff801c021d,
+	0x120000000820003c,
+	0x21000400001c0404,
+	0x12000000381c003c,
+/* 0x0278: rcp_L4 */
+	0x203ff800001c0410,
+	0x08acb4aca010aca0,
+	0xb3281c00001c101d,
+	0x120000000c00003c,
+	0x223ff800001c0404,
+	0xe4c03c007f9c0002,
+	0x120000001c1c003c,
+/* 0x02b0: rcp_L5 */
+	0xb4601c00001c021d,
+	0x120000000c00003c,
+	0x08aca0bca080aca0,
+	0x213ff800001c0404,
+	0x12000000081c003c,
+/* 0x02d8: rcp_L6 */
+	0xc400021a801c0001,
+	0x740000001b1fc00e,
+/* 0x02e8: rcp_L2 */
+	0x85800000005c3c02,
+/* 0x02f0: rcp_L3 */
+	0xb3181c00001c0c1d,
+	0x12000000d000003c,
+	0x08a010aca010a0a0,
+	0xc00000058a1c0409,
+	0x204007ffff9c041c,
+	0x401ff800001c1c1d,
+	0xe4c03c00001c001a,
+	0xe5400c00031c3816,
+	0x84000000021c1412,
+	0x745fc000001fc002,
+	0x08b8a0b8aca8a0a0,
+	0xcc000000029c1016,
+	0xcac88000001c1415,
+	0xcc001000029c1002,
+	0xe5400000001c2c02,
+	0xe5410000031c3c1a,
+	0x741ff800001fc026,
+	0xe4c03c007f9c0022,
+	0x08a4a4a4a4a4a4a4,
+	0xdb802000001c1812,
+	0xdb800000021c0002,
+	0xdb802000001c1812,
+	0xdb800000021c0002,
+	0xdb802000001c1812,
+	0xdb800000021c0002,
+	0xdb802000001c1812,
+	0x08a0a0a0a0a0a080,
+	0xdb800000021c0002,
+	0x48000001ff9c0809,
+	0xe0800000019c0812,
+	0xc00000058a1c040d,
+	0xe0800000021c0c0e,
+	0x407fffffff9c0c09,
+	0xb3101c03ff1c081d,
+	0x0880aca0a0aca0ac,
+	0x120000000c20003c,
+	0xc24000000a1c1011,
+	0xe0800000009c1006,
+	0x120000003c1c003c,
+/* 0x0428: rcp_L7 */
+	0xb3681c03ff9c0c1d,
+	0x120000001420003c,
+	0x20400000001c0404,
+	0x088ca0948caca094,
+	0xe4c03c007f9c0002,
+	0x403ff800001c0405,
+	0x12000000201c003c,
+/* 0x0460: rcp_L8 */
+	0xb3501c00001c0c1d,
+	0x204007ffff9c0404,
+	0x741fe8000003c01e,
+	0x741ff0000023c01e,
+	0x080000b8b8b0a094,
+	0x40000800001c0405,
+	0xe4c03c007f9c001a,
+	0xe4000000031c0002,
+/* 0x04a0: rcp_end */
+	0x19000000001c003c,
+/* 0x04a8: gk110_rsq_f64 */
 	0x19000000001c003c,
 };
 
@@ -77,5 +167,5 @@ uint64_t gk110_builtin_offsets[] = {
 	0x0000000000000000,
 	0x00000000000000f0,
 	0x0000000000000218,
-	0x0000000000000218,
+	0x00000000000004a8,
 };
-- 
2.12.0



More information about the mesa-dev mailing list