fpga: sim: Add PkgComplex, PkgMath, and PkgRandom

PkgComplex adds functions for doing complex arithmetic in SystemVerilog
simulation.

PkgMath provides mathematical operations and constants that aren't
built into SystemVerilog, such as a constant for pi and the function
round().

PkgRandom adds randomization functions beyond what standard Verilog
supports but that don't require any special licenses or simulators.

1 files changed, 230 insertions, 0 deletions

diff --git a/fpga/usrp3/sim/packages/PkgComplex.sv b/fpga/usrp3/sim/packages/PkgComplex.sv
new file mode 100644
index 000000000..9fc957b97
--- /dev/null
+++ b/fpga/usrp3/sim/packages/PkgComplex.sv
@@ -0,0 +1,230 @@
+//
+// Copyright 2021 Ettus Research, A National Instruments Brand
+//
+// SPDX-License-Identifier: LGPL-3.0-or-later
+//
+// Package: PkgComplex
+//
+// Description:
+//
+//   A package for doing complex arithmetic in SystemVerilog simulations.
+//   Fixed-point operations are implemented such that results clip to the range
+//   [-1.0, 1.0) and are rounded to the nearest ULP (half ULP is rounded away
+//   from zero, following SystemVerilog rounding behavior).
+//
+
+
+package PkgComplex;
+
+  //---------------------------------------------------------------------------
+  // Type Definitions
+  //---------------------------------------------------------------------------
+
+  // Define a signed 16-bit fixed-point type with 15 fractional bits (Q0.15).
+  typedef bit signed [15:0] s16_t;
+
+  // Signed complex 16-bit data type, the standard type used by UHD and RFNoC.
+  typedef struct packed {
+    s16_t re;
+    s16_t im;
+  } sc16_t;
+
+  // Complex floating point data type.
+  typedef struct {
+    real re;
+    real im;
+  } complex_t;
+
+  // Maximum and minimum allowed by the s16 type.
+  localparam s16_t MAX_S16 = 16'h7FFF;
+  localparam s16_t MIN_S16 = 16'h8000;
+
+
+  //---------------------------------------------------------------------------
+  // Conversion Functions
+  //---------------------------------------------------------------------------
+
+  // Create an sc16 value from two s16 values.
+  function sc16_t build_sc16(s16_t x = 0, s16_t y = 0);
+    sc16_t val;
+    val.re = x;
+    val.im = y;
+    return val;
+  endfunction : build_sc16
+
+  // Create a complex value from two real values.
+  function complex_t build_complex(real x = 0.0, real y = 0.0);
+    complex_t val;
+    val.re = x;
+    val.im = y;
+    return val;
+  endfunction : build_complex
+
+  // Convert s16 to real.
+  function real s16_to_real(s16_t x);
+    return real'(x) / (2.0**15);
+  endfunction : s16_to_real
+
+  // Convert real to s16.
+  function s16_t real_to_s16(real x);
+    real val;
+    val = x * (2.0**15);
+    val = (val > MAX_S16) ? MAX_S16 : val;
+    val = (val < MIN_S16) ? MIN_S16 : val;
+    return s16_t'(val);
+  endfunction : real_to_s16
+
+  // Convert complex to sc16.
+  function sc16_t complex_to_sc16(complex_t x);
+    sc16_t val;
+    val.re = real_to_s16(x.re);
+    val.im = real_to_s16(x.im);
+    return val;
+  endfunction : complex_to_sc16
+
+  // Convert sc16 to complex.
+  function complex_t sc16_to_complex(sc16_t x);
+    complex_t val;
+    val.re = s16_to_real(x.re);
+    val.im = s16_to_real(x.im);
+    return val;
+  endfunction : sc16_to_complex
+
+  // Convert polar coordinates to a complex number. The phase should be in
+  // radians.
+  function complex_t polar_to_complex(real mag, real phase);
+    complex_t val;
+    val.re = mag * $cos(phase);
+    val.im = mag * $sin(phase);
+    return val;
+  endfunction : polar_to_complex
+
+  // Convert polar coordinates to an sc16 complex number. The phase should be
+  // in radians.
+  function sc16_t polar_to_sc16(real mag, real phase);
+    return complex_to_sc16(polar_to_complex(mag, phase));
+  endfunction : polar_to_sc16
+
+
+  //---------------------------------------------------------------------------
+  // Floating Point Complex Arithmetic
+  //---------------------------------------------------------------------------
+
+  // Add two complex numbers: x + y
+  function complex_t add(complex_t x, complex_t y);
+    complex_t val;
+    val.re = x.re + y.re;
+    val.im = x.im + y.im;
+    return val;
+  endfunction : add
+
+  // Subtract two complex numbers: x - y
+  function complex_t sub(complex_t x, complex_t y);
+    complex_t val;
+    val.re = x.re - y.re;
+    val.im = x.im - y.im;
+    return val;
+  endfunction : sub
+
+  // Multiply two complex numbers: x * y
+  function complex_t mul(complex_t x, complex_t y);
+    complex_t val;
+    val.re = x.re*y.re - x.im*y.im;
+    val.im = x.re*y.im + x.im*y.re;
+    return val;
+  endfunction : mul
+
+  // Divide two complex numbers: x / y
+  function complex_t div(complex_t x, complex_t y);
+    complex_t z;
+    z.re = (x.re*y.re + x.im*y.im) / (y.re*y.re + y.im*y.im);
+    z.im = (x.im*y.re + x.re*y.im) / (y.re*y.re + y.im*y.im);
+    return z;
+  endfunction : div
+
+  // Compute the exponential: e^x
+  function complex_t exp(complex_t x);
+    complex_t val;
+    // exp(a+jb) = exp(a)*exp(jb) = exp(a)*(cos(b) + j*sin(b))
+    val.re = $exp(x.re)*$cos(x.im);
+    val.im = $exp(x.re)*$sin(x.im);
+    return val;
+  endfunction : exp
+
+  // Compute the sine: sin(x)
+  function complex_t sin(complex_t x);
+    complex_t val;
+    val.re = $sin(x.re)*$cosh(x.im);
+    val.im = $cos(x.re)*$sinh(x.im);
+    return val;
+  endfunction : sin
+
+  // Compute the cosine: cos(x)
+  function complex_t cos(complex_t x);
+    complex_t val;
+    val.re =      $cos(x.re)*$cosh(x.im);
+    val.im = -1.0*$sin(x.re)*$sinh(x.im);
+    return val;
+  endfunction : cos
+
+  // Compute the magnitude/modulus/absolute value: |x|
+  function real mag(complex_t x);
+    return $sqrt(x.re*x.re + x.im*x.im);
+  endfunction : mag
+
+  // Compute the phase/argument: arg(x)
+  function real arg(complex_t x);
+    return $atan2(x.im, x.re);
+  endfunction : arg
+
+
+  //---------------------------------------------------------------------------
+  // Fixed-Point Complex Arithmetic
+  //---------------------------------------------------------------------------
+
+  // Add two complex numbers: x + y
+  function sc16_t add_sc16(sc16_t x, sc16_t y);
+    return complex_to_sc16(add(sc16_to_complex(x), sc16_to_complex(y)));
+  endfunction : add_sc16
+
+  // Subtract two complex numbers: x - y
+  function sc16_t sub_sc16(sc16_t x, sc16_t y);
+    return complex_to_sc16(sub(sc16_to_complex(x), sc16_to_complex(y)));
+  endfunction : sub_sc16
+
+  // Multiply two complex numbers: x * y
+  function sc16_t mul_sc16(sc16_t x, sc16_t y);
+    return complex_to_sc16(mul(sc16_to_complex(x), sc16_to_complex(y)));
+  endfunction : mul_sc16
+
+  // Divide two complex numbers: x / y
+  function sc16_t div_sc16(sc16_t x, sc16_t y);
+    return complex_to_sc16(div(sc16_to_complex(x), sc16_to_complex(y)));
+  endfunction : div_sc16
+
+  // Compute the exponential: e^x
+  function sc16_t exp_sc16(sc16_t x);
+    return complex_to_sc16(exp(sc16_to_complex(x)));
+  endfunction : exp_sc16
+
+  // Compute the sine: sin(x)
+  function sc16_t sin_sc16(sc16_t x);
+    return complex_to_sc16(sin(sc16_to_complex(x)));
+  endfunction : sin_sc16
+
+  // compute the cosine: cos(x)
+  function sc16_t cos_sc16(sc16_t x);
+    return complex_to_sc16(cos(sc16_to_complex(x)));
+  endfunction : cos_sc16
+
+  // Compute the magnitude/modulus/absolute value: |x|
+  function s16_t mag_sc16(sc16_t x);
+    return real_to_s16(mag(sc16_to_complex(x)));
+  endfunction : mag_sc16
+
+  // Compute the phase/argument: arg(x)
+  function s16_t arg_sc16(sc16_t x);
+    return real_to_s16(arg(sc16_to_complex(x)));
+  endfunction : arg_sc16
+
+endpackage : PkgComplex