From f668b9cb68e92ddf46fe471ca2f08cf18b4bbe50 Mon Sep 17 00:00:00 2001 From: william-dawson Date: Mon, 3 Feb 2020 12:55:42 +0900 Subject: [PATCH] Memory leak fix in the dense block branch. (#136) * Memory leak fix in the dense block branch. * Another missed destructor in the Dense Branch * Another fix. * Some further cleanups to avoid memory leaks. This also removes the need to use the preprocessor on this module (can use built-in include instead). * Optimization. * Remove need for explicit transpose in dense branch * Bug fix and better testing. --- Source/Fortran/DMatrixModule.F90 | 165 ++++++++++++++---- Source/Fortran/SMatrixAlgebraModule.F90 | 7 +- .../Fortran/dense_includes/ComposeMatrix.f90 | 1 + .../dense_includes/ConstructMatrixSFromD.f90 | 4 +- .../Fortran/sparse_includes/DenseBranch.f90 | 19 +- UnitTests/test_matrix.py | 84 +++++++++ 6 files changed, 221 insertions(+), 59 deletions(-) diff --git a/Source/Fortran/DMatrixModule.F90 b/Source/Fortran/DMatrixModule.F90 index b7d8efa5..40dfebf3 100644 --- a/Source/Fortran/DMatrixModule.F90 +++ b/Source/Fortran/DMatrixModule.F90 @@ -5,9 +5,10 @@ !! performance. MODULE DMatrixModule USE DataTypesModule, ONLY : NTREAL, NTCOMPLEX - USE SMatrixModule, ONLY : Matrix_lsr, Matrix_lsc + USE SMatrixModule, ONLY : Matrix_lsr, Matrix_lsc, & + & ConstructMatrixFromTripletList USE TripletListModule, ONLY : TripletList_r, TripletList_c, & - & AppendToTripletList, ConstructTripletList + & AppendToTripletList, ConstructTripletList, DestructTripletList USE TripletModule, ONLY : Triplet_r, Triplet_c IMPLICIT NONE PRIVATE @@ -101,6 +102,7 @@ PURE SUBROUTINE ConstructEmptyMatrixSup_ldr(this, rows, columns) !> Columns of the matrix INTEGER, INTENT(IN) :: columns + CALL DestructMatrix(this) this = ConstructEmptyMatrix_ldr(rows, columns) END SUBROUTINE ConstructEmptyMatrixSup_ldr !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -126,7 +128,7 @@ PURE SUBROUTINE ConstructMatrixDFromS_ldr(sparse_matrix, dense_matrix) !! Helper Variables TYPE(Triplet_r) :: temporary -#include "dense_includes/ConstructMatrixDFromS.f90" + INCLUDE "dense_includes/ConstructMatrixDFromS.f90" END SUBROUTINE ConstructMatrixDFromS_ldr !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -143,9 +145,7 @@ PURE SUBROUTINE ConstructMatrixSFromD_ldr(dense_matrix, sparse_matrix, & TYPE(Triplet_r) :: temporary TYPE(TripletList_r) :: temporary_list -#define SMTYPE Matrix_lsr -#include "dense_includes/ConstructMatrixSFromD.f90" -#undef SMTYPE + INCLUDE "dense_includes/ConstructMatrixSFromD.f90" END SUBROUTINE ConstructMatrixSFromD_ldr !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -156,7 +156,7 @@ PURE SUBROUTINE CopyMatrix_ldr(matA, matB) !> matB = matA TYPE(Matrix_ldr), INTENT(INOUT) :: matB -#include "dense_includes/CopyMatrix.f90" + INCLUDE "dense_includes/CopyMatrix.f90" END SUBROUTINE CopyMatrix_ldr !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -209,7 +209,7 @@ PURE SUBROUTINE TransposeMatrix_ldr(matA, matAT) !> matAT = matA^T. TYPE(Matrix_ldr), INTENT(INOUT) :: matAT -#include "dense_includes/TransposeMatrix.f90" + INCLUDE "dense_includes/TransposeMatrix.f90" END SUBROUTINE TransposeMatrix_ldr !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -227,7 +227,7 @@ PURE SUBROUTINE ComposeMatrix_ldr(mat_array, block_rows, block_columns, & !> The composed matrix. TYPE(Matrix_ldr), INTENT(INOUT) :: out_matrix -#include "dense_includes/ComposeMatrix.f90" + INCLUDE "dense_includes/ComposeMatrix.f90" END SUBROUTINE ComposeMatrix_ldr !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -247,21 +247,26 @@ PURE SUBROUTINE SplitMatrix_ldr(this, block_rows, block_columns, & !> Specifies the size of the columns. INTEGER, DIMENSION(:), INTENT(IN), OPTIONAL :: block_size_column_in -#include "dense_includes/SplitMatrix.f90" + INCLUDE "dense_includes/SplitMatrix.f90" END SUBROUTINE SplitMatrix_ldr !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> A wrapper for multiplying two dense matrices. - SUBROUTINE MultiplyMatrix_ldr(MatA,MatB,MatC) + SUBROUTINE MultiplyMatrix_ldr(MatA, MatB, MatC, IsATransposed_in, & + & IsBTransposed_in) !> The first matrix. TYPE(Matrix_ldr), INTENT(IN) :: MatA !> The second matrix. TYPE(Matrix_ldr), INTENT(IN) :: MatB !> MatC = MatA*MatB. TYPE(Matrix_ldr), INTENT(INOUT) :: MatC + !> True if A is already transposed. + LOGICAL, OPTIONAL, INTENT(IN) :: IsATransposed_in + !> True if B is already transposed. + LOGICAL, OPTIONAL, INTENT(IN) :: IsBTransposed_in !! Local variables - CHARACTER, PARAMETER :: TRANSA = 'N' - CHARACTER, PARAMETER :: TRANSB = 'N' + CHARACTER :: TRANSA + CHARACTER :: TRANSB INTEGER :: M INTEGER :: N INTEGER :: K @@ -271,16 +276,55 @@ SUBROUTINE MultiplyMatrix_ldr(MatA,MatB,MatC) DOUBLE PRECISION, PARAMETER :: BETA = 0.0 INTEGER :: LDC - MatC = Matrix_ldr(MatA%rows,MatB%columns) + !! Optional Parameters + TRANSA = 'N' + IF (PRESENT(IsATransposed_in)) THEN + IF (IsATransposed_in) THEN + TRANSA = 'T' + END IF + END IF + TRANSB = 'N' + IF (PRESENT(IsBTransposed_in)) THEN + IF (IsBTransposed_in) THEN + TRANSB = 'T' + END IF + END IF !! Setup Lapack - M = MatA%rows - N = MatB%columns - K = MatA%columns - LDA = M - LDB = K + IF (TRANSA .EQ. 'T') THEN + M = MatA%columns + ELSE + M = MatA%rows + END IF + + IF (TRANSB .EQ. 'T') THEN + N = MatB%rows + ELSE + N = MatB%columns + END IF + + IF (TRANSA .EQ. 'T') THEN + K = MatA%rows + ELSE + K = MatA%columns + END IF + + IF (TRANSA .EQ. 'T') THEN + LDA = K + ELSE + LDA = M + END IF + + IF (TRANSB .EQ. 'T') THEN + LDB = N + ELSE + LDB = K + END IF + LDC = M + !! Multiply + CALL ConstructEmptyMatrix(MatC, M, N) CALL DGEMM(TRANSA, TRANSB, M, N, K, ALPHA, MatA%data, LDA, MatB%data, & & LDB, BETA, MatC%data, LDC) @@ -295,6 +339,7 @@ PURE SUBROUTINE ConstructEmptyMatrixSup_ldc(this, rows, columns) !> The number of columns o the matrix. INTEGER, INTENT(IN) :: columns + CALL DestructMatrix(this) this = ConstructEmptyMatrix_ldc(rows, columns) END SUBROUTINE ConstructEmptyMatrixSup_ldc !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -320,7 +365,7 @@ PURE SUBROUTINE ConstructMatrixDFromS_ldc(sparse_matrix, dense_matrix) !! Helper Variables TYPE(Triplet_c) :: temporary -#include "dense_includes/ConstructMatrixDFromS.f90" + INCLUDE "dense_includes/ConstructMatrixDFromS.f90" END SUBROUTINE ConstructMatrixDFromS_ldc !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -337,9 +382,7 @@ PURE SUBROUTINE ConstructMatrixSFromD_ldc(dense_matrix, sparse_matrix, & TYPE(Triplet_c) :: temporary TYPE(TripletList_c) :: temporary_list -#define SMTYPE Matrix_lsc -#include "dense_includes/ConstructMatrixSFromD.f90" -#undef SMTYPE + INCLUDE "dense_includes/ConstructMatrixSFromD.f90" END SUBROUTINE ConstructMatrixSFromD_ldc !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -350,7 +393,7 @@ PURE SUBROUTINE CopyMatrix_ldc(matA, matB) !> matB = matA TYPE(Matrix_ldc), INTENT(INOUT) :: matB -#include "dense_includes/CopyMatrix.f90" + INCLUDE "dense_includes/CopyMatrix.f90" END SUBROUTINE CopyMatrix_ldc !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -387,12 +430,14 @@ FUNCTION MatrixNorm_ldc(this) RESULT(norm) REAL(NTREAL) :: norm !! Local Variables INTEGER :: II, JJ + COMPLEX(NTCOMPLEX) :: val, conjval norm = 0 DO II =1, this%rows DO JJ = 1, this%columns - norm = norm + & - & REAL(this%data(II,JJ)*CONJG(this%data(II,JJ)),KIND=NTREAL) + val = this%data(II,JJ) + conjval = CONJG(val) + norm = norm + REAL(val*conjval,KIND=NTREAL) END DO END DO END FUNCTION MatrixNorm_ldc @@ -404,7 +449,7 @@ PURE SUBROUTINE TransposeMatrix_ldc(matA, matAT) !> matAT = matA^T. TYPE(Matrix_ldc), INTENT(INOUT) :: matAT -#include "dense_includes/TransposeMatrix.f90" + INCLUDE "dense_includes/TransposeMatrix.f90" END SUBROUTINE TransposeMatrix_ldc !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -422,7 +467,7 @@ PURE SUBROUTINE ComposeMatrix_ldc(mat_array, block_rows, block_columns, & !> The composed matrix. TYPE(Matrix_ldc), INTENT(INOUT) :: out_matrix -#include "dense_includes/ComposeMatrix.f90" + INCLUDE "dense_includes/ComposeMatrix.f90" END SUBROUTINE ComposeMatrix_ldc !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -442,21 +487,26 @@ PURE SUBROUTINE SplitMatrix_ldc(this, block_rows, block_columns, & !> Specifies the size of the columns. INTEGER, DIMENSION(:), INTENT(IN), OPTIONAL :: block_size_column_in -#include "dense_includes/SplitMatrix.f90" + INCLUDE "dense_includes/SplitMatrix.f90" END SUBROUTINE SplitMatrix_ldc !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !> A wrapper for multiplying two dense matrices. - SUBROUTINE MultiplyMatrix_ldc(MatA,MatB,MatC) + SUBROUTINE MultiplyMatrix_ldc(MatA, MatB, MatC, IsATransposed_in, & + & IsBTransposed_in) !> The first matrix. TYPE(Matrix_ldc), INTENT(IN) :: MatA !> The second matrix. TYPE(Matrix_ldc), INTENT(IN) :: MatB !> MatC = MatA*MatB. TYPE(Matrix_ldc), INTENT(INOUT) :: MatC + !> True if A is already transposed. + LOGICAL, OPTIONAL, INTENT(IN) :: IsATransposed_in + !> True if B is already transposed. + LOGICAL, OPTIONAL, INTENT(IN) :: IsBTransposed_in !! Local variables - CHARACTER, PARAMETER :: TRANSA = 'N' - CHARACTER, PARAMETER :: TRANSB = 'N' + CHARACTER :: TRANSA + CHARACTER :: TRANSB INTEGER :: M INTEGER :: N INTEGER :: K @@ -466,16 +516,55 @@ SUBROUTINE MultiplyMatrix_ldc(MatA,MatB,MatC) COMPLEX*16, PARAMETER :: BETA = 0.0 INTEGER :: LDC - MatC = Matrix_ldc(MatA%rows,MatB%columns) + !! Optional Parameters + TRANSA = 'N' + IF (PRESENT(IsATransposed_in)) THEN + IF (IsATransposed_in) THEN + TRANSA = 'T' + END IF + END IF + TRANSB = 'N' + IF (PRESENT(IsBTransposed_in)) THEN + IF (IsBTransposed_in) THEN + TRANSB = 'T' + END IF + END IF !! Setup Lapack - M = MatA%rows - N = MatB%columns - K = MatA%columns - LDA = M - LDB = K + IF (TRANSA .EQ. 'T') THEN + M = MatA%columns + ELSE + M = MatA%rows + END IF + + IF (TRANSB .EQ. 'T') THEN + N = MatB%rows + ELSE + N = MatB%columns + END IF + + IF (TRANSA .EQ. 'T') THEN + K = MatA%rows + ELSE + K = MatA%columns + END IF + + IF (TRANSA .EQ. 'T') THEN + LDA = K + ELSE + LDA = M + END IF + + IF (TRANSB .EQ. 'T') THEN + LDB = N + ELSE + LDB = K + END IF + LDC = M + !! Multiply + CALL ConstructEmptyMatrix(MatC, M, N) CALL ZGEMM(TRANSA, TRANSB, M, N, K, ALPHA, MatA%data, LDA, MatB%data, & & LDB, BETA, MatC%data, LDC) diff --git a/Source/Fortran/SMatrixAlgebraModule.F90 b/Source/Fortran/SMatrixAlgebraModule.F90 index d8f7c498..b5f0f929 100644 --- a/Source/Fortran/SMatrixAlgebraModule.F90 +++ b/Source/Fortran/SMatrixAlgebraModule.F90 @@ -3,7 +3,8 @@ MODULE SMatrixAlgebraModule USE DataTypesModule, ONLY : NTREAL, NTCOMPLEX USE DMatrixModule, ONLY : Matrix_ldr, Matrix_ldc, ConstructMatrixDFromS, & - & ConstructMatrixSFromD, MultiplyMatrix, DestructMatrix + & ConstructMatrixSFromD, CopyMatrix, MultiplyMatrix, TransposeMatrix, & + & DestructMatrix USE MatrixMemoryPoolModule, ONLY : MatrixMemoryPool_lr, MatrixMemoryPool_lc, & & DestructMatrixMemoryPool, CheckMemoryPoolValidity, SetPoolSparsity, & & ConstructMatrixMemoryPool @@ -420,8 +421,6 @@ SUBROUTINE DenseBranch_lsr(matA, matB, matC, IsATransposed, IsBTransposed, & !> Threshold for flushing values. REAL(NTREAL), INTENT(IN) :: threshold !! Local Data - TYPE(Matrix_lsr) :: untransposedMatA - TYPE(Matrix_lsr) :: untransposedMatB TYPE(Matrix_ldr) :: DenseA TYPE(Matrix_ldr) :: DenseB TYPE(Matrix_ldr) :: DenseC @@ -447,8 +446,6 @@ SUBROUTINE DenseBranch_lsc(matA, matB, matC, IsATransposed, IsBTransposed, & !> Threshold for flushing values. REAL(NTREAL), INTENT(IN) :: threshold !! Local Data - TYPE(Matrix_lsc) :: untransposedMatA - TYPE(Matrix_lsc) :: untransposedMatB TYPE(Matrix_ldc) :: DenseA TYPE(Matrix_ldc) :: DenseB TYPE(Matrix_ldc) :: DenseC diff --git a/Source/Fortran/dense_includes/ComposeMatrix.f90 b/Source/Fortran/dense_includes/ComposeMatrix.f90 index 9d678f28..9b7ffe3b 100644 --- a/Source/Fortran/dense_includes/ComposeMatrix.f90 +++ b/Source/Fortran/dense_includes/ComposeMatrix.f90 @@ -22,6 +22,7 @@ !! Allocate Memory CALL ConstructEmptyMatrix(out_matrix, out_columns, out_rows) + !! Copy DO JJ = 1, block_columns DO II = 1, block_rows out_matrix%data(row_offsets(II):row_offsets(II+1)-1, & diff --git a/Source/Fortran/dense_includes/ConstructMatrixSFromD.f90 b/Source/Fortran/dense_includes/ConstructMatrixSFromD.f90 index 90a37259..92f18392 100644 --- a/Source/Fortran/dense_includes/ConstructMatrixSFromD.f90 +++ b/Source/Fortran/dense_includes/ConstructMatrixSFromD.f90 @@ -32,4 +32,6 @@ END DO END IF - sparse_matrix = SMTYPE(temporary_list, rows, columns) + CALL ConstructMatrixFromTripletList(sparse_matrix, temporary_list, & + & rows, columns) + CALL DestructTripletList(temporary_list) diff --git a/Source/Fortran/sparse_includes/DenseBranch.f90 b/Source/Fortran/sparse_includes/DenseBranch.f90 index ac30b404..66c17ffe 100644 --- a/Source/Fortran/sparse_includes/DenseBranch.f90 +++ b/Source/Fortran/sparse_includes/DenseBranch.f90 @@ -1,21 +1,10 @@ - !! Handle Transposed Case - IF (IsATransposed) THEN - CALL TransposeMatrix(matA,untransposedMatA) - ELSE - CALL CopyMatrix(matA,untransposedMatA) - END IF - IF (IsBTransposed) THEN - CALL TransposeMatrix(matB,untransposedMatB) - ELSE - CALL CopyMatrix(matB,untransposedMatB) - END IF - !! Convert Forward - CALL ConstructMatrixDFromS(untransposedMatA, DenseA) - CALL ConstructMatrixDFromS(untransposedMatB, DenseB) + CALL ConstructMatrixDFromS(matA, DenseA) + CALL ConstructMatrixDFromS(matB, DenseB) !! Multiply - CALL MultiplyMatrix(DenseA, DenseB, DenseC) + CALL MultiplyMatrix(DenseA, DenseB, DenseC, & + & IsATransposed_in = IsATransposed, IsBTransposed_in = IsBTransposed) !! Convert Back CALL ConstructMatrixSFromD(DenseC, matC, threshold) diff --git a/UnitTests/test_matrix.py b/UnitTests/test_matrix.py index f7815dcb..585daa34 100644 --- a/UnitTests/test_matrix.py +++ b/UnitTests/test_matrix.py @@ -247,6 +247,90 @@ def test_multiply(self): ResultMat = mmread(self.file3) self._compare_mat(CheckMat, ResultMat) + def test_multiply_nt(self): + '''Test routines to multiply two matrices.''' + from random import uniform + for param in self.parameters: + matrix1 = param.create_matrix(complex=self.complex) + matrix2 = param.create_matrix(complex=self.complex).H + mmwrite(self.file1, matrix1) + mmwrite(self.file2, matrix2.T) + alpha = uniform(1.0, 2.0) + beta = 0.0 + if abs(beta) > 0.0: + CheckMat = alpha * matrix1.dot(matrix2) + beta * matrix1 + else: + CheckMat = alpha * matrix1.dot(matrix2) + + ntmatrix1 = self.SMatrix(self.file1) + ntmatrix2 = self.SMatrix(self.file2) + ntmatrix3 = self.SMatrix(ntmatrix2.GetRows(), + ntmatrix1.GetRows()) + memory_pool = self.MatrixMemoryPool(ntmatrix2.GetRows(), + ntmatrix1.GetRows()) + ntmatrix3.Gemm(ntmatrix1, ntmatrix2, False, True, alpha, beta, + 0.0, memory_pool) + ntmatrix3.WriteToMatrixMarket(self.file3) + + ResultMat = mmread(self.file3) + self._compare_mat(CheckMat, ResultMat) + + def test_multiply_tn(self): + '''Test routines to multiply two matrices.''' + from random import uniform + for param in self.parameters: + matrix1 = param.create_matrix(complex=self.complex) + matrix2 = param.create_matrix(complex=self.complex).H + mmwrite(self.file1, matrix1.T) + mmwrite(self.file2, matrix2) + alpha = uniform(1.0, 2.0) + beta = 0.0 + if abs(beta) > 0.0: + CheckMat = alpha * matrix1.dot(matrix2) + beta * matrix1 + else: + CheckMat = alpha * matrix1.dot(matrix2) + + ntmatrix1 = self.SMatrix(self.file1) + ntmatrix2 = self.SMatrix(self.file2) + ntmatrix3 = self.SMatrix(ntmatrix2.GetColumns(), + ntmatrix1.GetColumns()) + memory_pool = self.MatrixMemoryPool(ntmatrix2.GetColumns(), + ntmatrix1.GetColumns()) + ntmatrix3.Gemm(ntmatrix1, ntmatrix2, True, False, alpha, beta, + 0.0, memory_pool) + ntmatrix3.WriteToMatrixMarket(self.file3) + + ResultMat = mmread(self.file3) + self._compare_mat(CheckMat, ResultMat) + + def test_multiply_tt(self): + '''Test routines to multiply two matrices.''' + from random import uniform + for param in self.parameters: + matrix1 = param.create_matrix(complex=self.complex) + matrix2 = param.create_matrix(complex=self.complex).H + mmwrite(self.file1, matrix1.T) + mmwrite(self.file2, matrix2.T) + alpha = uniform(1.0, 2.0) + beta = 0.0 + if abs(beta) > 0.0: + CheckMat = alpha * matrix1.dot(matrix2) + beta * matrix1 + else: + CheckMat = alpha * matrix1.dot(matrix2) + + ntmatrix1 = self.SMatrix(self.file1) + ntmatrix2 = self.SMatrix(self.file2) + ntmatrix3 = self.SMatrix(ntmatrix2.GetRows(), + ntmatrix1.GetColumns()) + memory_pool = self.MatrixMemoryPool(ntmatrix2.GetRows(), + ntmatrix1.GetColumns()) + ntmatrix3.Gemm(ntmatrix1, ntmatrix2, True, True, alpha, beta, + 0.0, memory_pool) + ntmatrix3.WriteToMatrixMarket(self.file3) + + ResultMat = mmread(self.file3) + self._compare_mat(CheckMat, ResultMat) + def test_multiply_zero(self): '''Test routines to multiply two matrices where one is zero.''' from random import uniform