test/aoj/2171/main.test.cpp
Depends on
Code
#define PROBLEM "http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=2171"
#define ERROR 1e-8
#include <iomanip>
#include <iostream>
#include <vector>
#include "Mylib/Graph/ShortestPath/warshall_floyd_for_matrix_graph.cpp"
#include "Mylib/IO/input_vector.cpp"
#include "Mylib/LinearAlgebra/simultaneous_linear_equations_float.cpp"
namespace hl = haar_lib;
int main() {
std::cin.tie(0);
std::ios::sync_with_stdio(false);
int n, s, t;
while (std::cin >> n >> s >> t, n) {
--s, --t;
auto q = hl::input_vector<int>(n);
auto g = hl::input_vector<int>(n, n);
auto dist = hl::warshall_floyd_for_matrix<int, 0>(g);
if (not dist[s][t]) {
std::cout << "impossible" << std::endl;
continue;
}
std::vector<std::vector<double>> a(n, std::vector<double>(n));
std::vector<double> b(n);
for (int i = 0; i < n; ++i) {
if (i == t) {
a[i][i] = 1;
b[i] = 0;
continue;
}
if (q[i]) {
int d = *dist[i][t];
int k = 0;
for (int j = 0; j < n; ++j) {
if (g[i][j] > 0 and g[i][j] + *dist[j][t] == d) {
++k;
}
}
for (int j = 0; j < n; ++j) {
if (g[i][j] > 0 and g[i][j] + *dist[j][t] == d) {
a[i][j] = -1;
b[i] += g[i][j];
}
}
a[i][i] += k;
} else {
int k = 0;
for (int j = 0; j < n; ++j) {
if (g[i][j] > 0) {
++k;
}
}
for (int j = 0; j < n; ++j) {
if (g[i][j] > 0) {
a[i][j] = -1;
b[i] += g[i][j];
}
}
a[i][i] += k;
}
}
auto res = hl::float_simultaneous_linear_equations(a, b, ERROR);
double ans = (*res).solution[s];
std::cout << std::setprecision(12) << std::fixed << ans << std::endl;
}
return 0;
}
#line 1 "test/aoj/2171/main.test.cpp"
#define PROBLEM "http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=2171"
#define ERROR 1e-8
#include <iomanip>
#include <iostream>
#include <vector>
#line 2 "Mylib/Graph/ShortestPath/warshall_floyd_for_matrix_graph.cpp"
#include <optional>
#line 4 "Mylib/Graph/ShortestPath/warshall_floyd_for_matrix_graph.cpp"
namespace haar_lib {
namespace warshall_floyd_for_matrix_impl {
template <typename T>
struct result {
std::vector<std::vector<std::optional<T>>> dist;
bool has_negative_cycle;
const auto& operator[](int i) const { return dist[i]; }
};
} // namespace warshall_floyd_for_matrix_impl
template <typename T, T INVALID>
auto warshall_floyd_for_matrix(const std::vector<std::vector<T>>& g) {
const int n = g.size();
auto dist = std::vector(n, std::vector<std::optional<T>>(n));
for (int i = 0; i < n; ++i) dist[i][i] = 0;
for (int i = 0; i < n; ++i) {
for (int j = 0; j < n; ++j) {
if (g[i][j] != INVALID) {
dist[i][j] = g[i][j];
}
}
}
for (int k = 0; k < n; ++k) {
for (int i = 0; i < n; ++i) {
for (int j = 0; j < n; ++j) {
if (dist[i][k] and dist[k][j]) {
if (not dist[i][j]) {
dist[i][j] = *dist[i][k] + *dist[k][j];
} else {
dist[i][j] = std::min(*dist[i][j], *dist[i][k] + *dist[k][j]);
}
}
}
}
}
bool has_negative_cycle = false;
for (int i = 0; i < n; ++i)
if (*dist[i][i] < 0) has_negative_cycle = true;
return warshall_floyd_for_matrix_impl::result<T>{dist, has_negative_cycle};
}
} // namespace haar_lib
#line 4 "Mylib/IO/input_vector.cpp"
namespace haar_lib {
template <typename T>
std::vector<T> input_vector(int N) {
std::vector<T> ret(N);
for (int i = 0; i < N; ++i) std::cin >> ret[i];
return ret;
}
template <typename T>
std::vector<std::vector<T>> input_vector(int N, int M) {
std::vector<std::vector<T>> ret(N);
for (int i = 0; i < N; ++i) ret[i] = input_vector<T>(M);
return ret;
}
} // namespace haar_lib
#line 3 "Mylib/LinearAlgebra/simultaneous_linear_equations_float.cpp"
#include <utility>
#line 5 "Mylib/LinearAlgebra/simultaneous_linear_equations_float.cpp"
namespace haar_lib {
namespace float_simultaneous_linear_equations_impl {
template <typename T>
struct result {
int rank, dim;
std::vector<T> solution;
};
} // namespace float_simultaneous_linear_equations_impl
template <typename T>
auto float_simultaneous_linear_equations(std::vector<std::vector<T>> a, std::vector<T> b, T eps) {
using result = float_simultaneous_linear_equations_impl::result<T>;
std::optional<result> ret;
const int n = a.size(), m = a[0].size();
int rank = 0;
for (int j = 0; j < m; ++j) {
int pivot = -1;
double M = eps;
for (int i = rank; i < n; ++i) {
if (std::abs(a[i][j]) > M) {
M = std::abs(a[i][j]);
pivot = i;
}
}
if (pivot == -1) continue;
std::swap(a[pivot], a[rank]);
std::swap(b[pivot], b[rank]);
{
double d = a[rank][j];
for (int k = 0; k < m; ++k) a[rank][k] /= d;
b[rank] /= d;
}
for (int i = 0; i < n; ++i) {
if (i == rank or std::abs(a[i][j]) <= eps) continue;
double d = a[i][j];
for (int k = 0; k < m; ++k) {
a[i][k] -= a[rank][k] * d;
}
b[i] -= b[rank] * d;
}
++rank;
}
for (int i = rank; i < n; ++i) {
if (std::abs(b[i]) > eps) {
return ret;
}
}
const int dim = m - rank;
std::vector<T> solution(m);
for (int i = 0; i < rank; ++i) solution[i] = b[i];
ret = result({rank, dim, solution});
return ret;
}
} // namespace haar_lib
#line 10 "test/aoj/2171/main.test.cpp"
namespace hl = haar_lib;
int main() {
std::cin.tie(0);
std::ios::sync_with_stdio(false);
int n, s, t;
while (std::cin >> n >> s >> t, n) {
--s, --t;
auto q = hl::input_vector<int>(n);
auto g = hl::input_vector<int>(n, n);
auto dist = hl::warshall_floyd_for_matrix<int, 0>(g);
if (not dist[s][t]) {
std::cout << "impossible" << std::endl;
continue;
}
std::vector<std::vector<double>> a(n, std::vector<double>(n));
std::vector<double> b(n);
for (int i = 0; i < n; ++i) {
if (i == t) {
a[i][i] = 1;
b[i] = 0;
continue;
}
if (q[i]) {
int d = *dist[i][t];
int k = 0;
for (int j = 0; j < n; ++j) {
if (g[i][j] > 0 and g[i][j] + *dist[j][t] == d) {
++k;
}
}
for (int j = 0; j < n; ++j) {
if (g[i][j] > 0 and g[i][j] + *dist[j][t] == d) {
a[i][j] = -1;
b[i] += g[i][j];
}
}
a[i][i] += k;
} else {
int k = 0;
for (int j = 0; j < n; ++j) {
if (g[i][j] > 0) {
++k;
}
}
for (int j = 0; j < n; ++j) {
if (g[i][j] > 0) {
a[i][j] = -1;
b[i] += g[i][j];
}
}
a[i][i] += k;
}
}
auto res = hl::float_simultaneous_linear_equations(a, b, ERROR);
double ans = (*res).solution[s];
std::cout << std::setprecision(12) << std::fixed << ans << std::endl;
}
return 0;
}
Back to top page