equal
deleted
inserted
replaced
643 double v1 = 2 * u1 - 1; |
643 double v1 = 2 * u1 - 1; |
644 double v2 = 2 * u2 - 1; |
644 double v2 = 2 * u2 - 1; |
645 double w = v1 * v1 + v2 * v2; |
645 double w = v1 * v1 + v2 * v2; |
646 if (w <= 1.0) |
646 if (w <= 1.0) |
647 { // Got good pair |
647 { // Got good pair |
648 double y = sqrt ((-2 * log (w)) / w); |
648 double y = std::sqrt ((-2 * std::log (w)) / w); |
649 m_next = mean + v2 * y * sqrt (variance); |
649 m_next = mean + v2 * y * std::sqrt (variance); |
650 // if next is in bounds, it is valid |
650 // if next is in bounds, it is valid |
651 m_nextValid = fabs (m_next - mean) <= bound; |
651 m_nextValid = std::fabs (m_next - mean) <= bound; |
652 double x1 = mean + v1 * y * sqrt (variance); |
652 double x1 = mean + v1 * y * std::sqrt (variance); |
653 // if x1 is in bounds, return it |
653 // if x1 is in bounds, return it |
654 if (fabs (x1 - mean) <= bound) |
654 if (std::fabs (x1 - mean) <= bound) |
655 { |
655 { |
656 return x1; |
656 return x1; |
657 } |
657 } |
658 // otherwise try and return m_next if it is valid |
658 // otherwise try and return m_next if it is valid |
659 else if (m_nextValid) |
659 else if (m_nextValid) |
861 return GetValue (1.0 + alpha, beta) * std::pow (u, 1.0 / alpha); |
861 return GetValue (1.0 + alpha, beta) * std::pow (u, 1.0 / alpha); |
862 } |
862 } |
863 |
863 |
864 double x, v, u; |
864 double x, v, u; |
865 double d = alpha - 1.0 / 3.0; |
865 double d = alpha - 1.0 / 3.0; |
866 double c = (1.0 / 3.0) / sqrt (d); |
866 double c = (1.0 / 3.0) / std::sqrt (d); |
867 |
867 |
868 while (1) |
868 while (1) |
869 { |
869 { |
870 do |
870 do |
871 { |
871 { |
888 } |
888 } |
889 if (u < 1 - 0.0331 * x * x * x * x) |
889 if (u < 1 - 0.0331 * x * x * x * x) |
890 { |
890 { |
891 break; |
891 break; |
892 } |
892 } |
893 if (log (u) < 0.5 * x * x + d * (1 - v + log (v))) |
893 if (std::log (u) < 0.5 * x * x + d * (1 - v + std::log (v))) |
894 { |
894 { |
895 break; |
895 break; |
896 } |
896 } |
897 } |
897 } |
898 |
898 |
938 double v1 = 2 * u1 - 1; |
938 double v1 = 2 * u1 - 1; |
939 double v2 = 2 * u2 - 1; |
939 double v2 = 2 * u2 - 1; |
940 double w = v1 * v1 + v2 * v2; |
940 double w = v1 * v1 + v2 * v2; |
941 if (w <= 1.0) |
941 if (w <= 1.0) |
942 { // Got good pair |
942 { // Got good pair |
943 double y = sqrt ((-2 * log (w)) / w); |
943 double y = std::sqrt ((-2 * std::log (w)) / w); |
944 m_next = mean + v2 * y * sqrt (variance); |
944 m_next = mean + v2 * y * std::sqrt (variance); |
945 // if next is in bounds, it is valid |
945 // if next is in bounds, it is valid |
946 m_nextValid = fabs (m_next - mean) <= bound; |
946 m_nextValid = std::fabs (m_next - mean) <= bound; |
947 double x1 = mean + v1 * y * sqrt (variance); |
947 double x1 = mean + v1 * y * std::sqrt (variance); |
948 // if x1 is in bounds, return it |
948 // if x1 is in bounds, return it |
949 if (fabs (x1 - mean) <= bound) |
949 if (std::fabs (x1 - mean) <= bound) |
950 { |
950 { |
951 return x1; |
951 return x1; |
952 } |
952 } |
953 // otherwise try and return m_next if it is valid |
953 // otherwise try and return m_next if it is valid |
954 else if (m_nextValid) |
954 else if (m_nextValid) |
1125 } |
1125 } |
1126 |
1126 |
1127 // Calculate the triangular random variable. |
1127 // Calculate the triangular random variable. |
1128 if (u <= (mode - min) / (max - min) ) |
1128 if (u <= (mode - min) / (max - min) ) |
1129 { |
1129 { |
1130 return min + sqrt (u * (max - min) * (mode - min) ); |
1130 return min + std::sqrt (u * (max - min) * (mode - min) ); |
1131 } |
1131 } |
1132 else |
1132 else |
1133 { |
1133 { |
1134 return max - sqrt ( (1 - u) * (max - min) * (max - mode) ); |
1134 return max - std::sqrt ( (1 - u) * (max - min) * (max - mode) ); |
1135 } |
1135 } |
1136 } |
1136 } |
1137 |
1137 |
1138 uint32_t |
1138 uint32_t |
1139 TriangularRandomVariable::GetInteger (uint32_t mean, uint32_t min, uint32_t max) |
1139 TriangularRandomVariable::GetInteger (uint32_t mean, uint32_t min, uint32_t max) |