How to generate a n-sized random float array that sums up to 0.0?

Question

Consider that I need a n-sized vector where each element is defined between [-1,1]. The element a[i] is a float generated by -1 + 2*rand(). I need a elegant way to ensure that the sum of the elements of my array is equal to zero.

I've found two possible solutions:

The first one is this matlab function https://www.mathworks.com/matlabcentral/fileexchange/9700-random-vectors-with-fixed-sum. It has also a implementation in R, however it is too much work to implement it on C, since this function is used for a 2d array.

The second one is provided in this thread here: Generate random values with fixed sum in C++. Essentially, the idea is to generate n numbers with a normal distribution then normalize them to with my sum. (I have implemented it using python bellow) for a vector with sum up to 1.0. It works for every sum value except for zero.

import random as rd

mySum = 1;
randomVector = []
randomSum = 0

for i in range(7):
    randomNumber = -1 + 2*rd.random()
    randomVector.append(randomNumber)
    randomSum  += randomNumber

coef = mySum/randomSum
myNewList = [j * coef for j in randomVector]
newsum = sum(myNewList)

So, is there a way to do that using C or C++? If you know a already implemented function it would be awesome. Thanks.

Answer 1

I figured out a solution to your problem. This is not perfect since its randomness is limited by the range requirement.

The strategy is:

1. Define a function able to generate a random float in a customizable range. No need to reinvent the wheel: I borrowed it from https://stackoverflow.com/a/44105089/11336762
2. Malloc array (I omit pointer check in my example) and initialize the seed. In my example I just used current time but it can be improved
3. For every element to be generated, pre-calculate random range. Given the i-th sum, make sure that the next sum is NEVER out of range: if the sum is positive, the range needs to be (-1,1-sum); if it is negative it the range needs to be (-1-sum,1)
4. Do this until (n-1)th element. Last element must be directly assigned as the sum with the sign changed.

    #include<stdio.h>
    #include<stdlib.h>
    #include<time.h>

    float float_rand( float min, float max )
    {
        float scale = rand() / (float) RAND_MAX; /* [0, 1.0] */
        return min + scale * ( max - min );      /* [min, max] */
    }

    void main( int argc, char *argv[] )
    {
        if( argc == 2 )
        {
            int i, n = atoi ( argv[1] );
            float *outArr = malloc( n * sizeof( float ) );
            float sum = 0;

            printf( "Input value: %d\n\n", n );

            /* Initialize seed */
            srand ( time( NULL ) );

            for( i=0; i<n-1; i++ )
            {
                /* Limit random generation range in order to make sure the next sum is  *
                 * not outside (-1,1) range.                                            */
                float min = (sum<0? -1-sum : -1);
                float max = (sum>0? 1-sum : 1);

                outArr[i] = float_rand( min, max );
                sum += outArr[i];
            }

            /* Set last array element */
            outArr[n-1] = -sum;

            /* Print results */
            sum=0;
            for( i=0; i<n; i++ )
            {
                sum += outArr[i];
                printf( "  outArr[%d]=%f \t(sum=%f)\n", i, outArr[i], sum );
            }

            free( outArr );
        }  
        else
        {
          printf( "Only a parameter allowed (integer N)\n" );
        }
    }

I tried it, and it works also when n=1. In case of n=0 a sanity check should be added to my example.

Some output examples:

N=1:

Input value: 1

  outArr[0]=-0.000000   (sum=-0.000000)

N=4

Input value: 4

  outArr[0]=-0.804071   (sum=-0.804071)
  outArr[1]=0.810685    (sum=0.006614)
  outArr[2]=-0.353444   (sum=-0.346830)
  outArr[3]=0.346830    (sum=0.000000)

N=8:

Input value: 8

  outArr[0]=-0.791314   (sum=-0.791314)
  outArr[1]=0.800182    (sum=0.008867)
  outArr[2]=-0.571293   (sum=-0.562426)
  outArr[3]=0.293300    (sum=-0.269126)
  outArr[4]=-0.082886   (sum=-0.352012)
  outArr[5]=0.818639    (sum=0.466628)
  outArr[6]=-0.301473   (sum=0.165155)
  outArr[7]=-0.165155   (sum=0.000000)

Answer 2

Thank you guys again for the help.

So, based on the idea of Cryostasys I developed the following C code to solve my problem:

#include <stdio.h>      /* printf, scanf, puts, NULL */
#include <stdlib.h>     /* srand, rand */
#include <time.h>       /* time */
#include <math.h>

int main()
{
    int arraySize = 10; //input value 
    double createdArray[arraySize]; //output value

    double randomPositiveVector[arraySize];
    double randomNegativeVector[arraySize];
    double positiveSum = 0.;
    double negativeSum = 0.;

    srand(time(NULL)); //seed for random generation

    for(int i = 0; i < arraySize; ++i)
    {
        double randomNumber = -1.+2.*rand()/((double) RAND_MAX); //random in [-1.0,1.0]
        printf("%f\n",randomNumber);
        if(randomNumber >=0)
        {
            randomPositiveVector[i] = randomNumber;
            positiveSum += randomNumber;
        }
        else
        {
            randomNegativeVector[i] = randomNumber;
            negativeSum += randomNumber;
        }
    }
    if(positiveSum == 0. || negativeSum == 0.) printf("ERROR\n");

    double positiveCoefficient =  1.0/positiveSum;
    double negativeCoefficient = -1.0/negativeSum;
    for(int i = 0; i < arraySize; ++i)
    {
        randomPositiveVector[i] = positiveCoefficient * randomPositiveVector[i];
        randomNegativeVector[i] = negativeCoefficient * randomNegativeVector[i];
        if(fabs(randomPositiveVector[i]) > 1e-6) //near to zero 
        {
            createdArray[i] = randomPositiveVector[i];
        }
        else
        {
            createdArray[i] = randomNegativeVector[i];
        }
    }

    for(int i = 0; i < arraySize; ++i)
    {
        printf("createdArray[%d] = %9f\n",i,createdArray[i]);

    }

    return(0);
}

Please note that the randomness of the values generated is decreased, as mentioned in the comments of the question. Also, the kind of random distribution is determined by the function that you use to generate the randomNumber above. In this case, I've used rand() from stdlib.h which is based on giving a seed to the function and it is going to generate a pseudo-random number. You could use a different option, for instance, drand48() from stdlib.h as well.

Nevertheless, it is required that at least one positive and one negative value is generated in order to this code work. One verification step was added to the code, and if it reaches this condition one should run again the code or do something about.

Output example (arraySize = 10):

createdArray[0] = -0.013824
createdArray[1] =  0.359639
createdArray[2] = -0.005851
createdArray[3] =  0.126829
createdArray[4] = -0.334745
createdArray[5] = -0.473096
createdArray[6] = -0.172484
createdArray[7] =  0.249523
createdArray[8] =  0.262370
createdArray[9] =  0.001640

Answer 3

One option is to generate some samples and then scale their values around the average. In C++ it would be something like the following

#include <iostream>
#include <iomanip>
#include <random>
#include <algorithm>
#include <cmath>

int main()
{
    std::random_device rd;
    std::seed_seq ss{rd(), rd(), rd(), rd()};
    std::mt19937 gen{ss};

    const int samples = 9;

    // Generates the samples in [0, 2]
    std::uniform_real_distribution dist(0.0, std::nextafter(2.0, 3.0));
    std::vector<double> nums(samples);
    double sum = 0.0;
    for ( auto & i : nums )
    {
        i = dist(gen);
        sum += i;
    }
    double average = sum / samples;
    double k = 1.0 / std::max(average, 2.0 - average);

    // Transform the values (apart from the last) to meet the requirements
    sum = 0.0;
    for ( size_t i = 0; i < nums.size() - 1; ++i )
    {
        nums[i] = (nums[i] - average) * k;
        sum += nums[i];
    };
    // This trick (to ensure the needed precision) only works if the sum
    // is always evaluated in the same order
    nums.back() = 0.0 - sum;

    sum = 0.0;
    for ( size_t i = 0; i < nums.size(); ++i )
    {
        sum += nums[i];
        std::cout << std::setw(10) << std::fixed << nums[i] << '\n';
    }
    if (sum != 0.0)
        std::cout << "Failed.\n";
}

Testable here.