Method Name: makeRandomDD
Author: Andrew Aposhian
Language: C++
This method must be called on an instance of a DenseArray of doubles. See the docs for DenseArray for details about the DenseArray class template.
Description/Purpose: The purpose of this method is to set a DenseArray to contain entries set to random real floating point values in the range [min, max], such that the matrix is diagonally dominant. The min and max arguments have default values of 0.0 and 1.0 respectively.
Input:
min
: a double representing the minimum value to which an entry in the matrix can be setmax
: a double representing the maximum value to which an entry in the matrix can be set. This must be strictly greater thanmin
.
Output: None
Usage/Example: The example below shows instantiating a DenseArray of doubles, calling makeRandomDD on that DenseArray instance with the default arguments, then printing the array.
std::cout << "Random DD Matrix: " << std::endl;
DenseArray<double>* ddA = new DenseArray<double>(5, 5);
ddA->makeRandomDD(-5, 5);
ddA->print();
Output from lines above:
Random DD Matrix:
-4.96466 2.24876 -0.278311 0.0989411 -2.29522
2.23385 -4.99332 0.0213465 -1.34278 1.39068
0.0131428 -4.01333 -4.99231 0.566554 -0.394448
-4.16298 0.473458 0.293169 4.99812 0.0512707
-0.397526 1.11989 -1.12894 -2.08432 4.9713
Implementation/Code: See DenseArray.ipp on GitHub
template <typename T>
void DenseArray<T>::makeRandomDD(const double min, const double max) {
assertSquare();
assertProperMinMax(min, max);
const double ABS_OF_MIN = abs(min);
const double ABS_OF_MAX = abs(max);
const double MAX_ABS = (ABS_OF_MIN > ABS_OF_MAX) ? ABS_OF_MIN : ABS_OF_MAX;
double minAbs;
if (min <= 0 && max >= 0) {
minAbs = 0;
}
else {
minAbs = (ABS_OF_MIN > ABS_OF_MAX) ? ABS_OF_MAX : ABS_OF_MIN;
}
const double MIN_ABS = minAbs;
// Assert value range is valid
if (MIN_ABS * (colDim() - 1) >= MAX_ABS) {
std::cout << "Invalid value range to make diagonally dominant array with n = ";
std::cout << colDim() << " columns: " << std::endl;
std::cout << "(min absolute val=" << MIN_ABS << ")";
std::cout << " * (n - 1) >= (max absolute val=" << MAX_ABS << ")" << std::endl;
throw std::exception();
}
initializeEmptyArray();
const double OFF_DIAG_COUNT = colDim() - 1;
// For each row
for (unsigned int i = 0; i < rowDim(); i++) {
// Get off-diagonal values for row
std::vector<double> rowVals(OFF_DIAG_COUNT);
double absRowSum = 0;
for (unsigned int j = 0; j < OFF_DIAG_COUNT; j++) {
// Find bounds for value
double minSumOfRest = absRowSum + ((OFF_DIAG_COUNT - j - 1) * MIN_ABS);
double maxAbsOfVal = MAX_ABS - minSumOfRest;
double minBound;
if (min >= 0 || maxAbsOfVal > ABS_OF_MIN) {
minBound = min;
}
else {
minBound = -maxAbsOfVal;
}
double maxBound;
if (max <= 0 || maxAbsOfVal > ABS_OF_MAX) {
maxBound = max;
}
else {
maxBound = maxAbsOfVal;
}
// Get a random value
double randValue;
do {
randValue = getRandDouble(minBound, maxBound);
// Continue to loop to look for a different random value
// if the new absolute row sum equals MAX_ABS.
// This guarantees strict diagonal dominance.
// >= instead of == just for extra safety with
// floating point arithmetic.
} while (abs(randValue) + absRowSum >= MAX_ABS);
// Save the random value in rowVals
rowVals[j] = randValue;
// Update the absolute row sum
absRowSum += abs(randValue);
}
// Randomly permute rowVals
std::random_shuffle(rowVals.begin(), rowVals.end());
// Set array values in off-diagonal of row
for (unsigned int j = 0; j < colDim(); j++) {
if (j == i) {
// Skip the diagonal
continue;
}
set(i, j, rowVals.back());
rowVals.pop_back();
}
// Get sign of diagonal value:
// ABS_OF_MIN > absRowSum and not ABS_OF_MAX > absRowSum
// implies the random value must be negative.
// ABS_OF_MAX > absRowSum and not ABS_OF_MIN > absRowSum
// implies the random value must be positive.
// If ABS_OF_MIN > absRowSum and ABS_OF_MAX > absRowSum, then
// the random value can be positive or negative, so we randomly
// select whether it will be positive or negative.
bool positive;
if (ABS_OF_MIN > absRowSum && ABS_OF_MAX > absRowSum) {
// Randomly select whether value will be positive or negative
double signSelector = getRandDouble(-1.0, 1.0);
positive = (signSelector > 0) ? true : false;
}
else if (ABS_OF_MIN > absRowSum) {
positive = false;
}
else if (ABS_OF_MAX > absRowSum) {
positive = true;
}
else {
std::cout << "min and max are equal" << std::endl;
throw std::exception();
}
// Find a random value between min and max with an absolute
// value greater than absRowSum.
double randMaxAbs;
do {
if (positive) {
randMaxAbs = getRandDouble(absRowSum, ABS_OF_MAX);
}
else {
randMaxAbs = getRandDouble(-absRowSum, -ABS_OF_MIN);
}
// Continue to loop to look for a random value
// if randMaxAbs value has an absolute value equal
// to absRowSum. This guarantees strict diagonal dominance.
// The probablity of this looping is low.
// <= instead of == just for extra safety with
// floating point arithmetic.
} while (abs(randMaxAbs) <= absRowSum);
// Set array value on diagonal to the random value
set(i, i, randMaxAbs);
}
}
Last Modified: February/2019