81 lines
2.8 KiB
C++
81 lines
2.8 KiB
C++
#include <cmath>
|
|
#include <iostream>
|
|
#include <stdexcept>
|
|
|
|
struct Circle {
|
|
double x;
|
|
double y;
|
|
double radius;
|
|
|
|
void display() {
|
|
std::cout << "centre(" << x << ", " << y << ")" << ", radius = " << radius << std::endl;
|
|
}
|
|
};
|
|
|
|
double square(const double& value) {
|
|
return value * value;
|
|
}
|
|
|
|
/**
|
|
* The parameters s1, s2, s3 indicate whether the result circle is
|
|
* internally tangent to the corresponding argument circle (-1), or
|
|
* externally tangent to the corresponding argument circle (+1).
|
|
*/
|
|
Circle solve_apollonius(const Circle& circle1, const Circle& circle2, const Circle& circle3,
|
|
const double& s1, const double& s2, const double& s3) {
|
|
if ( abs(s1) != 1 || abs(s2) != 1 || abs(s3) != 1 ) {
|
|
throw std::invalid_argument("Parameters s1, s2, s3 must be +1 or -1");
|
|
}
|
|
|
|
const double v11 = 2 * circle2.x - 2 * circle1.x;
|
|
const double v12 = 2 * circle2.y - 2 * circle1.y;
|
|
const double v13 = square(circle1.x) - square(circle2.x) + square(circle1.y) - square(circle2.y)
|
|
- square(circle1.radius) + square(circle2.radius);
|
|
const double v14 = 2 * s2 * circle2.radius - 2 * s1 * circle1.radius;
|
|
|
|
const double v21 = 2 * circle3.x - 2 * circle2.x;
|
|
const double v22 = 2 * circle3.y - 2 * circle2.y;
|
|
const double v23 = square(circle2.x) - square(circle3.x) + square(circle2.y) - square(circle3.y)
|
|
- square(circle2.radius) + square(circle3.radius);
|
|
const double v24 = 2 * s3 * circle3.radius - 2 * s2 * circle2.radius;
|
|
|
|
const double w12 = v12 / v11;
|
|
const double w13 = v13 / v11;
|
|
const double w14 = v14 / v11;
|
|
|
|
const double w22 = v22 / v21 - w12;
|
|
const double w23 = v23 / v21 - w13;
|
|
const double w24 = v24 / v21 - w14;
|
|
|
|
const double p = -w23 / w22;
|
|
const double q = w24 / w22;
|
|
const double m = -w12 * p - w13;
|
|
const double n = w14 - w12 * q;
|
|
|
|
const double a = square(n) + square(q) - 1;
|
|
const double b = 2 * m * n - 2 * n * circle1.x + 2 * p * q - 2 * q * circle1.y + 2 * s1 * circle1.radius;
|
|
const double c = square(circle1.x) + square(m) - 2 * m * circle1.x + square(p)
|
|
+ square(circle1.y) - 2 * p * circle1.y - square(circle1.radius);
|
|
|
|
const double discriminant = b * b - 4 * a * c;
|
|
const double root = ( -b - sqrt(discriminant) ) / ( 2 * a );
|
|
|
|
const double centre_x = m + n * root;
|
|
const double centre_y = p + q * root;
|
|
|
|
return Circle(centre_x, centre_y, root);
|
|
}
|
|
|
|
int main() {
|
|
Circle circle1(0.0, 0.0, 1.0);
|
|
Circle circle2(4.0, 0.0, 1.0);
|
|
Circle circle3(2.0, 4.0, 2.0);
|
|
|
|
std::cout << "The three initial circles are:" << std::endl;
|
|
std::cout << " Circle 1: "; circle1.display();
|
|
std::cout << " Circle 2: "; circle2.display();
|
|
std::cout << " Circle 3: "; circle3.display();
|
|
std::cout << std::endl;
|
|
std::cout << "The internal circle is: "; solve_apollonius(circle1, circle2, circle3, -1.0, -1.0, -1.0).display();
|
|
std::cout << "The external circle is: "; solve_apollonius(circle1, circle2, circle3, 1.0, 1.0, 1.0).display();
|
|
}
|