lcompilers
diff --git a/integration_tests/test_random.py b/integration_tests/test_random.py
@@ -30,14 +30,41 @@ def test_randint():
     print(ri2)
     assert ri2 >= -50 and ri2 <= 76
 
+def test_uniform():
+    r: f64
+    r = random.uniform(5., 76.)
+    print(r)
+    r = random.uniform(-50., 76.)
+    print(r)
+
 def test_paretovariate():
     r: f64
     r = random.paretovariate(2.0)
     print(r)
     r = random.paretovariate(-5.6)
     print(r)
 
-test_random()
-test_randrange()
-test_randint()
-test_paretovariate()
+def test_expovariate():
+    r: f64
+    r = random.expovariate(2.0)
+    print(r)
+    r = random.expovariate(-5.6)
+    print(r)
+
+def test_weibullvariate():
+    r: f64
+    r = random.weibullvariate(2.0, 3.0)
+    print(r)
+    r = random.weibullvariate(-5.6, 1.2)
+    print(r)
+
+def check():
+    test_random()
+    test_randrange()
+    test_randint()
+    test_uniform()
+    test_paretovariate()
+    test_expovariate()
+    test_weibullvariate()
+
+check()
diff --git a/src/runtime/impure/lfortran_intrinsics.c b/src/runtime/impure/lfortran_intrinsics.c
@@ -30,9 +30,9 @@ LFORTRAN_API void _lfortran_random_number(int n, double *v)
     }
 }
 
-LFORTRAN_API float _lfortran_random_float()
+LFORTRAN_API double _lfortran_random()
 {
-    return ((float) rand() / (float) RAND_MAX);
+    return (rand() / (double) RAND_MAX);
 }
 
 LFORTRAN_API int _lfortran_randrange(int lower, int upper)

diff --git a/src/runtime/impure/lfortran_intrinsics.h b/src/runtime/impure/lfortran_intrinsics.h
@@ -31,7 +31,7 @@ typedef double _Complex double_complex_t;
 
 LFORTRAN_API double _lfortran_sum(int n, double *v);
 LFORTRAN_API void _lfortran_random_number(int n, double *v);
-LFORTRAN_API float _lfortran_random_float();
+LFORTRAN_API double _lfortran_random();
 LFORTRAN_API int _lfortran_randrange(int lower, int upper);
 LFORTRAN_API int _lfortran_random_int(int lower, int upper);
 LFORTRAN_API void _lfortran_printf(const char* format, ...);

diff --git a/src/runtime/random.py b/src/runtime/random.py
@@ -1,14 +1,35 @@
 from ltypes import i32, f64, ccall
 
+e: f64 = 2.718281828459045235360287471352662497757
+eps: f64 = 1e-16
+
+#: TODO: Call `log` from C directly until we fix the multiple import issue
+def _log(x: f64) -> f64:
+    return _lfortran_dlog(x)
+
+@ccall
+def _lfortran_dlog(x: f64) -> f64:
+    pass
+
+def _exp(x: f64) -> f64:
+    return e**x
+
+def _sqrt(x: f64) -> f64:
+    return x**(1/2)
+
+def _abs(x: f64) -> f64:
+    if x < 0.0:
+        return -x
+    return x
 
 def random() -> f64:
     """
     Returns a random floating point number in the range [0.0, 1.0)
     """
-    return _lfortran_random_float()
+    return _lfortran_random()
 
 @ccall
-def _lfortran_random_float() -> f64:
+def _lfortran_random() -> f64:
     pass
 
 def randrange(lower: i32, upper: i32) -> i32:
@@ -31,10 +52,32 @@ def randint(lower: i32, upper: i32) -> i32:
 def _lfortran_random_int(lower: i32, upper: i32) -> i32:
     pass
 
+def uniform(a: f64, b: f64) -> f64:
+    """
+    Get a random number in the range [a, b) or [a, b] depending on rounding.
+    """
+    return a + (b - a) * random()
+
 def paretovariate(alpha: f64) -> f64:
     """
     Return a random number from a Pareto distribution with parameter `alpha`.
     """
     u: f64
     u = 1.0 - random()
     return u ** (-1.0 / alpha)
+
+def expovariate(l: f64) -> f64:
+    """
+    Return a random number from an exponential distribution with parameter
+    `l` (lambda).
+    """
+    assert _abs(l) > eps
+    return -_log(1.0 - random()) / l
+
+def weibullvariate(alpha: f64, beta: f64) -> f64:
+    """
+    Return a random number from a Weibull distribution with parameters `alpha`
+    and `beta`.
+    """
+    assert _abs(beta) > eps
+    return alpha * (-_log(1.0 - random())) ** (1.0 / beta)