Profiler Assembly and FloatingPoint Introduction to CUDA Programming

Profiler, Assembly, and Floating-Point Introduction to CUDA Programming Andreas Moshovos Winter 2009 Some material from: Wen-Mei Hwu and David Kirk NVIDIA Robert Strzodka, Dominik Göddeke, NVISION 08 presentation http: //www. mathematik. uni-dortmund. de/~goeddeke/pubs/NVISION 08 -long. pdf

The CUDA Profiler • Both GUI and command-line • Non-GUI control: – CUDA_PROFILE • set to 1 or 0 to enable or disable the profiler – CUDA_PROFILE_LOG • set to the name of the log file (will default to. /cuda_profile. log) – CUDA_PROFILE_CSV • set to 1 or 0 to enable or disable a comma separated version of the log – CUDA_PROFILE_CONFIG • specify a config file with up to 4 signals

Profiler Signals

Profiler Counters • Grid size X, Y • Block size X, Y, Z • Dyn smem per block: – Dynamic shared memory • Sta smem per block: – static shared memory • Reg per thread • Mem transfer dir – Direction: 0 host to device, 1 device to host • Mem transfer size – bytes

Interpreting Profiler Counters • Values represent events within a thread warp • Only targets one multiprocessor – Values will not correspond to the total number of warps launched for a particular kernel. – Launch enough thread blocks to ensure that the target multiprocessor is given a consistent percentage of the total work. • Values are best used to identify relative performance differences between non-optimized and optimized code – e. g. , make the number of non-coalesced loads go from some non-zero value to zero

CUDA Visual Profiler • Helps measure and find potential performance problem – GPU and CPU timing for all kernel invocations and memcpys – Time stamps • Access to hardware performance counters

Assembly

PTX: Assembly for NVIDIA GPUs • Parallel Thread e. Xecution • Virtual Assembly – Translated to actual machine code at runtime – Allows for different hardware implementations • Might enable additional optimizations – E. g. , %clock register to time blocks of code

Code Generation Flow • C/C++ Application ASM-level Library Programmer C/C++ Compiler PTX Code Parallel Thread e. Xecution (PTX) – Virtual Machine and ISA – Programming model – Execution resources and state • ISA – Instruction Set Architecture – Variable declarations – Instructions and operands PTX Code • Translator is an optimizing compiler – Translates PTX to Target code – Program install time PTX to Target • Translator Driver implements VM runtime – Coupled with Translator C G 80 … Target code GPU

How to See the PTX code • nvcc –keep – Produces. ptx and. cubin • nvcc --opencc-options -LIST: source=on

PTX Example CUDA float 4 me = gx[gtid]; me. x += me. y * me. z; PTX ld. global. v 4. f 32 {$f 1, $f 3, $f 5, $f 7}, [$r 9+0]; # 174 me. x += me. y * me. z; mad. f 32 $f 1, $f 5, $f 3, $f 1; Registers are virtual – The actual hardware registers are hidden from PTX

PTX Syntax Example

Another Example: CUDA Function • CUDA __device__ void interaction( float 4 b 0, float 4 b 1, float 3 *accel) { r. x = b 1. x - b 0. x; r. y = b 1. y - b 0. y; r. z = b 1. z - b 0. z; float dist. Sqr = r. x * r. x + r. y * r. y + r. z * r. z; float s = 1. 0 f/sqrt(dist. Sqr); accel->x += r. x * s; accel->y += r. y * s; accel->z += r. z * s; } • PTX sub. f 32 mul. f 32 add. f 32 rsqrt. f 32 mad. f 32 mov. f 32 $f 18, $f 15; $f 19, $f 3, $f 16; $f 20, $f 5, $f 17; $f 21, $f 18; $f 22, $f 19; $f 23, $f 20; $f 24, $f 21, $f 22; $f 25, $f 23, $f 24; $f 26, $f 25; $f 13, $f 18, $f 26, $f 13; $f 14, $f 13; $f 11, $f 19, $f 26, $f 11; $f 12, $f 11; $f 9, $f 20, $f 26, $f 9; $f 10, $f 9;

PTX Data types

Predicated Execution If (cond) – p = Evaluate cond – Branch not true After – Then Code • After: Then Code After Code – p = Evaluate cond – (p) Then Code – After Code

PTX Predicated Execution

Variable Declaration

Parameterized Variable Names • How to create 100 register “variables” • . reg. b 32 %r<100> • Declares %r 0 - %r 99

Addresses as Operands The value of x The value of tbl[12] The base address of tlb

Compiling a loop that calls a function - again • CUDA – sx is shared – mx, accel are local for (i = 0; i < K; i++) { if (i != thread. Idx. x) { interaction( sx[i], mx, &accel ); } } • PTX mov. s 32 $r 12, 0; $Lt_0_26: setp. eq. u 32 $p 1, $r 12, $r 5; @$p 1 bra $Lt_0_27; mul. lo. u 32 $r 13, $r 12, 16; add. u 32 $r 14, $r 13, $r 1; ld. shared. f 32 $f 15, [$r 14+0]; ld. shared. f 32 $f 16, [$r 14+4]; ld. shared. f 32 $f 17, [$r 14+8]; [func body from previous slide inlined here] $Lt_0_27: add. s 32 $r 12, 1; mov. s 32 $r 15, 128; setp. ne. s 32 $p 2, $r 15; @$p 2 bra $Lt_0_26;

Yet Another Example: SAXPY code cvt. u 32. u 16 $blockid, %ctaid. x; // Calculate i from thread/block IDs cvt. u 32. u 16 $blocksize, %ntid. x; cvt. u 32. u 16 $tid, %tid. x; mad 24. lo. u 32 $i, $blockid, $blocksize, $tid; ld. param. u 32 $n, [N]; // Nothing to do if n ≤ i setp. le. u 32 $p 1, $n, $i; @$p 1 bra $L_finish; mul. lo. u 32 $offset, $i, 4; // Load y[i] ld. param. u 32 $yaddr, [Y]; add. u 32 $yaddr, $offset; ld. global. f 32 $y_i, [$yaddr+0]; ld. param. u 32 $xaddr, [X]; // Load x[i] add. u 32 $xaddr, $offset; ld. global. f 32 $x_i, [$xaddr+0]; ld. param. f 32 $alpha, [ALPHA]; // Compute and store alpha*x[i] + y[i] mad. f 32 $y_i, $alpha, $x_i, $y_i; st. global. f 32 [$yaddr+0], $y_i; $L_finish: exit;

The %clock register • Real time clock cycle counter • How to read: – mov. u 32 $r 1, %clock; • Can be used to time code • It measures real time not just time spent executing this thread – If a thread is blocks time still elapses

PTX Reference • Please Read the PTX ISA specification – Posted under the handouts section

Occupancy Calculator • http: //developer. download. nvidia. com/compute/cuda/C UDA_Occupancy_calculator. xls • GPU Occupancy – – Active warps / max warps Threads/block Registers/thread Shared memory/block • Nvcc –cubin – code { name = my_kernel lmem = 0 smem = 24 reg = 5 bar = 0 bincode { �} const { �} }

Occupancy Calculator Example

Floating Point Considerations

Comparison of FP Capabilities G 80 SSE IBM Altivec Cell SPE Precision IEEE 754 Rounding modes for FADD and FMUL Round to nearest and round to zero All 4 IEEE, round to nearest, zero, inf Round to nearest only Round to zero/truncate only Denormal handling Flush to zero Supported, 1000’s of cycles Flush to zero Na. N support Yes Yes No Overflow and Infinity support Yes, only clamps to max norm Yes No, infinity Flags No Yes Some Square root Software only Hardware Software only Division Software only Hardware Software only Reciprocal estimate accuracy 24 bit 12 bit Reciprocal sqrt estimate accuracy 23 bit 12 bit log 2(x) and 2^x estimates accuracy 23 bit No 12 bit No

IEEE Floating Point Representation • A floating point binary number consists of three parts: – sign (S), exponent (E), and mantissa (M). – Each (S, E, M) pattern uniquely identifies a floating point number. • For each bit pattern, its IEEE floating-point value is derived as: – value = (-1)S * M * {2 E}, where 1. 0 ≤ M < 10. 0 B • The interpretation of S is simple: S=0 results in a positive number and S=1 a negative number.

Normalized Representation • Specifying that 1. 0 B ≤ M < 10. 0 B makes the mantissa value for each floating point number unique. – For example, the only one mantissa value allowed for 0. 5 D is M =1. 0 • 0. 5 D = 1. 0 B * 2 -1 – Neither 10. 0 B * 2 -2 nor 0. 1 B * 2 0 qualifies • Because all mantissa values are of the form 1. XX…, one can omit the “ 1. ” part in the representation. – The mantissa value of 0. 5 D in a 2 -bit mantissa is 00, which is derived by omitting “ 1. ” from 1. 00.

Exponent Representation • In an n-bits exponent representation, 2 n-1 -1 is added to its 2's complement representation to form its excess representation. – See Table for a 3 -bit exponent representation • A simple unsigned integer comparator can be used to compare the magnitude of two FP numbers • Symmetric range for +/exponents (111 reserved) 2’s complement Actual decimal Excess-3 000 0 011 001 1 100 010 2 101 011 3 110 100 (reserved pattern) 111 101 -3 000 110 -2 001 111 -1 010 E = represented E - BIAS

A Hypothetical 5 -bit Floating Point Representation • Assume 1 -bit S, 2 -bit E, and 2 -bit M – 0. 5 D = 1. 00 B * 2 -1 – 0. 5 D = 0 00 00, – where • S = 0, • E = 00 • M = (1. )00 2’s complement Actual decimal Excess-1 00 0 01 01 1 10 10 (reserved pattern) 11 11 -1 00

Representable Numbers • The representable numbers of a given format is the set of all numbers that can be exactly represented in the format. • • See Table for representable numbers of an unsigned 3 -bit integer format 00 0 0 00 1 1 01 2 0 01 3 1 10 4 0 10 5 1 11 6 0 -1 0 1 2 3 4 5 6 7 8 9 11 7 1

Hypothetical 5 -bit FP: Representable Numbers Abrupt underflow No-zero E M 00 00 2 -1 -(2 -1) 0 0 01 2 -1+1*2 -3 -(2 -1+1*2 -3) 0 0 1*2 -2 -1*2 -2 10 2 -1+2*2 -3 -(2 -1+2*2 -3) 0 0 2*2 -2 -2*2 -2 11 2 -1+3*2 -3 -(2 -1+3*2 -3) 0 0 3*2 -2 -3*2 -2 00 20 -(20) 01 20+1*2 -2 -(20+1*2 -2) 10 20+2*2 -2 -(20+2*2 -2) 11 20+3*2 -2 -(20+3*2 -2) 00 21 -(21) 01 21+1*2 -1 -(21+1*2 -1) 10 21+2*2 -1 -(21+2*2 -1) 11 21+3*2 -1 -(21+3*2 -1) 01 10 11 S=0 Gradual underflow S=1 Reserved pattern S=0 S=1

Flush To Zero • Treat all bit patterns with E=0 as 0. 0 – This takes away several representable numbers near zero and lump them all into 0. 0 – For a representation with large M, a large number of representable numbers will be removed. 0 1 2 3 4

Hypothetical 5 -bit FP: Representable Numbers Abrupt underflow No-zero E M 00 00 2 -1 -(2 -1) 0 0 01 2 -1+1*2 -3 -(2 -1+1*2 -3) 0 0 1*2 -2 -1*2 -2 10 2 -1+2*2 -3 -(2 -1+2*2 -3) 0 0 2*2 -2 -2*2 -2 11 2 -1+3*2 -3 -(2 -1+3*2 -3) 0 0 3*2 -2 -3*2 -2 00 20 -(20) 01 20+1*2 -2 -(20+1*2 -2) 10 20+2*2 -2 -(20+2*2 -2) 11 20+3*2 -2 -(20+3*2 -2) 00 21 -(21) 01 21+1*2 -1 -(21+1*2 -1) 10 21+2*2 -1 -(21+2*2 -1) 11 21+3*2 -1 -(21+3*2 -1) 01 10 11 S=0 Gradual underflow S=1 Reserved pattern S=0 S=1

Denormalized Numbers • The actual method adopted by the IEEE standard is called denormalized numbers or gradual underflow. – The method relaxes the normalization requirement for numbers very close to 0. – whenever E=0, the mantissa is no longer assumed to be of the form 1. XX. Rather, it is assumed to be 0. XX. In general, if the n-bit exponent is 0, the value is • 0. M * 2 - 2 ^(n-1) + 2 0 1 2 3

Hypothetical 5 -bit FP: Representable Numbers Abrupt underflow No-zero E M 00 00 2 -1 -(2 -1) 0 0 01 2 -1+1*2 -3 -(2 -1+1*2 -3) 0 0 1*2 -2 -1*2 -2 10 2 -1+2*2 -3 -(2 -1+2*2 -3) 0 0 2*2 -2 -2*2 -2 11 2 -1+3*2 -3 -(2 -1+3*2 -3) 0 0 3*2 -2 -3*2 -2 00 20 -(20) 01 20+1*2 -2 -(20+1*2 -2) 10 20+2*2 -2 -(20+2*2 -2) 11 20+3*2 -2 -(20+3*2 -2) 00 21 -(21) 01 21+1*2 -1 -(21+1*2 -1) 10 21+2*2 -1 -(21+2*2 -1) 11 21+3*2 -1 -(21+3*2 -1) 01 10 11 S=0 Gradual underflow S=1 Reserved pattern S=0 S=1

Floating Point Numbers • As the exponent gets larger • The distance between two representable numbers increases

Arithmetic Instruction Throughput • int and float add, shift, min, max and float mul, mad: 4 cycles per warp – int multiply (*) is by default 32 -bit • requires multiple cycles / warp – Use __mul 24() / __umul 24() intrinsics for 4 -cycle 24 -bit int multiply – For G 80, for G 20 should be OK • Integer divide and modulo are expensive – Compiler will convert literal power-of-2 divides to shifts – Be explicit in cases where compiler can’t tell that divisor is a power of 2 – Useful trick: foo % n == foo & (n-1) if n is a power of 2

Arithmetic Instruction Throughput • Reciprocal, reciprocal square root, sin/cos, log, exp: 16 cycles per warp – These are the versions prefixed with “__” – Examples: __rcp(), __sin(), __exp() • Other functions are combinations of the above – y / x == rcp(x) * y == 20 cycles per warp – sqrt(x) == rcp(rsqrt(x)) == 32 cycles per warp

Runtime Math Library • There are two types of runtime math operations – __func(): direct mapping to hardware ISA • Fast but low accuracy (see prog. guide for details) • Examples: __sin(x), __exp(x), __pow(x, y) – func() : compile to multiple instructions • Slower but higher accuracy (5 ulp, units in the least place, or less) • Examples: sin(x), exp(x), pow(x, y) • The -use_fast_math compiler option forces every func() to compile to __func()

• Make float-safe! G 20 has your doubleprogram precision support – G 80 is single-precision only – Double precision has additional performance cost • Only one unit per multiprocessor – Careless use of double or undeclared types may run more slowly on G 80+ • Important to be float-safe (be explicit whenever you want single precision) to avoid using double precision where it is not needed – Add ‘f’ specifier on float literals: • foo = bar * 0. 123; • foo = bar * 0. 123 f; // double assumed // float explicit – Use float version of standard library functions • foo = sin(bar); • foo = sinf(bar); // double assumed // single precision explicit

Deviations from IEEE-754 • Addition and Multiplication are IEEE 754 compliant – Maximum 0. 5 ulp (units in the least place) error • However, often combined into multiply-add (FMAD) – Intermediate result is truncated • • Division is non-compliant (2 ulp) Not all rounding modes are supported Denormalized numbers are not supported No mechanism to detect floating-point exceptions

Units in the Last Place Error • If the result of a FP computation is: – 3. 12 x 10^-2 = 0. 0312 • But the answer when computed to infinite precision is: – -0. 0312159 • Then ulp is: – 0. 0314 – 0. 0312 = 0. 159 • For binary representations the maximum ulp is 0. 5 – Round to nearest number