mirror of
https://github.com/Z3Prover/z3
synced 2025-04-27 10:55:50 +00:00
268 lines
8 KiB
C++
268 lines
8 KiB
C++
/*++
|
|
Copyright (c) 2017 Microsoft Corporation
|
|
|
|
Module Name:
|
|
|
|
<name>
|
|
|
|
Abstract:
|
|
|
|
<abstract>
|
|
|
|
Author:
|
|
|
|
Lev Nachmanson (levnach)
|
|
|
|
Revision History:
|
|
|
|
|
|
--*/
|
|
#include <algorithm>
|
|
#include "math/lp/scaler.h"
|
|
#include "math/lp/numeric_pair.h"
|
|
namespace lp {
|
|
// for scaling an LP
|
|
template <typename T, typename X> T scaler<T, X>::right_side_balance() {
|
|
T ret = zero_of_type<T>();
|
|
unsigned i = m_A.row_count();
|
|
while (i--) {
|
|
T rs = abs(convert_struct<T, X>::convert(m_b[i]));
|
|
if (!is_zero<T>(rs)) {
|
|
numeric_traits<T>::log(rs);
|
|
ret += rs * rs;
|
|
}
|
|
}
|
|
return ret;
|
|
}
|
|
|
|
template <typename T, typename X> T scaler<T, X>::A_min() const {
|
|
T min = zero_of_type<T>();
|
|
for (unsigned i = 0; i < m_A.row_count(); i++) {
|
|
T t = m_A.get_min_abs_in_row(i);
|
|
min = i == 0 ? t : std::min(t, min);
|
|
}
|
|
return min;
|
|
}
|
|
|
|
template <typename T, typename X> T scaler<T, X>::A_max() const {
|
|
T max = zero_of_type<T>();
|
|
for (unsigned i = 0; i < m_A.row_count(); i++) {
|
|
T t = m_A.get_max_abs_in_row(i);
|
|
max = i == 0? t : std::max(t, max);
|
|
}
|
|
return max;
|
|
}
|
|
|
|
template <typename T, typename X> T scaler<T, X>::get_A_ratio() const {
|
|
T min = A_min();
|
|
T max = A_max();
|
|
lp_assert(!m_settings.abs_val_is_smaller_than_zero_tolerance(min));
|
|
T ratio = max / min;
|
|
return ratio;
|
|
}
|
|
|
|
template <typename T, typename X> T scaler<T, X>::get_max_ratio_on_rows() const {
|
|
T ret = T(1);
|
|
unsigned i = m_A.row_count();
|
|
while (i--) {
|
|
T den = m_A.get_min_abs_in_row(i);
|
|
lp_assert(!m_settings.abs_val_is_smaller_than_zero_tolerance(den));
|
|
T t = m_A.get_max_abs_in_row(i)/ den;
|
|
if (t > ret)
|
|
ret = t;
|
|
}
|
|
return ret;
|
|
}
|
|
|
|
template <typename T, typename X> T scaler<T, X>::get_max_ratio_on_columns() const {
|
|
T ret = T(1);
|
|
unsigned i = m_A.column_count();
|
|
while (i--) {
|
|
T den = m_A.get_min_abs_in_column(i);
|
|
if (m_settings.abs_val_is_smaller_than_zero_tolerance(den))
|
|
continue; // got a zero column
|
|
T t = m_A.get_max_abs_in_column(i)/den;
|
|
if (t > ret)
|
|
ret = t;
|
|
}
|
|
return ret;
|
|
}
|
|
|
|
template <typename T, typename X> void scaler<T, X>::scale_rows_with_geometric_mean() {
|
|
unsigned i = m_A.row_count();
|
|
while (i--) {
|
|
T max = m_A.get_max_abs_in_row(i);
|
|
T min = m_A.get_min_abs_in_row(i);
|
|
lp_assert(max > zero_of_type<T>() && min > zero_of_type<T>());
|
|
if (is_zero(max) || is_zero(min))
|
|
continue;
|
|
T gm = T(sqrt(numeric_traits<T>::get_double(max*min)));
|
|
if (m_settings.is_eps_small_general(gm, 0.01)) {
|
|
continue;
|
|
}
|
|
m_A.multiply_row(i, one_of_type<T>() / gm);
|
|
m_b[i] /= gm;
|
|
}
|
|
}
|
|
|
|
template <typename T, typename X> void scaler<T, X>::scale_columns_with_geometric_mean() {
|
|
unsigned i = m_A.column_count();
|
|
while (i--) {
|
|
T max = m_A.get_max_abs_in_column(i);
|
|
T min = m_A.get_min_abs_in_column(i);
|
|
T den = T(sqrt(numeric_traits<T>::get_double(max*min)));
|
|
if (m_settings.is_eps_small_general(den, 0.01))
|
|
continue; // got a zero column
|
|
T gm = T(1)/ den;
|
|
T cs = m_column_scale[i] * gm;
|
|
if (m_settings.is_eps_small_general(cs, 0.1))
|
|
continue;
|
|
m_A.multiply_column(i, gm);
|
|
m_column_scale[i] = cs;
|
|
}
|
|
}
|
|
|
|
template <typename T, typename X> void scaler<T, X>::scale_once_for_ratio() {
|
|
T max_ratio_on_rows = get_max_ratio_on_rows();
|
|
T max_ratio_on_columns = get_max_ratio_on_columns();
|
|
bool scale_rows_first = max_ratio_on_rows > max_ratio_on_columns;
|
|
// if max_ratio_on_columns is the largest then the rows are in worse shape than columns
|
|
if (scale_rows_first) {
|
|
scale_rows_with_geometric_mean();
|
|
scale_columns_with_geometric_mean();
|
|
} else {
|
|
scale_columns_with_geometric_mean();
|
|
scale_rows_with_geometric_mean();
|
|
}
|
|
}
|
|
|
|
template <typename T, typename X> bool scaler<T, X>::scale_with_ratio() {
|
|
T ratio = get_A_ratio();
|
|
// The ratio is greater than or equal to one. We would like to diminish it and bring it as close to 1 as possible
|
|
unsigned reps = m_settings.reps_in_scaler;
|
|
do {
|
|
scale_once_for_ratio();
|
|
T new_r = get_A_ratio();
|
|
if (new_r >= T(0.9) * ratio)
|
|
break;
|
|
} while (reps--);
|
|
|
|
bring_rows_and_columns_maximums_to_one();
|
|
return true;
|
|
}
|
|
|
|
template <typename T, typename X> void scaler<T, X>::bring_row_maximums_to_one() {
|
|
unsigned i = m_A.row_count();
|
|
while (i--) {
|
|
T t = m_A.get_max_abs_in_row(i);
|
|
if (m_settings.abs_val_is_smaller_than_zero_tolerance(t)) continue;
|
|
m_A.multiply_row(i, one_of_type<T>() / t);
|
|
m_b[i] /= t;
|
|
}
|
|
}
|
|
|
|
template <typename T, typename X> void scaler<T, X>::bring_column_maximums_to_one() {
|
|
unsigned i = m_A.column_count();
|
|
while (i--) {
|
|
T max = m_A.get_max_abs_in_column(i);
|
|
if (m_settings.abs_val_is_smaller_than_zero_tolerance(max)) continue;
|
|
T t = T(1) / max;
|
|
m_A.multiply_column(i, t);
|
|
m_column_scale[i] *= t;
|
|
}
|
|
}
|
|
|
|
template <typename T, typename X> void scaler<T, X>::bring_rows_and_columns_maximums_to_one() {
|
|
if (get_max_ratio_on_rows() > get_max_ratio_on_columns()) {
|
|
bring_row_maximums_to_one();
|
|
bring_column_maximums_to_one();
|
|
} else {
|
|
bring_column_maximums_to_one();
|
|
bring_row_maximums_to_one();
|
|
}
|
|
}
|
|
|
|
template <typename T, typename X> bool scaler<T, X>::scale_with_log_balance() {
|
|
T balance = get_balance();
|
|
T balance_before_scaling = balance;
|
|
// todo : analyze the scale order : rows-columns, or columns-rows. Iterate if needed
|
|
for (int i = 0; i < 10; i++) {
|
|
scale_rows();
|
|
scale_columns();
|
|
T nb = get_balance();
|
|
if (nb < T(0.9) * balance) {
|
|
balance = nb;
|
|
} else {
|
|
balance = nb;
|
|
break;
|
|
}
|
|
}
|
|
return balance <= balance_before_scaling;
|
|
}
|
|
// Returns true if and only if the scaling was successful.
|
|
// It is the caller responsibility to restore the matrix
|
|
template <typename T, typename X> bool scaler<T, X>::scale() {
|
|
if (numeric_traits<T>::precise()) return true;
|
|
if (m_settings.scale_with_ratio)
|
|
return scale_with_ratio();
|
|
return scale_with_log_balance();
|
|
}
|
|
|
|
template <typename T, typename X> void scaler<T, X>::scale_rows() {
|
|
for (unsigned i = 0; i < m_A.row_count(); i++)
|
|
scale_row(i);
|
|
}
|
|
|
|
template <typename T, typename X> void scaler<T, X>::scale_row(unsigned i) {
|
|
T row_max = std::max(m_A.get_max_abs_in_row(i), abs(convert_struct<T, X>::convert(m_b[i])));
|
|
T alpha = numeric_traits<T>::one();
|
|
if (numeric_traits<T>::is_zero(row_max)) {
|
|
return;
|
|
}
|
|
if (row_max < m_scaling_minimum) {
|
|
do {
|
|
alpha *= T(2);
|
|
row_max *= T(2);
|
|
} while (row_max < m_scaling_minimum);
|
|
m_A.multiply_row(i, alpha);
|
|
m_b[i] *= alpha;
|
|
} else if (row_max > m_scaling_maximum) {
|
|
do {
|
|
alpha /= T(2);
|
|
row_max /= T(2);
|
|
} while (row_max > m_scaling_maximum);
|
|
m_A.multiply_row(i, alpha);
|
|
m_b[i] *= alpha;
|
|
}
|
|
}
|
|
|
|
template <typename T, typename X> void scaler<T, X>::scale_column(unsigned i) {
|
|
T column_max = m_A.get_max_abs_in_column(i);
|
|
T alpha = numeric_traits<T>::one();
|
|
|
|
if (numeric_traits<T>::is_zero(column_max)){
|
|
return; // the column has zeros only
|
|
}
|
|
if (column_max < m_scaling_minimum) {
|
|
do {
|
|
alpha *= T(2);
|
|
column_max *= T(2);
|
|
} while (column_max < m_scaling_minimum);
|
|
} else if (column_max > m_scaling_maximum) {
|
|
do {
|
|
alpha /= T(2);
|
|
column_max /= T(2);
|
|
} while (column_max > m_scaling_maximum);
|
|
} else {
|
|
return;
|
|
}
|
|
m_A.multiply_column(i, alpha);
|
|
m_column_scale[i] = alpha;
|
|
}
|
|
|
|
template <typename T, typename X> void scaler<T, X>::scale_columns() {
|
|
for (unsigned i = 0; i < m_A.column_count(); i++) {
|
|
scale_column(i);
|
|
}
|
|
}
|
|
}
|