(no commit message)
[libreriscv.git] / ztrans_proposal.mdwn
1 # Zftrans - transcendental operations
2
3 With thanks to:
4
5 * Jacob Lifshay
6 * Dan Petroski
7 * Mitch Alsup
8 * Allen Baum
9 * Andrew Waterman
10 * Luis Vitorio Cargnini
11
12 [[!toc levels=2]]
13
14 See:
15
16 * <http://bugs.libre-riscv.org/show_bug.cgi?id=127>
17 * <https://www.khronos.org/registry/spir-v/specs/unified1/OpenCL.ExtendedInstructionSet.100.html>
18 * Discussion: <http://lists.libre-riscv.org/pipermail/libre-riscv-dev/2019-August/002342.html>
19 * [[rv_major_opcode_1010011]] for opcode listing.
20 * [[zfpacc_proposal]] for accuracy settings proposal
21
22 Extension subsets:
23
24 * **Zftrans**: standard transcendentals (best suited to 3D)
25 * **ZftransExt**: extra functions (useful, not generally needed for 3D,
26 can be synthesised using Ztrans)
27 * **Ztrigpi**: trig. xxx-pi sinpi cospi tanpi
28 * **Ztrignpi**: trig non-xxx-pi sin cos tan
29 * **Zarctrigpi**: arc-trig. a-xxx-pi: atan2pi asinpi acospi
30 * **Zarctrignpi**: arc-trig. non-a-xxx-pi: atan2, asin, acos
31 * **Zfhyp**: hyperbolic/inverse-hyperbolic. sinh, cosh, tanh, asinh,
32 acosh, atanh (can be synthesised - see below)
33 * **ZftransAdv**: much more complex to implement in hardware
34 * **Zfrsqrt**: Reciprocal square-root.
35
36 Minimum recommended requirements for 3D: Zftrans, Ztrignpi,
37 Zarctrignpi, with Ztrigpi and Zarctrigpi as augmentations.
38
39 Minimum recommended requirements for Mobile-Embedded 3D: Ztrignpi, Zftrans, with Ztrigpi as an augmentation.
40
41 # TODO:
42
43 * Decision on accuracy, moved to [[zfpacc_proposal]]
44 <http://lists.libre-riscv.org/pipermail/libre-riscv-dev/2019-August/002355.html>
45 * Errors **MUST** be repeatable.
46 * How about four Platform Specifications? 3DUNIX, UNIX, 3DEmbedded and Embedded?
47 <http://lists.libre-riscv.org/pipermail/libre-riscv-dev/2019-August/002361.html>
48 Accuracy requirements for dual (triple) purpose implementations must
49 meet the higher standard.
50 * Reciprocal Square-root is in its own separate extension (Zfrsqrt) as
51 it is desirable on its own by other implementors. This to be evaluated.
52
53 # Requirements <a name="requirements"></a>
54
55 This proposal is designed to meet a wide range of extremely diverse needs,
56 allowing implementors from all of them to benefit from the tools and hardware
57 cost reductions associated with common standards adoption.
58
59 **There are *four* different, disparate platform's needs (two new)**:
60
61 * 3D Embedded Platform (new)
62 * Embedded Platform
63 * 3D UNIX Platform (new)
64 * UNIX Platform
65
66 **The use-cases are**:
67
68 * 3D GPUs
69 * Numerical Computation
70 * (Potentially) A.I. / Machine-learning (1)
71
72 (1) although approximations suffice in this field, making it more likely
73 to use a custom extension. High-end ML would inherently definitely
74 be excluded.
75
76 **The power and die-area requirements vary from**:
77
78 * Ultra-low-power (smartwatches where GPU power budgets are in milliwatts)
79 * Mobile-Embedded (good performance with high efficiency for battery life)
80 * Desktop Computing
81 * Server / HPC (2)
82
83 (2) Supercomputing is left out of the requirements as it is traditionally
84 covered by Supercomputer Vectorisation Standards (such as RVV).
85
86 **The software requirements are**:
87
88 * Full public integration into GNU math libraries (libm)
89 * Full public integration into well-known Numerical Computation systems (numpy)
90 * Full public integration into upstream GNU and LLVM Compiler toolchains
91 * Full public integration into Khronos OpenCL SPIR-V compatible Compilers
92 seeking public Certification and Endorsement from the Khronos Group
93 under their Trademarked Certification Programme.
94
95 **The "contra"-requirements are**:
96
97 * NOT for use with RVV (RISC-V Vector Extension). These are *scalar* opcodes.
98 Ultra Low Power Embedded platforms (smart watches) are sufficiently
99 resource constrained that Vectorisation (of any kind) is likely to be
100 unnecessary and inappropriate.
101 * The requirements are **not** for the purposes of developing a full custom
102 proprietary GPU with proprietary firmware driven by *hardware* centric
103 optimised design decisions as a priority over collaboration.
104 * A full custom proprietary GPU ASIC Manufacturer *may* benefit from
105 this proposal however the fact that they typically develop proprietary
106 software that is not shared with the rest of the community likely to
107 use this proposal means that they have completely different needs.
108 * This proposal is for *sharing* of effort in reducing development costs
109
110 # Requirements Analysis <a name="requirements_analysis"></a>
111
112 **Platforms**:
113
114 3D Embedded will require significantly less accuracy and will need to make
115 power budget and die area compromises that other platforms (including Embedded)
116 will not need to make.
117
118 3D UNIX Platform has to be performance-price-competitive: subtly-reduced
119 accuracy in FP32 is acceptable where, conversely, in the UNIX Platform,
120 IEEE754 compliance is a hard requirement that would compromise power
121 and efficiency on a 3D UNIX Platform.
122
123 Even in the Embedded platform, IEEE754 interoperability is beneficial,
124 where if it was a hard requirement the 3D Embedded platform would be severely
125 compromised in its ability to meet the demanding power budgets of that market.
126
127 Thus, learning from the lessons of
128 [SIMD considered harmful](https://www.sigarch.org/simd-instructions-considered-harmful/)
129 this proposal works in conjunction with the [[zfpacc_proposal]], so as
130 not to overburden the OP32 ISA space with extra "reduced-accuracy" opcodes.
131
132 **Use-cases**:
133
134 There really is little else in the way of suitable markets. 3D GPUs
135 have extremely competitive power-efficiency and power-budget requirements
136 that are completely at odds with the other market at the other end of
137 the spectrum: Numerical Computation.
138
139 Interoperability in Numerical Computation is absolutely critical: it
140 implies (correlates directly with) IEEE754 compliance. However full
141 IEEE754 compliance automatically and inherently penalises a GPU on
142 performance and die area, where accuracy is simply just not necessary.
143
144 To meet the needs of both markets, the two new platforms have to be created,
145 and [[zfpacc_proposal]] is a critical dependency. Runtime selection of
146 FP accuracy allows an implementation to be "Hybrid" - cover UNIX IEEE754
147 compliance *and* 3D performance in a single ASIC.
148
149 **Power and die-area requirements**:
150
151 This is where the conflicts really start to hit home.
152
153 A "Numerical High performance only" proposal (suitable for Server / HPC
154 only) would customise and target the Extension based on a quantitative
155 analysis of the value of certain opcodes *for HPC only*. It would
156 conclude, reasonably and rationally, that it is worthwhile adding opcodes
157 to RVV as parallel Vector operations, and that further discussion of
158 the matter is pointless.
159
160 A "Proprietary GPU effort" (even one that was intended for publication
161 of its API through, for example, a public libre-licensed Vulkan SPIR-V
162 Compiler) would conclude, reasonably and rationally, that, likewise, the
163 opcodes were best suited to be added to RVV, and, further, that their
164 requirements conflict with the HPC world, due to the reduced accuracy.
165 This on the basis that the silicon die area required for IEEE754 is far
166 greater than that needed for reduced-accuracy, and thus their product
167 would be completely unacceptable in the market if it had to meet IEEE754,
168 unnecessarily.
169
170 An "Embedded 3D" GPU has radically different performance, power
171 and die-area requirements (and may even target SoftCores in FPGA).
172 Sharing of the silicon to cover multi-function uses (CORDIC for example)
173 is absolutely essential in order to keep cost and power down, and high
174 performance simply is not. Multi-cycle FSMs instead of pipelines may
175 be considered acceptable, and so on. Subsets of functionality are
176 also essential.
177
178 An "Embedded Numerical" platform has requirements that are separate and
179 distinct from all of the above!
180
181 Mobile Computing needs (tablets, smartphones) again pull in a different
182 direction: high performance, reasonable accuracy, but efficiency is
183 critical. Screen sizes are not at the 4K range: they are within the
184 800x600 range at the low end (320x240 at the extreme budget end), and
185 only the high-performance smartphones and tablets provide 1080p (1920x1080).
186 With lower resolution, accuracy compromises are possible which the Desktop
187 market (4k and soon to be above) would find unacceptable.
188
189 Meeting these disparate markets may be achieved, again, through
190 [[zfpacc_proposal]], by subdividing into four platforms, yet, in addition
191 to that, subdividing the extension into subsets that best suit the different
192 market areas.
193
194 **Software requirements**:
195
196 A "custom" extension is developed in near-complete isolation from the
197 rest of the RISC-V Community. Cost savings to the Corporation are
198 large, with no direct beneficial feedback to (or impact on) the rest
199 of the RISC-V ecosystem.
200
201 However given that 3D revolves around Standards - DirectX, Vulkan, OpenGL,
202 OpenCL - users have much more influence than first appears. Compliance
203 with these standards is critical as the userbase (Games writers,
204 scientific applications) expects not to have to rewrite extremely large
205 and costly codebases to conform with *non-standards-compliant* hardware.
206
207 Therefore, compliance with public APIs (Vulkan, OpenCL, OpenGL, DirectX)
208 is paramount, and compliance with Trademarked Standards is critical.
209 Any deviation from Trademarked Standards means that an implementation
210 may not be sold and also make a claim of being, for example, "Vulkan
211 compatible".
212
213 This in turn reinforces and makes a hard requirement a need for public
214 compliance with such standards, over-and-above what would otherwise be
215 set by a RISC-V Standards Development Process, including both the
216 software compliance and the knock-on implications that has for hardware.
217
218 **Collaboration**:
219
220 The case for collaboration on any Extension is already well-known.
221 In this particular case, the precedent for inclusion of Transcendentals
222 in other ISAs, both from Graphics and High-performance Computing, has
223 these primitives well-established in high-profile software libraries and
224 compilers in both GPU and HPC Computer Science divisions. Collaboration
225 and shared public compliance with those standards brooks no argument.
226
227 The combined requirements of collaboration and multi accuracy requirements
228 mean that *overall this proposal is categorically and wholly unsuited
229 to relegation of "custom" status*.
230
231 # Quantitative Analysis <a name="analysis"></a>
232
233 This is extremely challenging. Normally, an Extension would require full,
234 comprehensive and detailed analysis of every single instruction, for every
235 single possible use-case, in every single market. The amount of silicon
236 area required would be balanced against the benefits of introducing extra
237 opcodes, as well as a full market analysis performed to see which divisions
238 of Computer Science benefit from the introduction of the instruction,
239 in each and every case.
240
241 With 34 instructions, four possible Platforms, and sub-categories of
242 implementations even within each Platform, over 136 separate and distinct
243 analyses is not a practical proposition.
244
245 A little more intelligence has to be applied to the problem space,
246 to reduce it down to manageable levels.
247
248 Fortunately, the subdivision by Platform, in combination with the
249 identification of only two primary markets (Numerical Computation and
250 3D), means that the logical reasoning applies *uniformly* and broadly
251 across *groups* of instructions rather than individually, making it a primarily
252 hardware-centric and accuracy-centric decision-making process.
253
254 In addition, hardware algorithms such as CORDIC can cover such a wide
255 range of operations (simply by changing the input parameters) that the
256 normal argument of compromising and excluding certain opcodes because they
257 would significantly increase the silicon area is knocked down.
258
259 However, CORDIC, whilst space-efficient, and thus well-suited to
260 Embedded, is an old iterative algorithm not well-suited to High-Performance
261 Computing or Mid to High-end GPUs, where commercially-competitive
262 FP32 pipeline lengths are only around 5 stages.
263
264 Not only that, but some operations such as LOG1P, which would normally
265 be excluded from one market (due to there being an alternative macro-op
266 fused sequence replacing it) are required for other markets due to
267 the higher accuracy obtainable at the lower range of input values when
268 compared to LOG(1+P).
269
270 (Thus we start to see why "proprietary" markets are excluded from this
271 proposal, because "proprietary" markets would make *hardware*-driven
272 optimisation decisions that would be completely inappropriate for a
273 common standard).
274
275 ATAN and ATAN2 is another example area in which one market's needs
276 conflict directly with another: the only viable solution, without compromising
277 one market to the detriment of the other, is to provide both opcodes
278 and let implementors make the call as to which (or both) to optimise,
279 at the *hardware* level.
280
281 Likewise it is well-known that loops involving "0 to 2 times pi", often
282 done in subdivisions of powers of two, are costly to do because they
283 involve floating-point multiplication by PI in each and every loop.
284 3D GPUs solved this by providing SINPI variants which range from 0 to 1
285 and perform the multiply *inside* the hardware itself. In the case of
286 CORDIC, it turns out that the multiply by PI is not even needed (is a
287 loop invariant magic constant).
288
289 However, some markets may not wish to *use* CORDIC, for reasons mentioned
290 above, and, again, one market would be penalised if SINPI was prioritised
291 over SIN, or vice-versa.
292
293 In essence, then, even when only the two primary markets (3D and
294 Numerical Computation) have been identified, this still leaves two
295 (three) diametrically-opposed *accuracy* sub-markets as the prime
296 conflict drivers:
297
298 * Embedded Ultra Low Power
299 * IEEE754 compliance
300 * Khronos Vulkan compliance
301
302 Thus the best that can be done is to use Quantitative Analysis to work
303 out which "subsets" - sub-Extensions - to include, provide an additional
304 "accuracy" extension, be as "inclusive" as possible, and thus allow
305 implementors to decide what to add to their implementation, and how best
306 to optimise them.
307
308 This approach *only* works due to the uniformity of the function space,
309 and is **not** an appropriate methodology for use in other Extensions
310 with huge (non-uniform) market diversity even with similarly large
311 numbers of potential opcodes. BitManip is the perfect counter-example.
312
313 # Proposed Opcodes vs Khronos OpenCL vs IEEE754-2019<a name="khronos_equiv"></a>
314
315 This list shows the (direct) equivalence between proposed opcodes,
316 their Khronos OpenCL equivalents, and their IEEE754-2019 equivalents.
317 98% of the opcodes in this proposal that are in the IEEE754-2019 standard
318 are present in the Khronos Extended Instruction Set.
319
320 For RISCV opcode encodings see
321 [[rv_major_opcode_1010011]]
322
323 See
324 <https://www.khronos.org/registry/spir-v/specs/unified1/OpenCL.ExtendedInstructionSet.100.html>
325 and <https://ieeexplore.ieee.org/document/8766229>
326
327 * Special FP16 opcodes are *not* being proposed, except by indirect / inherent
328 use of the "fmt" field that is already present in the RISC-V Specification.
329 * "Native" opcodes are *not* being proposed: implementors will be expected
330 to use the (equivalent) proposed opcode covering the same function.
331 * "Fast" opcodes are *not* being proposed, because the Khronos Specification
332 fast\_length, fast\_normalise and fast\_distance OpenCL opcodes require
333 vectors (or can be done as scalar operations using other RISC-V instructions).
334
335 The OpenCL FP32 opcodes are **direct** equivalents to the proposed opcodes.
336 Deviation from conformance with the Khronos Specification - including the
337 Khronos Specification accuracy requirements - is not an option, as it
338 results in non-compliance, and the vendor may not use the Trademarked words
339 "Vulkan" etc. in conjunction with their product.
340
341 IEEE754-2019 Table 9.1 lists "additional mathematical operations".
342 Interestingly the only functions missing when compared to OpenCL are
343 compound, exp2m1, exp10m1, log2p1, log10p1, pown (integer power) and powr.
344
345 [[!table data="""
346 opcode | OpenCL FP32 | OpenCL FP16 | OpenCL native | OpenCL fast | IEEE754 |
347 FSIN | sin | half\_sin | native\_sin | NONE | sin |
348 FCOS | cos | half\_cos | native\_cos | NONE | cos |
349 FTAN | tan | half\_tan | native\_tan | NONE | tan |
350 NONE (1) | sincos | NONE | NONE | NONE | NONE |
351 FASIN | asin | NONE | NONE | NONE | asin |
352 FACOS | acos | NONE | NONE | NONE | acos |
353 FATAN | atan | NONE | NONE | NONE | atan |
354 FSINPI | sinpi | NONE | NONE | NONE | sinPi |
355 FCOSPI | cospi | NONE | NONE | NONE | cosPi |
356 FTANPI | tanpi | NONE | NONE | NONE | tanPi |
357 FASINPI | asinpi | NONE | NONE | NONE | asinPi |
358 FACOSPI | acospi | NONE | NONE | NONE | acosPi |
359 FATANPI | atanpi | NONE | NONE | NONE | atanPi |
360 FSINH | sinh | NONE | NONE | NONE | sinh |
361 FCOSH | cosh | NONE | NONE | NONE | cosh |
362 FTANH | tanh | NONE | NONE | NONE | tanh |
363 FASINH | asinh | NONE | NONE | NONE | asinh |
364 FACOSH | acosh | NONE | NONE | NONE | acosh |
365 FATANH | atanh | NONE | NONE | NONE | atanh |
366 FATAN2 | atan2 | NONE | NONE | NONE | atan2 |
367 FATAN2PI | atan2pi | NONE | NONE | NONE | atan2pi |
368 FRSQRT | rsqrt | half\_rsqrt | native\_rsqrt | NONE | rSqrt |
369 FCBRT | cbrt | NONE | NONE | NONE | NONE (2) |
370 FEXP2 | exp2 | half\_exp2 | native\_exp2 | NONE | exp2 |
371 FLOG2 | log2 | half\_log2 | native\_log2 | NONE | log2 |
372 FEXPM1 | expm1 | NONE | NONE | NONE | expm1 |
373 FLOG1P | log1p | NONE | NONE | NONE | logp1 |
374 FEXP | exp | half\_exp | native\_exp | NONE | exp |
375 FLOG | log | half\_log | native\_log | NONE | log |
376 FEXP10 | exp10 | half\_exp10 | native\_exp10 | NONE | exp10 |
377 FLOG10 | log10 | half\_log10 | native\_log10 | NONE | log10 |
378 FPOW | pow | NONE | NONE | NONE | pow |
379 FPOWN | pown | NONE | NONE | NONE | pown |
380 FPOWR | powr | NONE | NONE | NONE | powr |
381 FROOTN | rootn | NONE | NONE | NONE | rootn |
382 FHYPOT | hypot | NONE | NONE | NONE | hypot |
383 FRECIP | NONE | half\_recip | native\_recip | NONE | NONE (3) |
384 NONE | NONE | NONE | NONE | NONE | compound |
385 NONE | NONE | NONE | NONE | NONE | exp2m1 |
386 NONE | NONE | NONE | NONE | NONE | exp10m1 |
387 NONE | NONE | NONE | NONE | NONE | log2p1 |
388 NONE | NONE | NONE | NONE | NONE | log10p1 |
389 """]]
390
391 Note (1) FSINCOS is macro-op fused (see below).
392
393 Note (2) synthesised in IEEE754-2019 as "pown(x, 3)"
394
395 Note (3) synthesised in IEEE754-2019 using "1.0 / x"
396
397 ## List of 2-arg opcodes
398
399 [[!table data="""
400 opcode | Description | pseudocode | Extension |
401 FATAN2 | atan2 arc tangent | rd = atan2(rs2, rs1) | Zarctrignpi |
402 FATAN2PI | atan2 arc tangent / pi | rd = atan2(rs2, rs1) / pi | Zarctrigpi |
403 FPOW | x power of y | rd = pow(rs1, rs2) | ZftransAdv |
404 FPOWN | x power of n (n int) | rd = pow(rs1, rs2) | ZftransAdv |
405 FPOWR | x power of y (x +ve) | rd = exp(rs1 log(rs2)) | ZftransAdv |
406 FROOTN | x power 1/n (n integer)| rd = pow(rs1, 1/rs2) | ZftransAdv |
407 FHYPOT | hypotenuse | rd = sqrt(rs1^2 + rs2^2) | ZftransAdv |
408 """]]
409
410 ## List of 1-arg transcendental opcodes
411
412 [[!table data="""
413 opcode | Description | pseudocode | Extension |
414 FRSQRT | Reciprocal Square-root | rd = sqrt(rs1) | Zfrsqrt |
415 FCBRT | Cube Root | rd = pow(rs1, 1.0 / 3) | ZftransAdv |
416 FRECIP | Reciprocal | rd = 1.0 / rs1 | Zftrans |
417 FEXP2 | power-of-2 | rd = pow(2, rs1) | Zftrans |
418 FLOG2 | log2 | rd = log(2. rs1) | Zftrans |
419 FEXPM1 | exponential minus 1 | rd = pow(e, rs1) - 1.0 | ZftransExt |
420 FLOG1P | log plus 1 | rd = log(e, 1 + rs1) | ZftransExt |
421 FEXP | exponential | rd = pow(e, rs1) | ZftransExt |
422 FLOG | natural log (base e) | rd = log(e, rs1) | ZftransExt |
423 FEXP10 | power-of-10 | rd = pow(10, rs1) | ZftransExt |
424 FLOG10 | log base 10 | rd = log(10, rs1) | ZftransExt |
425 """]]
426
427 ## List of 1-arg trigonometric opcodes
428
429 [[!table data="""
430 opcode | Description | pseudo-code | Extension |
431 FSIN | sin (radians) | rd = sin(rs1) | Ztrignpi |
432 FCOS | cos (radians) | rd = cos(rs1) | Ztrignpi |
433 FTAN | tan (radians) | rd = tan(rs1) | Ztrignpi |
434 FASIN | arcsin (radians) | rd = asin(rs1) | Zarctrignpi |
435 FACOS | arccos (radians) | rd = acos(rs1) | Zarctrignpi |
436 FATAN | arctan (radians) | rd = atan(rs1) | Zarctrignpi |
437 FSINPI | sin times pi | rd = sin(pi * rs1) | Ztrigpi |
438 FCOSPI | cos times pi | rd = cos(pi * rs1) | Ztrigpi |
439 FTANPI | tan times pi | rd = tan(pi * rs1) | Ztrigpi |
440 FASINPI | arcsin / pi | rd = asin(rs1) / pi | Zarctrigpi |
441 FACOSPI | arccos / pi | rd = acos(rs1) / pi | Zarctrigpi |
442 FATANPI | arctan / pi | rd = atan(rs1) / pi | Zarctrigpi |
443 FSINH | hyperbolic sin (radians) | rd = sinh(rs1) | Zfhyp |
444 FCOSH | hyperbolic cos (radians) | rd = cosh(rs1) | Zfhyp |
445 FTANH | hyperbolic tan (radians) | rd = tanh(rs1) | Zfhyp |
446 FASINH | inverse hyperbolic sin | rd = asinh(rs1) | Zfhyp |
447 FACOSH | inverse hyperbolic cos | rd = acosh(rs1) | Zfhyp |
448 FATANH | inverse hyperbolic tan | rd = atanh(rs1) | Zfhyp |
449 """]]
450
451 # Subsets
452
453 The full set is based on the Khronos OpenCL opcodes. If implemented
454 entirely it would be too much for both Embedded and also 3D.
455
456 The subsets are organised by hardware complexity, need (3D, HPC), however
457 due to synthesis producing inaccurate results at the range limits,
458 the less common subsets are still required for IEEE754 HPC.
459
460 MALI Midgard, an embedded / mobile 3D GPU, for example only has the
461 following opcodes:
462
463 E8 - fatan_pt2
464 F0 - frcp (reciprocal)
465 F2 - frsqrt (inverse square root, 1/sqrt(x))
466 F3 - fsqrt (square root)
467 F4 - fexp2 (2^x)
468 F5 - flog2
469 F6 - fsin1pi
470 F7 - fcos1pi
471 F9 - fatan_pt1
472
473 These in FP32 and FP16 only: no FP32 hardware, at all.
474
475 Vivante Embedded/Mobile 3D (etnaviv <https://github.com/laanwj/etna_viv/blob/master/rnndb/isa.xml>) only has the following:
476
477 sin, cos2pi
478 cos, sin2pi
479 log2, exp
480 sqrt and rsqrt
481 recip.
482
483 It also has fast variants of some of these, as a CSR Mode.
484
485 AMD's R600 GPU (R600\_Instruction\_Set\_Architecture.pdf) and the
486 RDNA ISA (RDNA\_Shader\_ISA\_5August2019.pdf, Table 22, Section 6.3) have:
487
488 COS2PI (appx)
489 EXP2
490 LOG (IEEE754)
491 RECIP
492 RSQRT
493 SQRT
494 SIN2PI (appx)
495
496 AMD RDNA has F16 and F32 variants of all the above, and also has F64
497 variants of SQRT, RSQRT and RECIP. It is interesting that even the
498 modern high-end AMD GPU does not have TAN or ATAN, where MALI Midgard
499 does.
500
501 Also a general point, that customised optimised hardware targetting
502 FP32 3D with less accuracy simply can neither be used for IEEE754 nor
503 for FP64 (except as a starting point for hardware or software driven
504 Newton Raphson or other iterative method).
505
506 Also in cost/area sensitive applications even the extra ROM lookup tables
507 for certain algorithms may be too costly.
508
509 These wildly differing and incompatible driving factors lead to the
510 subset subdivisions, below.
511
512 ## Transcendental Subsets
513
514 ### Zftrans
515
516 LOG2 EXP2 RECIP RSQRT
517
518 Zftrans contains the minimum standard transcendentals best suited to
519 3D. They are also the minimum subset for synthesising log10, exp10,
520 exp1m, log1p, the hyperbolic trigonometric functions sinh and so on.
521
522 They are therefore considered "base" (essential) transcendentals.
523
524 ### ZftransExt
525
526 LOG, EXP, EXP10, LOG10, LOGP1, EXP1M
527
528 These are extra transcendental functions that are useful, not generally
529 needed for 3D, however for Numerical Computation they may be useful.
530
531 Although they can be synthesised using Ztrans (LOG2 multiplied
532 by a constant), there is both a performance penalty as well as an
533 accuracy penalty towards the limits, which for IEEE754 compliance is
534 unacceptable. In particular, LOG(1+rs1) in hardware may give much better
535 accuracy at the lower end (very small rs1) than LOG(rs1).
536
537 Their forced inclusion would be inappropriate as it would penalise
538 embedded systems with tight power and area budgets. However if they
539 were completely excluded the HPC applications would be penalised on
540 performance and accuracy.
541
542 Therefore they are their own subset extension.
543
544 ### Zfhyp
545
546 SINH, COSH, TANH, ASINH, ACOSH, ATANH
547
548 These are the hyperbolic/inverse-hyperbolic functions. Their use in 3D is limited.
549
550 They can all be synthesised using LOG, SQRT and so on, so depend
551 on Zftrans. However, once again, at the limits of the range, IEEE754
552 compliance becomes impossible, and thus a hardware implementation may
553 be required.
554
555 HPC and high-end GPUs are likely markets for these.
556
557 ### ZftransAdv
558
559 CBRT, POW, POWN, POWR, ROOTN
560
561 These are simply much more complex to implement in hardware, and typically
562 will only be put into HPC applications.
563
564 * **Zfrsqrt**: Reciprocal square-root.
565
566 ## Trigonometric subsets
567
568 ### Ztrigpi vs Ztrignpi
569
570 * **Ztrigpi**: SINPI COSPI TANPI
571 * **Ztrignpi**: SIN COS TAN
572
573 Ztrignpi are the basic trigonometric functions through which all others
574 could be synthesised, and they are typically the base trigonometrics
575 provided by GPUs for 3D, warranting their own subset.
576
577 In the case of the Ztrigpi subset, these are commonly used in for loops
578 with a power of two number of subdivisions, and the cost of multiplying
579 by PI inside each loop (or cumulative addition, resulting in cumulative
580 errors) is not acceptable.
581
582 In for example CORDIC the multiplication by PI may be moved outside of
583 the hardware algorithm as a loop invariant, with no power or area penalty.
584
585 Again, therefore, if SINPI (etc.) were excluded, programmers would be penalised by being forced to divide by PI in some circumstances. Likewise if SIN were excluded, programmers would be penaslised by being forced to *multiply* by PI in some circumstances.
586
587 Thus again, a slightly different application of the same general argument applies to give Ztrignpi and
588 Ztrigpi as subsets. 3D GPUs will almost certainly provide both.
589
590 ### Zarctrigpi and Zarctrignpi
591
592 * **Zarctrigpi**: ATAN2PI ASINPI ACOSPI
593 * **Zarctrignpi**: ATAN2 ACOS ASIN
594
595 These are extra trigonometric functions that are useful in some
596 applications, but even for 3D GPUs, particularly embedded and mobile class
597 GPUs, they are not so common and so are typically synthesised, there.
598
599 Although they can be synthesised using Ztrigpi and Ztrignpi, there is,
600 once again, both a performance penalty as well as an accuracy penalty
601 towards the limits, which for IEEE754 compliance is unacceptable, yet
602 is acceptable for 3D.
603
604 Therefore they are their own subset extensions.
605
606 # Synthesis, Pseudo-code ops and macro-ops
607
608 The pseudo-ops are best left up to the compiler rather than being actual
609 pseudo-ops, by allocating one scalar FP register for use as a constant
610 (loop invariant) set to "1.0" at the beginning of a function or other
611 suitable code block.
612
613 * FSINCOS - fused macro-op between FSIN and FCOS (issued in that order).
614 * FSINCOSPI - fused macro-op between FSINPI and FCOSPI (issued in that order).
615
616 FATANPI example pseudo-code:
617
618 lui t0, 0x3F800 // upper bits of f32 1.0
619 fmv.x.s ft0, t0
620 fatan2pi.s rd, rs1, ft0
621
622 Hyperbolic function example (obviates need for Zfhyp except for
623 high-performance or correctly-rounding):
624
625 ASINH( x ) = ln( x + SQRT(x**2+1))
626
627 # Evaluation and commentary
628
629 This section will move later to discussion.
630
631 ## Reciprocal
632
633 Used to be an alias. Some implementors may wish to implement divide as
634 y times recip(x).
635
636 Others may have shared hardware for recip and divide, others may not.
637
638 To avoid penalising one implementor over another, recip stays.
639
640 ## To evaluate: should LOG be replaced with LOG1P (and EXP with EXPM1)?
641
642 RISC principle says "exclude LOG because it's covered by LOGP1 plus an ADD".
643 Research needed to ensure that implementors are not compromised by such
644 a decision
645 <http://lists.libre-riscv.org/pipermail/libre-riscv-dev/2019-August/002358.html>
646
647 > > correctly-rounded LOG will return different results than LOGP1 and ADD.
648 > > Likewise for EXP and EXPM1
649
650 > ok, they stay in as real opcodes, then.
651
652 ## ATAN / ATAN2 commentary
653
654 Discussion starts here:
655 <http://lists.libre-riscv.org/pipermail/libre-riscv-dev/2019-August/002470.html>
656
657 from Mitch Alsup:
658
659 would like to point out that the general implementations of ATAN2 do a
660 bunch of special case checks and then simply call ATAN.
661
662 double ATAN2( double y, double x )
663 { // IEEE 754-2008 quality ATAN2
664
665 // deal with NANs
666 if( ISNAN( x ) ) return x;
667 if( ISNAN( y ) ) return y;
668
669 // deal with infinities
670 if( x == +∞ && |y|== +∞ ) return copysign( π/4, y );
671 if( x == +∞ ) return copysign( 0.0, y );
672 if( x == -∞ && |y|== +∞ ) return copysign( 3π/4, y );
673 if( x == -∞ ) return copysign( π, y );
674 if( |y|== +∞ ) return copysign( π/2, y );
675
676 // deal with signed zeros
677 if( x == 0.0 && y != 0.0 ) return copysign( π/2, y );
678 if( x >=+0.0 && y == 0.0 ) return copysign( 0.0, y );
679 if( x <=-0.0 && y == 0.0 ) return copysign( π, y );
680
681 // calculate ATAN2 textbook style
682 if( x > 0.0 ) return ATAN( |y / x| );
683 if( x < 0.0 ) return π - ATAN( |y / x| );
684 }
685
686
687 Yet the proposed encoding makes ATAN2 the primitive and has ATAN invent
688 a constant and then call/use ATAN2.
689
690 When one considers an implementation of ATAN, one must consider several
691 ranges of evaluation::
692
693 x [ -∞, -1.0]:: ATAN( x ) = -π/2 + ATAN( 1/x );
694 x (-1.0, +1.0]:: ATAN( x ) = + ATAN( x );
695 x [ 1.0, +∞]:: ATAN( x ) = +π/2 - ATAN( 1/x );
696
697 I should point out that the add/sub of π/2 can not lose significance
698 since the result of ATAN(1/x) is bounded 0..π/2
699
700 The bottom line is that I think you are choosing to make too many of
701 these into OpCodes, making the hardware function/calculation unit (and
702 sequencer) more complicated that necessary.
703
704 --------------------------------------------------------
705
706 We therefore I think have a case for bringing back ATAN and including ATAN2.
707
708 The reason is that whilst a microcode-like GPU-centric platform would do ATAN2 in terms of ATAN, a UNIX-centric platform would do it the other way round.
709
710 (that is the hypothesis, to be evaluated for correctness. feedback requested).
711
712 This because we cannot compromise or prioritise one platfrom's
713 speed/accuracy over another. That is not reasonable or desirable, to
714 penalise one implementor over another.
715
716 Thus, all implementors, to keep interoperability, must both have both
717 opcodes and may choose, at the architectural and routing level, which
718 one to implement in terms of the other.
719
720 Allowing implementors to choose to add either opcode and let traps sort it
721 out leaves an uncertainty in the software developer's mind: they cannot
722 trust the hardware, available from many vendors, to be performant right
723 across the board.
724
725 Standards are a pig.
726
727 ---
728
729 I might suggest that if there were a way for a calculation to be performed
730 and the result of that calculation chained to a subsequent calculation
731 such that the precision of the result-becomes-operand is wider than
732 what will fit in a register, then you can dramatically reduce the count
733 of instructions in this category while retaining
734
735 acceptable accuracy:
736
737 z = x / y
738
739 can be calculated as::
740
741 z = x * (1/y)
742
743 Where 1/y has about 26-to-32 bits of fraction. No, it's not IEEE 754-2008
744 accurate, but GPUs want speed and
745
746 1/y is fully pipelined (F32) while x/y cannot be (at reasonable area). It
747 is also not "that inaccurate" displaying 0.625-to-0.52 ULP.
748
749 Given that one has the ability to carry (and process) more fraction bits,
750 one can then do high precision multiplies of π or other transcendental
751 radixes.
752
753 And GPUs have been doing this almost since the dawn of 3D.
754
755 // calculate ATAN2 high performance style
756 // Note: at this point x != y
757 //
758 if( x > 0.0 )
759 {
760 if( y < 0.0 && |y| < |x| ) return - π/2 - ATAN( x / y );
761 if( y < 0.0 && |y| > |x| ) return + ATAN( y / x );
762 if( y > 0.0 && |y| < |x| ) return + ATAN( y / x );
763 if( y > 0.0 && |y| > |x| ) return + π/2 - ATAN( x / y );
764 }
765 if( x < 0.0 )
766 {
767 if( y < 0.0 && |y| < |x| ) return + π/2 + ATAN( x / y );
768 if( y < 0.0 && |y| > |x| ) return + π - ATAN( y / x );
769 if( y > 0.0 && |y| < |x| ) return + π - ATAN( y / x );
770 if( y > 0.0 && |y| > |x| ) return +3π/2 + ATAN( x / y );
771 }
772
773 This way the adds and subtracts from the constant are not in a precision
774 precarious position.