Powered by Hyde
template <class T> struct Fraction {
T a, b;
Fraction (T a = T(0), T b = T(1)): a(a), b(b) {
this->Normalize();
}
Fraction (const Fraction& f): a(f.a), b(f.b) {}
Fraction (double r, int factor = 4) {
this->b = 1;
for (int i = 0; i < factor; ++i) {
this->b *= T(10);
}
this->a = T(this->b * r + 0.5);
this->Normalize();
}
Fraction (int a = 0, int b = 1): Fraction(T(a), T(b)) {}
void Normalize () {
T d = __gcd (a, b);
a /= d, b /= d;
if (b < 0) {
a = -a, b = -b;
}
}
Fraction operator - () const {
return Fraction (-this->a, this->b);
}
Fraction& operator += (const Fraction& rhs) {
T a = this->a * rhs.b + this->b * rhs.a, b = this->b * rhs.b;
this->a = a, this->b = b;
this->Normalize();
return *this;
}
Fraction& operator += (const T& rhs) {
this->a += this->b * rhs;
return *this;
}
Fraction& operator -= (const Fraction& rhs) {
T a = this->a * rhs.b - this->b * rhs.a, b = this->b * rhs.b;
this->a = a, this->b = b;
this->Normalize();
return *this;
}
Fraction& operator -= (const T& rhs) {
this->a -= this->b * rhs;
return *this;
}
Fraction& operator *= (const Fraction& rhs) {
this->a *= rhs.a, this->b *= rhs.b;
this->Normalize();
return *this;
}
Fraction& operator *= (const T& rhs) {
this->a *= rhs;
this->Normalize();
return *this;
}
Fraction& operator /= (const Fraction& rhs) {
this->a *= rhs.b, this->b *= rhs.a;
this->Normalize();
return *this;
}
Fraction& operator /= (const T& rhs) {
this->b *= rhs;
this->Normalize();
return *this;
}
friend inline Fraction abs (const Fraction& f) {
return Fraction (abs(f.a), f.b);
}
friend inline Fraction floor (const Fraction& f) {
return f.a / f.b;
}
friend inline Fraction ceil (const Fraction& f) {
return f.b == 1 ? f.a : f.a / f.b + T(1);
}
friend inline Fraction sqrt (const Fraction& f) {
return Fraction(sqrt(f.a / double(f.b)));
}
Fraction& operator %= (const Fraction& rhs) {
T n = floor(*this / rhs);
*this -= rhs * n;
return *this;
}
friend inline Fraction operator % (Fraction lhs, const Fraction& rhs) {
lhs %= rhs;
return lhs;
}
bool operator == (const Fraction& rhs) const {
return this->a == rhs.a && this->b == rhs.b;
}
bool operator == (const T& rhs) const {
return this->a == rhs && this->b == T(1);
}
bool operator < (const Fraction& rhs) const {
return this->a * rhs.b < this->b * rhs.a;
}
bool operator < (const T& rhs) const {
return this->a < this->b * rhs;
}
friend ostream& operator << (ostream& os, Fraction <T>&& f) {
os << f.a;
if (f.b != T(1)) {
os << "/" << f.b;
}
return os;
}
friend ostream& operator << (ostream& os, Fraction <T>& f) {
os << f.a;
if (f.b != T(1)) {
os << "/" << f.b;
}
return os;
}
};template <class T> struct Point {
T x, y;
Point (T x, T y): x(x), y(y) {}
Point (const Point& p): x(p.x), y(p.y) {}
Point& operator += (const Point& p) {
x += p.x, y += p.y;
return *this;
}
Point& operator -= (const Point& p) {
x -= p.x, y -= p.y;
return *this;
}
Point& operator += (const T& c) {
x += c, y += c;
return *this;
}
Point& operator -= (const T& c) {
x -= c, y -= c;
return *this;
}
Point& operator *= (const T& c) {
x *= c, y *= c;
return *this;
}
Point& operator /= (const T& c) {
x /= c, y /= c;
return *this;
}
bool operator == (const Point& p) const {
return abs(x - p.x) <= EPS && abs(y - p.y) <= EPS;
}
bool operator < (const Point& p) const {
return x < p.x || (x == p.x && y < p.y);
}
friend ostream& operator << (ostream& os, Point <T>&& p) {
os << "(" << p.x << ", " << p.y << ")";
return os;
}
friend ostream& operator << (ostream& os, Point <T>& p) {
os << "(" << p.x << ", " << p.y << ")";
return os;
}
friend Point RotateCCW90 (Point p) {
return Point(-p.y, p.x);
}
friend Point RotateCW90 (Point p) {
return Point(p.y, -p.x);
}
// Rotate clockwise by angle t
friend Point RotateCCW (Point p, double t) {
return Point(p.x * cos(t) - p.y * sin(t), p.x * sin(t) + p.y * cos(t));
}
friend T dot (Point a, Point b) {
return a.x * b.x + a.y * b.y;
}
friend T cross (Point a, Point b) {
return a.x * b.y - a.y * b.x;
}
friend T dist2 (Point a, Point b) {
return dot(a - b, a - b);
}
// project point c on line (a, b)
friend Point ProjectPointLine (Point a, Point b, Point c) {
return a + (b - a) * dot(c - a, b - a) / dot(b - a, b - a);
}
// project point c onto line segment [a, b]
friend Point ProjectPointSegment(Point a, Point b, Point c) {
T r = dot(b - a, b - a);
if (abs(r) <= EPS) {
return a;
}
r = dot(c - a, b - a) / r;
if (r < T(0)) return a;
if (r > T(1)) return b;
return a + (b - a) * r;
}
// compute distance from c to segment [a, b]
friend T DistancePointSegment(Point a, Point b, Point c) {
return sqrt(dist2(c, ProjectOnSegment(a, b, c)));
}
// determine if lines (a, b) and (c, d) are parallel or collinear
friend bool LinesParallel(Point a, Point b, Point c, Point d) {
return abs(cross(b - a, c - d)) <= EPS;
}
friend bool LinesCollinear(Point a, Point b, Point c, Point d) {
return LinesParallel(a, b, c, d) && abs(cross(a - b, a - c)) <= EPS && abs(cross(c - d, c - a)) <= EPS;
}
// determine if line segment [a, b] intersects with line segment [c, d]
friend bool SegmentsIntersect(Point a, Point b, Point c, Point d) {
if (LinesCollinear(a, b, c, d)) {
if (dist2(a, c) <= EPS || dist2(a, d) <= EPS || dist2(b, c) <= EPS || dist2(b, d) <= EPS) {
return true;
} else if (dot(c - a, c - b) > T(0) && dot(d - a, d - b) > T(0) && dot(c - b, d - b) > T(0)) {
return false;
} else {
return true;
}
} else {
return cross(d - a, b - a) * cross(c - a, b - a) <= T(0) && cross(a - c, d - c) * cross(b - c, d - c) <= T(0);
}
}
// intersection point between lines (a, b) and (c, d)
friend Point LineLineIntersection(Point a, Point b, Point c, Point d) {
b = b - a, d = c - d, c = c - a;
assert(dot(b, b) >= EPS && dot(d, d) >= EPS);
return a + b * cross(c, d) / cross(b, d);
}
// compute center of circle given three points
friend Point CircleCenter(Point a, Point b, Point c) {
b = (a + b) / T(2);
c = (a + c) / T(2);
return LineLineIntersection(b, b + RotateCW90(a - b), c, c + RotateCW90(a - c));
}
// returns 1(0) if p is a strictly interior (exterior) point of v
friend bool PointInPolygon(const vector<Point>& v, Point p) {
bool ret = false;
for (auto q = v.begin(); q != v.end(); q++) {
auto r = (next(q) == v.end() ? v.begin() : next(q));
if ((q->y <= p.y && p.y < r->y || r->y <= p.y && p.y < q->y) &&
p.x < q->x + (r->x - q->x) * (p.y - q->y) / (r->y - q->y)) {
ret = !ret;
}
}
return ret;
}
// determine if point is on the boundary of a polygon
friend bool PointOnPolygon(const vector<Point>& v, Point p) {
for (auto q = v.begin(); q != v.end(); q++) {
auto r = (next(q) == v.end() ? v.begin() : next(q));
if (dist2(ProjectPointSegment(*q, *r, p), p) <= EPS) {
return true;
}
}
return false;
}
// compute intersection of line (a, b) with circle centered at c with radius r > 0
friend vector <Point> CircleLineIntersection(Point a, Point b, Point c, T r) {
vector <Point> ret;
b = b - a;
a = a - c;
T A = dot(b, b);
T B = dot(a, b);
T C = dot(a, a) - r * r;
T D = B * B - A * C;
if (D <= -EPS) {
return ret;
}
ret.push_back(c + a + b * (-B + sqrt(D + EPS)) / A);
if (D >= EPS) {
ret.push_back(c + a + b * (-B - sqrt(D)) / A);
}
return ret;
}
// compute intersection of circle centered at a with radius r with circle centered at b with radius R
friend vector<Point> CircleCircleIntersection(Point a, Point b, T r, T R) {
vector<Point> ret;
T d = sqrt(dist2(a, b));
if (d > r + R || d + min(r, R) < max(r, R)) {
return ret;
}
T x = (d * d - R * R + r * r) / (d * T(2));
T y = sqrt(r * r - x * x);
Point v = (b - a) / d;
ret.push_back(a + v * x + RotateCCW90(v) * y);
if (y > T(0)) {
ret.push_back(a + v * x - RotateCCW90(v) * y);
}
return ret;
}
friend T SignedArea(const vector<Point>& v) {
T area(0);
for (auto p = v.begin(); p != v.end(); p++) {
auto q = (next(p) == v.end() ? v.begin() : next(p));
area = area + (p->x * q->y - q->x * p->y);
}
return area / 2.0;
}
friend T Area(const vector<Point>& v) {
return abs(SignedArea(v));
}
friend Point Centroid(const vector<Point>& v) {
Point c(0, 0);
T scale = 6.0 * SignedArea(v);
for (auto p = v.begin(); p != v.end(); p++) {
auto q = (next(p) == v.end() ? v.begin() : next(p));
c = c + (*p + *q) * (p->x * q->y - q->x * p->y);
}
return c / scale;
}
// tests whether or not a given polygon (in CW or CCW order) is simple
friend bool IsSimple(const vector<Point>& v) {
for (auto p = v.begin(); p != v.end(); p++) {
for (auto r = next(p); r != v.end(); r++) {
auto q = (next(p) == v.end() ? v.begin() : next(p));
auto s = (next(r) == v.end() ? v.begin() : next(r));
if (p != s && q != r && SegmentsIntersect(*p, *q, *r, *s)) {
return false;
}
}
}
return true;
}
// area x 2 of triangle (a, b, c)
friend T TwiceArea (Point a, Point b, Point c) {
return cross(a, b) + cross (b, c) + cross (c, a);
}
friend T area2 (Point a, Point b, Point c) {
return cross(a, b) + cross(b, c) + cross(c, a);
}
friend void ConvexHull(vector<Point>& v) {
sort(v.begin(), v.end());
vector<Point> up, dn;
for (auto& p : v) {
while (up.size() > 1 && area2(up[up.size() - 2], up.back(), p) >= 0) {
up.pop_back();
}
while (dn.size() > 1 && area2(dn[dn.size() - 2], dn.back(), p) <= 0) {
dn.pop_back();
}
up.push_back(p);
dn.push_back(p);
}
v = dn;
v.pop_back();
reverse(dn.begin(), dn.end());
for (auto& p : dn) {
v.push_back(p);
}
}
};