Skip to content

riscv_fpu: static/FMA rounding modes, exact integral rounding, FMA underflow tininess #246

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
LekKit merged 10 commits into from
Jul 22, 2026
Merged

riscv_fpu: static/FMA rounding modes, exact integral rounding, FMA underflow tininess #246

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
10 commits
Select commit Hold shift + click to select a range
9c7f6e4
fpu_lib: round to integral values on the encoding, not via biasing adds
SolAstrius Jul 17, 2026
5f27427
riscv_fpu: honor the static rounding-mode field
SolAstrius Jul 17, 2026
8cf6ee1
riscv_fpu: honor the rounding mode in the FMA ops
SolAstrius Jul 17, 2026
b17bfc5
fpu_lib: set the FMA underflow flag by IEEE after-rounding tininess
SolAstrius Jul 17, 2026
1587da8
fpu_lib: clarify the integral-rounding edge cases
SolAstrius Jul 20, 2026
b3a3f3a
riscv_fpu: share static-rm handling between the OP-FP and FMA dispatches
SolAstrius Jul 20, 2026
c2379db
fpu_lib: odd-round the FMA32 underflow recompute
SolAstrius Jul 20, 2026
de1b6f6
fpu_lib: clang-format the touched declarations
SolAstrius Jul 20, 2026
bd053f0
fpu_lib: decide FMA32 tininess by an f64 threshold compare
SolAstrius Jul 20, 2026
4c1c294
riscv_fpu: keep the RMM preparation away from non-arithmetic ops
SolAstrius Jul 20, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
52 changes: 47 additions & 5 deletions src/cpu/riscv_fpu.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -81,8 +81,9 @@ static void riscv_prepare_rmm(rvvm_hart_t* vm, const uint32_t insn, const size_t
neg = fpu_signbit64(riscv_view_d(vm, rs1)) != fpu_signbit64(riscv_view_d(vm, rs2));
break;
default:
neg = fpu_signbit64(riscv_view_d(vm, rs1));
break;
// Only add/sub/mul/div need the directed synthesis: sqrt has no exact
// ties, and ops taking rm as an argument handle RMM natively
Comment thread
purplesyringa marked this conversation as resolved.
return;
}

// Round to positive/negative infinity based on the result sign
Expand All @@ -93,7 +94,33 @@ static void riscv_prepare_rmm(rvvm_hart_t* vm, const uint32_t insn, const size_t
}
}

slow_path void riscv_emulate_f_opc_op(rvvm_hart_t* vm, const uint32_t insn)
// funct3 is an rm field only on the rounding-capable OP-FP ops; on
// fsgnj/fmin/fmax/fcmp/fclass/fmv it encodes the operation itself. "Implicitly"
// rounding: ops that take rm as an argument (fcvt to integer, fround) consume
// the field themselves and are not listed here.
static forceinline bool riscv_f_op_is_implicitly_rounding(const uint32_t insn)
{
switch (insn & 0xFE000000UL) {
case 0x00000000UL: // fadd.s
case 0x02000000UL: // fadd.d
case 0x08000000UL: // fsub.s
case 0x0A000000UL: // fsub.d
case 0x10000000UL: // fmul.s
case 0x12000000UL: // fmul.d
case 0x18000000UL: // fdiv.s
case 0x1A000000UL: // fdiv.d
case 0x58000000UL: // fsqrt.s
case 0x5A000000UL: // fsqrt.d
case 0x40000000UL: // fcvt.s.d, fround.s (Zfa)
case 0x42000000UL: // fcvt.d.s, fround.d (Zfa)
case 0xD0000000UL: // fcvt.s.w[u]/l[u]
case 0xD2000000UL: // fcvt.d.w[u]/l[u]
return true;
}
return false;
}

static slow_path void riscv_emulate_f_opc_op_impl(rvvm_hart_t* vm, const uint32_t insn)
{
const size_t rds = bit_ext_u32(insn, 7, 5);
const uint32_t rm = bit_ext_u32(insn, 12, 3);
Expand All @@ -102,8 +129,9 @@ slow_path void riscv_emulate_f_opc_op(rvvm_hart_t* vm, const uint32_t insn)

if (likely(riscv_fpu_is_enabled(vm))) {

if (unlikely(vm->csr.fcsr >> 5 == 0x04)) {
// Handle RMM rounding
if (unlikely(((rm == 0x07) ? vm->csr.fcsr >> 5 : rm) == 0x04)) {
// Handle RMM rounding in the effective mode: a static rmm field
// behaves exactly like frm == RMM
riscv_prepare_rmm(vm, insn, rs1, rs2);
}

Expand Down Expand Up @@ -406,4 +434,18 @@ slow_path void riscv_emulate_f_opc_op(rvvm_hart_t* vm, const uint32_t insn)
riscv_illegal_insn(vm, insn);
}

slow_path void riscv_emulate_f_opc_op(rvvm_hart_t* vm, const uint32_t insn)
{
const uint32_t rm = bit_ext_u32(insn, 12, 3);
// A static rm field on an implicitly rounding op overrides the frm-tracked
// host mode around the op
if (unlikely(rm != 0x07) && riscv_fpu_rm_is_valid(rm) && riscv_f_op_is_implicitly_rounding(insn)) {
const uint32_t prev_rm = riscv_fpu_static_rm_enter(rm);
riscv_emulate_f_opc_op_impl(vm, insn);
riscv_fpu_static_rm_leave(prev_rm);
} else {
riscv_emulate_f_opc_op_impl(vm, insn);
}
}

#endif
51 changes: 51 additions & 0 deletions src/cpu/riscv_fpu.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,37 @@ static forceinline bool riscv_fpu_rm_is_valid(uint32_t rm)
return rm > 1;
}

// Host mode implementing a given rm: RMM has no host equivalent and runs in RNE,
// which differs only on exact halfway ties
static forceinline uint32_t riscv_fpu_host_rm(uint32_t rm)
{
return (rm == FPU_LIB_ROUND_MM) ? FPU_LIB_ROUND_NE : rm;
}

/*
* A rounding-capable op runs under its effective mode: the frm-tracked host mode
* when rm == DYN, or the static rm field applied around the op. Enter returns the
* mode to restore afterwards via leave, or DYN when nothing was changed.
*/
static forceinline uint32_t riscv_fpu_static_rm_enter(uint32_t rm)
{
if (unlikely(rm != 0x07)) {
// Always report a mode to restore: the op itself may change it further
// (the RMM preparation), and leave must undo that too
const uint32_t prev = fpu_get_rounding_mode();
fpu_set_rounding_mode(riscv_fpu_host_rm(rm));
return prev;
}
return 0x07;
}

static forceinline void riscv_fpu_static_rm_leave(uint32_t prev)
{
if (unlikely(prev != 0x07)) {
fpu_set_rounding_mode(prev);
}
}

// Bit-precise register read (raw low 32 bits, no NaN-box check) -- for fmv.x.w
static forceinline fpu_f32_t riscv_view_s(rvvm_hart_t* vm, size_t reg)
{
Expand Down Expand Up @@ -135,20 +166,25 @@ static forceinline void riscv_emulate_f_fmadd(rvvm_hart_t* vm, const uint32_t in
const size_t rs3 = insn >> 27;

if (likely(riscv_fpu_is_enabled(vm) && riscv_fpu_rm_is_valid(rm))) {
// A static rm field overrides the frm-tracked host mode, as in the OP-FP dispatch
const uint32_t prev_rm = riscv_fpu_static_rm_enter(rm);
switch (bit_ext_u32(insn, 25, 2)) {
case 0x0: // fmadd.s
riscv_emit_s(vm, rds,
fpu_fma32(riscv_read_s(vm, rs1), //
riscv_read_s(vm, rs2), //
riscv_read_s(vm, rs3)));
riscv_fpu_static_rm_leave(prev_rm);
return;
case 0x1: // fmadd.d
riscv_emit_d(vm, rds,
fpu_fma64(riscv_view_d(vm, rs1), //
riscv_view_d(vm, rs2), //
riscv_view_d(vm, rs3)));
riscv_fpu_static_rm_leave(prev_rm);
return;
}
riscv_fpu_static_rm_leave(prev_rm);
}

riscv_illegal_insn(vm, insn);
Expand All @@ -163,20 +199,25 @@ static forceinline void riscv_emulate_f_fmsub(rvvm_hart_t* vm, const uint32_t in
const size_t rs3 = insn >> 27;

if (likely(riscv_fpu_is_enabled(vm) && riscv_fpu_rm_is_valid(rm))) {
// A static rm field overrides the frm-tracked host mode, as in the OP-FP dispatch
const uint32_t prev_rm = riscv_fpu_static_rm_enter(rm);
switch (bit_ext_u32(insn, 25, 2)) {
case 0x0: // fmsub.s
riscv_emit_s(vm, rds,
fpu_fma32(riscv_read_s(vm, rs1), //
riscv_read_s(vm, rs2), //
fpu_neg32(riscv_read_s(vm, rs3))));
riscv_fpu_static_rm_leave(prev_rm);
return;
case 0x1: // fmsub.d
riscv_emit_d(vm, rds,
fpu_fma64(riscv_view_d(vm, rs1), //
riscv_view_d(vm, rs2), //
fpu_neg64(riscv_view_d(vm, rs3))));
riscv_fpu_static_rm_leave(prev_rm);
return;
}
riscv_fpu_static_rm_leave(prev_rm);
}

riscv_illegal_insn(vm, insn);
Expand All @@ -191,20 +232,25 @@ static forceinline void riscv_emulate_f_fnmsub(rvvm_hart_t* vm, const uint32_t i
const size_t rs3 = insn >> 27;

if (likely(riscv_fpu_is_enabled(vm) && riscv_fpu_rm_is_valid(rm))) {
// A static rm field overrides the frm-tracked host mode, as in the OP-FP dispatch
const uint32_t prev_rm = riscv_fpu_static_rm_enter(rm);
switch (bit_ext_u32(insn, 25, 2)) {
case 0x0: // fnmsub.s
riscv_emit_s(vm, rds,
fpu_fma32(fpu_neg32(riscv_read_s(vm, rs1)), //
riscv_read_s(vm, rs2), //
riscv_read_s(vm, rs3)));
riscv_fpu_static_rm_leave(prev_rm);
return;
case 0x1: // fnmsub.d
riscv_emit_d(vm, rds,
fpu_fma64(fpu_neg64(riscv_view_d(vm, rs1)), //
riscv_view_d(vm, rs2), //
riscv_view_d(vm, rs3)));
riscv_fpu_static_rm_leave(prev_rm);
return;
}
riscv_fpu_static_rm_leave(prev_rm);
}

riscv_illegal_insn(vm, insn);
Expand All @@ -219,21 +265,26 @@ static forceinline void riscv_emulate_f_fnmadd(rvvm_hart_t* vm, const uint32_t i
const size_t rs3 = insn >> 27;

if (likely(riscv_fpu_is_enabled(vm) && riscv_fpu_rm_is_valid(rm))) {
// A static rm field overrides the frm-tracked host mode, as in the OP-FP dispatch
const uint32_t prev_rm = riscv_fpu_static_rm_enter(rm);
switch (bit_ext_u32(insn, 25, 2)) {
case 0x0: // fnmadd.s = -(rs1*rs2) - rs3; negate operands so the single
// rounding sees the correctly-signed result (directed modes)
riscv_emit_s(vm, rds,
fpu_fma32(fpu_neg32(riscv_read_s(vm, rs1)), //
riscv_read_s(vm, rs2), //
fpu_neg32(riscv_read_s(vm, rs3))));
riscv_fpu_static_rm_leave(prev_rm);
return;
case 0x1: // fnmadd.d
riscv_emit_d(vm, rds,
fpu_fma64(fpu_neg64(riscv_view_d(vm, rs1)), //
riscv_view_d(vm, rs2), //
fpu_neg64(riscv_view_d(vm, rs3))));
riscv_fpu_static_rm_leave(prev_rm);
return;
}
riscv_fpu_static_rm_leave(prev_rm);
}

riscv_illegal_insn(vm, insn);
Expand Down
108 changes: 86 additions & 22 deletions src/util/fpu_lib.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -405,54 +405,118 @@ slow_path uint32_t fpu_fclass64(fpu_f64_t d)
return ret;
}

/*
* Round to an integral value directly on the encoding: no host FP arithmetic, so
* no spurious exception flags and no double rounding. Only DYN falls back to the
* tracked mode; an explicit RMM request is honored as-is.
*/
slow_path fpu_f32_t fpu_round_f32_internal(fpu_f32_t f, uint32_t mode)
{
uint32_t u = fpu_bit_f32_to_u32(f);
uint32_t s = u & FPU_LIB_FP32_SIGNEDFP_MASK;
if (unlikely(mode > FPU_LIB_ROUND_UP)) {
const uint32_t u = fpu_bit_f32_to_u32(f);
const uint32_t s = u & FPU_LIB_FP32_SIGNEDFP_MASK;
const int32_t e = fpu_exponent32(f);
bool away = false;
if (unlikely(mode > FPU_LIB_ROUND_MM)) {
mode = fpu_get_rounding_mode();
}
if (e >= 23 || !(u << 1)) {
return f; // Already integral: |f| >= 2^23, +/-0, inf, NaN
}
if (e < 0) {
// |f| in (0, 1) rounds to +/-0 or +/-1
switch (mode) {
case FPU_LIB_ROUND_NE:
// |f| in (0.5, 1) is nearest to 1; exactly 0.5 is a tie and goes
// to the even neighbour, which is 0
away = (e == -1) && (u & FPU_LIB_FP32_MANTISSA_MASK);
Comment thread
purplesyringa marked this conversation as resolved.
break;
case FPU_LIB_ROUND_MM:
away = (e == -1);
break;
case FPU_LIB_ROUND_DN:
away = !!s;
break;
case FPU_LIB_ROUND_UP:
away = !s;
break;
}
return fpu_bit_u32_to_f32(s | (away ? 0x3F800000U : 0)); // +/-1.0 or +/-0.0
}
// |f| in [1, 2^23): split the mantissa into integer part and fraction bits
const uint32_t frac = u & (FPU_LIB_FP32_MANTISSA_MASK >> e);
const uint32_t half = (1U << 22) >> e;
const uint32_t step = (1U << 23) >> e; // 1.0 at this exponent
switch (mode) {
case FPU_LIB_ROUND_NE:
// On a tie, round up iff the integer part is odd. For e == 0 the step
// bit is the exponent LSB rather than a mantissa bit, but the biased
// exponent of [1, 2) is odd (0x7F/0x3FF), matching its odd integer 1.
away = frac > half || (frac == half && (u & step));
Comment thread
purplesyringa marked this conversation as resolved.
break;
case FPU_LIB_ROUND_MM:
return fpu_add32(f, fpu_bit_u32_to_f32(0x3F000000U | s));
away = frac >= half;
break;
case FPU_LIB_ROUND_DN:
if (s && fpu_is_fractional32(f)) {
return fpu_sub32(f, fpu_bit_u32_to_f32(0x3F800000U));
}
away = s && frac;
break;
case FPU_LIB_ROUND_UP:
if (!s && fpu_is_fractional32(f)) {
return fpu_add32(f, fpu_bit_u32_to_f32(0x3F800000U));
}
away = !s && frac;
break;
}
return f;
// Adding step may carry from the mantissa into the exponent field: that only
// happens when the truncated mantissa wraps to zero, i.e. when rounding away
// lands exactly on the next power of two, where the carry is the intended
// encoding (the same trick as the classic nextafter bit-increment).
return fpu_bit_u32_to_f32((u - frac) + (away ? step : 0));
}

slow_path fpu_f64_t fpu_round_f64_internal(fpu_f64_t d, uint32_t mode)
{
uint64_t u = fpu_bit_f64_to_u64(d);
uint64_t s = u & FPU_LIB_FP64_SIGNEDFP_MASK;
if (unlikely(mode > FPU_LIB_ROUND_UP)) {
const uint64_t u = fpu_bit_f64_to_u64(d);
const uint64_t s = u & FPU_LIB_FP64_SIGNEDFP_MASK;
const int32_t e = fpu_exponent64(d);
bool away = false;
if (unlikely(mode > FPU_LIB_ROUND_MM)) {
mode = fpu_get_rounding_mode();
}
if (e >= 52 || !(u << 1)) {
return d; // Already integral: |d| >= 2^52, +/-0, inf, NaN
}
if (e < 0) {
switch (mode) {
case FPU_LIB_ROUND_NE:
away = (e == -1) && (u & FPU_LIB_FP64_MANTISSA_MASK);
break;
case FPU_LIB_ROUND_MM:
away = (e == -1);
break;
case FPU_LIB_ROUND_DN:
away = !!s;
break;
case FPU_LIB_ROUND_UP:
away = !s;
break;
}
return fpu_bit_u64_to_f64(s | (away ? 0x3FF0000000000000ULL : 0)); // +/-1.0 or +/-0.0
}
const uint64_t frac = u & (FPU_LIB_FP64_MANTISSA_MASK >> e);
const uint64_t half = (1ULL << 51) >> e;
const uint64_t step = (1ULL << 52) >> e; // 1.0 at this exponent
switch (mode) {
case FPU_LIB_ROUND_NE:
away = frac > half || (frac == half && (u & step));
break;
case FPU_LIB_ROUND_MM:
return fpu_add64(d, fpu_bit_u64_to_f64(0x3FE0000000000000ULL | s));
away = frac >= half;
break;
case FPU_LIB_ROUND_DN:
if (s && fpu_is_fractional64(d)) {
return fpu_sub64(d, fpu_bit_u64_to_f64(0x3FF0000000000000ULL));
}
away = s && frac;
break;
case FPU_LIB_ROUND_UP:
if (!s && fpu_is_fractional64(d)) {
return fpu_add64(d, fpu_bit_u64_to_f64(0x3FF0000000000000ULL));
}
away = !s && frac;
break;
}
return d;
return fpu_bit_u64_to_f64((u - frac) + (away ? step : 0));
Comment thread
purplesyringa marked this conversation as resolved.
}

#if defined(USE_SOFT_FPU_SQRT)
Expand Down
Loading
Loading