You should create a point struct this makes things a bit easier. It has an double x, double y, double z coordinate. This point will be immensely useful, as we can define a line segment as a pair<point,point>, a triangular section of a plane as a vector<point> of 3 points, and a sphere as a center + a radius. The point is the most basic unit in geometry, and it is a good idea to package it compactly in a struct.

In the next 5 chapters, we will show how to handle cases where an object falls onto another object, but these objects are simpler: points, lines (When I refer to a line, I mean a line segment), planes (I really refer to a triangular section of a plane), or spheres. Note that a cube or tetrahedron is just composed of a bunch of these simple shapes.

It turns out there are only 10 types of collisions you need to consider:

• A point falling onto a plane [CH 2]
• A plane falling onto a point [CH 2]
• A line falling onto another line [CH 3]
• A point falling onto sphere [CH 4]
• A sphere falling onto a point [CH 4]
• A line falling onto a sphere [CH 5]
• A sphere falling onto a line [CH 5]
• A plane falling onto a sphere [CH 6]
• A sphere falling onto a plane [CH 6]
• A sphere falling onto a sphere [CH 7].

The you can convince yourself that you don't have to consider other things. For example: a line falling onto a plane is encapsulated by a line falling onto a line if it falls onto an edge of the plane, or a point falling onto a plane if the line lands on the plane.

We have methods with parameters (center, object 1, object 2). object 1 and the center are in their own coordinate system and are rigid and will fall together. Then, we place object 1 with the center in object 2's coordinate system at a very tall height and make it fall. As object 1 falls, it's center will also fall that exact distance, until, object 1 collides with object 2. (The center has no collision). Then, the method will output the center's z coordinate in object 2's coordinate system. If object 1 falls through object 2, then the method returns -1.

For example, if the center is $$$(0, 0, 3)$$$, object 1 is a point $$$(0, 0, 2)$$$, and object 2 is the plane $$$(-1, -1, 0)$$$, $$$(-1, 1, 0)$$$, $$$(1, 0, 0)$$$. (this is a flat triangle on the x-y plane that encloses the origin), the method should return 1, since object 1 will fall onto the origin, and the center, always being $$$(3-2) = 1$$$ higher in the z-direction, would be at the point $$$(0, 0, 1)$$$, with a z coordinate of 1.

There are several helper methods to implement that will be helpful:

point add(point a, point b). This method returns a new point with the sum of the x, y, z coordinates of 2 points. This is adding 2 points $$$a+b$$$.

point sub(point a, point b). This method returns a new point with $$$x = a_x - b_x$$$, and similarly for y and z. This is subtracting 2 points $$$a-b$$$.

point mul(point a, double b). This method returns a new point with $$$x = a_x \cdot b$$$. This is multiplying $$$a \cdot b$$$.

point cross(point a, point b). This treats a and b like vectors, and returns their cross product. (Refer to https://en.wikipedia.org/wiki/Cross_product)

int pointLoc(point a, point b, point c). This method determines if you stand at a, look at b, whether c would be collinear, to the left of you, or to the right of you. (1 if to the left, 0 if collinear, -1 if to the right).

This is fairly simple. You can check the quantity $$$(c_y - a_y)(b_x - a_x) - (b_y - c_y)(c_x - a_x)$$$ whether it is positive, negative or zero. Positive means that C is to the left of AB. negative means C is to the right of AB, and zero means C is collinear with AB. I am not sure why this works. Perhaps it is an interesting property of the cross product, but you don't need to know to implement point location test.

double dist2D(point a, point b). Returns the distance of the projection of the points on the 2d plane, ignores the z coordinates. The distance is $$$\sqrt{(a_x - b_x)^2 + (a_y - b_y)^2}$$$.

We will be doing the point falls onto a plane case.

double pointplane(point center, point p, vector<point> pl). This method returns the z coordinate of the center after point p falls onto plane pl.

We first want to make sure that the point will actually hit the plane. To do this, we just want to check if the point lies in the plane's shadow, which is a triangle. We will use the pointLoc() method in Chapter 1. Let int a1 = pointLoc(pl[0], pl[1], p);, int a2 = pointLoc(pl[1], pl[2], p);, int a3 = pointLoc(pl[2], pl[0], p); If a1, a2, and a3 are $$$\geq 0$$$, or $$$\leq 0$$$, then the point lies in the triangle. You can convince yourself this is true, since if you start from pl[0], and drive around the triangle and then come back, if the point is inside the triangle, it will always be on your right side, or it will always be on your left side if you drive the other direction.

Now, we want to know where the point will hit the plane. We search "equation of a plane through 3 points" on google, and find this very helpful math stackexchange post: https://math.stackexchange.com/questions/2686606/equation-of-a-plane-passing-through-3-points. If we find the cross product: point coefs = cross(sub(pl[1], pl[0]), sub(pl[2], pl[0])), the plane can be represented as $$$coefs_x \cdot x + coefs_y \cdot y + coefs_z \cdot z = cst$$$, where $$$cst$$$ is a constant. We can find it, since pl[0] must lie on the plane, so $$$cst = coefs_x \cdot pl[0]_x + coefs_y \cdot pl[0]_y + coefs_z \cdot pl[0]_z$$$. When the point falls, its x coordinate and y coordinate stay the same as the start, and the z coordinate is what we want to find. Once the point falls, it will be a point on the plane, so it will satisfy the equation. $$$coefs_x \cdot p_x + coefs_y \cdot p_y + coefs_z \cdot z = cst$$$. We can use algebra and rearrange and find $$$z = \frac{cst - coefs_x \cdot p_x + coefs_y \cdot p_y}{coefs_z}$$$. Then, we can find the z-coordinate of the point. We can then find the center's z coordinate, $$$center_z - p_z + z$$$. If $$$coefs_z = 0$$$, the plane is vertical, it won't matter, since one of the point's adjacent edges will hit the plane's top edge instead.

The plane falling onto a point case is very similar.

Define: double planepoint(point center, point p, vector<point> pl). This method returns the z coordinate of the center after the plane pl falls onto point p.

We will find the point's z-coordinate in the plane's first coordinate system using the methods described above $$$(z)$$$, and we can find how much higher the center is compared to it $$$(center_z - z)$$$, and add that. Then, we can translate everything back to the point's coordinate system, and the center's z coordinate is $$$center_z - z + p_z$$$

We will be doing the line falls onto another line case.

double lineline(point center, pair<point, point> l1, pair<point,point> l2) This method returns the z coordinate of the center after line l1 falls onto l2.

We need to make sure the lines actually will fall on top of each other, so we need to make sure their shadows will intersect. If a line segment AB intersects another line segment CD in 2 dimensions, then if C is on one side of AB, D must be on the other. Similarly for CD, if A is on one side of CD, B must be on the other. It turns out that this is a necessary and sufficient condition, and we can use the pointLoc method we used earlier. So, our condition is (pointLoc(l1.first, l1.second, l2.first) != pointLoc(l1.first, l1.second, l2.second) && pointLoc(l2.first, l2.second, l1.first) != pointLoc(l2.first, l2.second, l1.second)).

We must now find the coordinates of where the 2 line segments will intersect. To do this, we can refer to Wikipedia: https://en.wikipedia.org/wiki/Line%E2%80%93line_intersection#Given_two_points_on_each_line, and yoink the coordinates. You can let x1 = l1.first.x, y1 = l1.first.y, x2 = l1.second.x, y2 = l1.second.y, x3 = l2.first.x, y3 = l2.first.y, x4 = l2.second.x, y4 = l2.second.y, then the x and y coordinates of the intersection are $$$x = \frac{(x_1 y_2 - y_1 x_2)(x_3 - x_4) - (x_1 - x_2)(x_3 y_4 - y_3 x_4)}{(x_1 - x_2)(y_3 - y_4) - (y_1 - y_2)(x_3 - x_4)}$$$, $$$y = \frac{(x_1 y_2 - y_1 x_2)(y_3 - y_4) - (y_1 - y_2)(x_3 y_4 - y_3 x_4)}{(x_1 - x_2)(y_3 - y_4) - (y_1 - y_2)(x_3 - x_4)}$$$. Then, we can use the dist2D method defined in Chapter 1 to find the intersection's z-coordinate.

We can find how far the intersection is along the second line to find its z coordinate using ratio = dist2D(l2.first, inter)/dist2D(l2.first, l2.second), and its z coordinate is inter.z = l2.first.z*(1-ratio) + l2.second.z*ratio. And we can find the z coordinate of the intersection is along the first line along the first coordinate system using dist2D and a similar technique (let ratio2 = dist2D(l1.first, inter)/dist2D(l1.first, l2.second), inter2.z = l1.first.z*(1-ratio2) + l1.second.z*ratio2), and our answer would be (inter.z + center.z - inter2.z).

We will now be doing the point-sphere/sphere-point case. The sc is the sphere's center.

double pointsphere(point center, point p, point sc, double radius)

double spherepoint(point center, point sc, double radius, point p){

To solve the point-sphere case, you can find the z-coordinate of the point once it is on the sphere using the pythagorean theorem: $$$(p_x - sc_x)^2 + (p_y - sc_y)^2 + (p_z - sc_z)^2 = radius^2$$$. Since the point is falling from the top, it will land at $$$z = \sqrt{radius^2 - (p_x - sc_x)^2 + (p_y - sc_y)^2} + sc_z$$$, so answer for that case is z + center.z - p.z. For the sphere-point case, the sphere center would be at $$$z = p_z + \sqrt{radius^2 - (p_x - sc_x)^2 + (p_y - sc_y)^2}$$$, and thus, the answer would be z + (center.z - sc.z)

This is probably the hardest case — a line falls onto a sphere, or a sphere falls onto a line.

double linesphere(point center, pair <point,point> l, point sc, double radius)

double sphereline(point center, point sc, double radius, pair <point,point> l)

We will take care of the line falls on sphere case first. The sphere falls on line case is very similar.

There are no good articles/resources about this online, so I had to improvise a solution myself.

Let l.first $$$ = p_1$$$, l.second $$$ = p_2$$$

Note that a point on the line can be represented using a parametric equation: $$$(x, y, z) = (p_{1x} + (p_{2x} - p_{1x})t, p_{1y} + (p_{2y}-p_{1y}) t , p_{1z} + (p_{2z}-p_{1z})t)$$$, with $$$0 \leq t \leq 1$$$

The answer when a point $$$(x, y, z)$$$ falls on to the sphere is $$$\sqrt{r^2 - (x - sc_x)^2 - (y - sc_y)^2} + sc_z + center_z - z$$$. This can be expressed in terms of t as $$$f(t) = \sqrt{r^2 - (p_{1x} + (p_{2x}-p_{1x})t - sc_x)^2 - (p_{1y} + (p_{2y}-p_{1y})t - sc_y)^2} + (sc_z + center_z) - (p_{1z} + (p_{2z} - p_{1z})t)$$$. Simplifying, we will get an expression where we have constants $$$a, b, c, d, e, f$$$, such that $$$f(t) = \sqrt{r^2 - (a+bt)^2 - (c+dt)^2} + (e+ft)$$$, and this can be further simplified to for constants $$$g, h, i$$$, $$$f(t) = \sqrt{g + ht + it^2} + (e+ft)$$$. We want to find the maximum $$$f(t)$$$, and to do that, we can take the derivative and setting it to 0, $$$f'(t) = \frac{2it + h}{2\sqrt{g+ht+it^2}} + f = 0$$$. Simplifying, you get $$$4i^2t^2 + 4iht + h^2 = 4f^2 (g + ht + it^2)$$$, and $$$(4i^2 - 4f^2g)t^2 + (4ih - 4f^2h)t + (h^2 - 4f^2i) = 0$$$ which is a quadratic in terms of $$$t$$$, which we can solve. Once we have the solution in $$$t$$$, and make sure it is between 0 and 1, we know that that t corresponds to a point, and that point is destined to be the first one to hit the sphere, so it will reduce to the point-sphere case. The sphere-line case is very similar, and after some math will also reduce to the sphere-point case.

double planesphere(point center, vector <point> pl, point sc, double radius) double sphereplane(point center, point sc, double radius, vector <point> pl)

The key observation here is that once the sphere touches the plane, it's center will be normal to the plane. Therefore, imagine a normal line to the plane running through the sphere's center. The two points that the line hits are the only 2 points you need to check. In code, we will find the normal vector $$$n$$$, and just consider the points $$$p_1, p_2 = sc \pm \lvert \frac{radius}{\lvert n \rvert} \rvert \cdot n$$$, where $$$\lvert n \rvert$$$ is the magnitude of n. This boils down to the plane-point case. We only need to consider the $$$p$$$ with higher z coordinate, and call planepoint(center, p, pl). Similar logic can be applied to the sphere-plane case.

The sphere-sphere case is very easy as well. double spheresphere(point center, point sc1, double rad1, point sc2, double rad2)

You should be able to convince yourself that instead of dropping a sphere onto another sphere, it is equivalent to drop a point onto the other sphere, except the other sphere has $$$rad2 = (rad2 + rad1)$$$. Essentially, you can call double pointsphere(center, sc1, sc2, rad1+rad2) and it should return the same answer.

https://en.wikipedia.org/wiki/Euler_angles is a very good resource.

It has a gif showing how z-x-z rotations work, the same as in the problem, and even better, it has a rotation matrix corresponding to z-x-z rotations which we can use to rotate points.

It is a good idea to create a rotation struct, which contains 3 angles, $$$\theta_1, \theta_2, \theta_3$$$. If $$$s_1 = sin(\theta_1)$$$, $$$c_2 = cos(\theta_1)$$$, and $$$s_2, c_2, s_3, c_3$$$ are defined similarly, the rotation matrix would be

Then, you can perform the rotation on a point $$$P$$$ by multiplying $$$M \cdot P$$$ to get a new point $$$P'$$$. This does the rotation about the origin.

Create a cloud struct. It has a

• int t (type),
• rotation r (rotation),
• a point c (center),
• a double s (size).

We also need to store the vertices, edges, faces, and spheres of each cloud. (if the cloud is a cube/tetrahedron)

We need to store how many of each simple object there are. (vertices, edges, faces, and spheres are all considerd simple).

A method should be made to initialize these.

A cube of size 1 has vertices coords $$$(\pm 0.5, \pm 0.5, \pm 0.5)$$$ ...

You should scale it up, and then translate to the correct x/y coordinates.

It has 12 edges, and it has 6 faces, but since our faces are triangular, we can split each of it's faces into 2 triangles. You can create this manually.

A tetrahedron of size 1 has vertices coords $$$(-\frac{0.5}{\sqrt{3}}, 0.5, -\frac{0.5}{\sqrt{6}}), (-\frac{0.5}{\sqrt{3}}, -0.5, -\frac{0.5}{\sqrt{6}}), (\frac{1}{\sqrt{3}}, 0, -\frac{0.5}{\sqrt{6}}), (0, 0, \frac{1.5}{\sqrt{6}})$$$ (https://www.qfbox.info/4d/tetrahedron)

It has 6 edges and it has 4 faces, these are simple to construct.

Then, you can create an intersection(cloud c1, cloud c2) method, and it returns the center of when cloud c1 fall onto cloud c2. To do so, we can split c1 and c2 into the shapes, and check these collisions:

• plane — point
• line — line
• point — plane
• sphere — point
• point — sphere
• sphere — line
• line — sphere
• sphere — plane
• plane — sphere
• sphere — sphere

Then, we can find how high the highest collision is.

Then, you can create a "ground" cloud which just consists of 2 triangular planes $$$(-20000, -20000, 0), (-20000, 20000, 0), (20000, -20000, 0)$$$, and $$$(-20000, 20000, 0), (20000, -20000, 0), (20000, 20000, 0)$$$.

The rest is very simple: it is a matter of initializing a cloud , then finding the highest point it intersects with all the fallen clouds. Then, we create a new cloud with that point. You can then find the highest point with

Clouds intersect when even a single point of the cloud touches another cloud. Therefore, your program may not count it when a cloud "skims" another cloud. Therefore, you should add some tolerances to some methods:

point-plane and plane-point: Make sure to have tolerances when the plane's z coord is < 1e-9 — the plane is vertical and you don't count it.

line-line: If the denominator is less than 1e-9, you can return -1.0.

sphere-point and point-sphere: Make sure to count intersections where the point "skims" the sphere, just barely misses the sphere.

sphere-line and line-sphere: First make the quadratic monic. (If the quadratic term is within 1e-9, treat it as a linear). If the discriminant is within 1e-9 of 0 (maybe negative), we can treat it as 0.

For pointLoc, make sure to have tolerances over the zero case as well. Check if that quantity is very close to 0 (within 1e-9 * dist(A, B) * dist(A, C)).

With some luck, you will get to test 33. That case is a bunch of cubes, and it can make your code TLE, since cubes have alot of line segments/faces/components in general. To optimize, note that if 2 objects are far away (> s1 + s2) in the x/y direction, there is no possibility of intersection, and you don't need to check them. Also, you can do the objects in a order from youngest (generally highest) to oldest (generally lowest), and keep track of the answer, and if it already intersects, if (ans > s1 + s2 + z), you can also skip it, since there is no possibility of it intersecting with an object that low. Be careful about this optimization with the ground — make it a special case where you don't use this optimization.