Problem A. Farthest Point · 數據結構與算法/leetcode/lintcode題解

# Problem A. Farthest Point（圓周上最遠整點） ### Source - [hihoCoder](http://hihocoder.com/contest/mstest2015sept2/problem/1) ### Problem 時間限制:5000ms 單點時限:1000ms 內存限制:256MB ### 描述 Given a circle on a two-dimentional plane. Output the **integral** point in or on the boundary of the circle which hasthe largest distance from the center. ### 輸入 One line with three floats which are all accurate to three decimal places,indicating the coordinates of the center x, y and the radius r. For 80% of the data: |x|,|y|<=1000, 1<=r<=1000 For 100% of the data: |x|,|y|<=100000, 1<=r<=100000 ### 輸出 One line with two integers separated by one space, indicating the answer. If there are multiple answers, print the one with the largest x-coordinate. If there are still multiple answers, print the one with the largesty-coordinate. #### 樣例輸入 ~~~ 1.000 1.000 5.000 ~~~ #### 樣例輸出 ~~~ 6 1 ~~~ ### 題解1 - 圓周枚舉其實自己最開始做這道題時用的就是枚舉，但是似乎忘記加圓心坐標了，一直 WA... 題目要求是返回最大的 x, 所以我們首先根據半徑范圍將 x 的整數解范圍求出來。然后求出可能的 y, 由于題中給出的解有3位小數，如果要精確求解的話，可以將圓方程兩邊同乘1000，然后判斷是否為整數。 ### Java ~~~ import java.io.*; import java.util.*; import java.util.Queue; class Point { long x; long y; Point(long x, long y) { this.x = x; this.y = y; } } public class Main { public static void main(String[] args) { Scanner in = new Scanner(System.in); double xd = in.nextDouble(), yd = in.nextDouble(), rd = in.nextDouble(); Point result = solve(xd, yd, rd); System.out.println(result.x + " " + result.y); } private static Point solve(double x0, double y0, double r) { // convert double to long(accurate) long xl0 = (long)(x0 * 1000), yl0 = (long)(y0 * 1000), rl0 = (long)(r * 1000); Point res = new Point(Long.MIN_VALUE, Long.MIN_VALUE); int lower_x = (int)Math.ceil(x0 - r), upper_x = (int)Math.floor(x0 + r); for (int i = upper_x; i >= lower_x; i--) { // circle above long y1l = yl0 + (long)(Math.sqrt(rl0*rl0 - (i*1000 - xl0)*(i*1000 - xl0)) + 0.5); if ((i*1000 - xl0)*(i*1000 - xl0) + (y1l - yl0)*(y1l - yl0) == rl0*rl0) { // ensure y1 is integer if (y1l % 1000 == 0) { res.x = i; res.y = y1l / 1000; return res; } } // circle below y1l = yl0 - (long)(Math.sqrt(rl0*rl0 - (i*1000 - xl0)*(i*1000 - xl0)) + 0.5); if ((i*1000 - xl0)*(i*1000 - xl0) + (y1l - yl0)*(y1l - yl0) == rl0*rl0) { // ensure y1 is integer if (y1l % 1000 == 0) { res.x = i; res.y = y1l / 1000; return res; } } } return res; } } ~~~ ### 源碼分析自右向左枚舉，先枚舉圓的上半部分，再枚舉圓的下半部分。注意1000的轉換。 ### 復雜度分析最壞情況下 O(R)O(R)O(R). ### 題解2 - 整數分解看似容易實則比較難的一道題，現場通過率非常低。我們仔細審下題，求圓周上的整點，有多個整點時輸出最大的 x 和最大的 y. 容易想到的方案是枚舉所有可能的 x 和 y, 然后代入等式測試是否相等，這個過不了大的 x 和 y. 如果用開方的方法必然有誤差，我用這種方法不知道貢獻了多少 WA, 淚流滿面... 作為在線測試，**更為合理的方案應為先暴力搜索拿到百分之八十的分數。** 從 Microsoft 和 Google APAC 在線測試的風格來看是偏向于程序設計競賽的，那么題目的考點自然就在競賽范圍之內，這道題看似是浮點型的數據，~~實際上考的卻是整數中數論的基礎。~~**注意題中的 accurate to three decimal places, 那么也就意味著我們對給定的數據同乘 10310^3103 后一定是整數！！**！這個關鍵的信息我在測試過程中也沒注意到，直到第二天早上醒來后突然就想到了！興奮地六點多就爬起來了。首先肯定是要寫出圓方程的，設圓心坐標為 (x0,y0)(x_0, y_0)(x0,y0), 半徑為 rrr, 那么我們有：(x?x0)2+(y?y0)2=r2(x - x_0)^2 + (y - y_0)^2 = r^2(x?x0)2+(y?y0)2=r2 設 m=103(x?x0)m = 10^3(x - x_0)m=103(x?x0), n=103(y?y0)n = 10^3(y - y_0)n=103(y?y0), R=103rR = 10^3rR=103r, 那么我們有新的圓方程：m2+n2=R2m^2 + n^2 = R^2m2+n2=R2其中 `m, n, R` 均為整數。接下來我們看看給出的數據范圍，x, y, r 均是 10610^6106 以內，那么圓方程兩邊同乘 10610^6106 （括號內的數即乘上 10310^3103）后數據在 101810^{18}1018 以內。我們來估算下整數的范圍，210≈1032^{10} \approx 10^3210≈103, Java 中 int 型為4個字節，最大為 231?1≈2?1092^{31} - 1 \approx 2 \cdot 10^9231?1≈2?109, long 型為8個字節，最大為 263?1≈23?10182^{63} - 1 \approx 2^3 \cdot 10^{18}263?1≈23?1018, 估算下來應該選用 long 保存 m, n, R. 接下來就是數論部分的推導了，先來一個簡單的推導，勾股數部分的推導不直觀。首先從已知部分出發，已知的只有勾股數方程和 m, n 均是整數，那么接下來肯定是要利用整數的理論無疑了。我們首先對以上圓方程移項開方，考慮到圓的對稱性，我們其實只需要考慮圓的八分之一即可。這里考慮`0 < m < r`部分，`m == 0`表示在點在軸上，最后單獨加上。 m=R2?n2=(R+n)(R?n)m = \sqrt{R^2 - n^2} = \sqrt{(R + n)(R - n)}m=√R2?n2=√(R+n)(R?n)由于 m 一定是整數，故根號內一定為完全平方數，由于排除了軸上的點，那么`-R < n < R`, 設`G = gcd(R + n, R - n)`, p2=(R+n)/Gp^2 = (R + n) / Gp2=(R+n)/G, q2=(R?n)/Gq^2 = (R - n) / Gq2=(R?n)/G, 于是我們有`m = Gpq`, `p > q`, 由于`G` 是`R + n` 和`R - n` 的最大公約數，故`p` 和`q`一定互質，且有：p2+q2=2R/Gp^2 + q^2 = 2R / Gp2+q2=2R/G由于`p`,`q` 都大于等于1，那么我們能推斷出`G` 一定是 `2R` 的約數！根據約數(素數)部分的基礎理論，我們可以在 O(2R)O(\sqrt{2R})O(√2R) 時間內找出所有約數。然后對以上等式進行縮放得到`p` 的范圍，枚舉求解，判斷`p^2` 和`q^2` 是否互質(最大公約數是否為1)。 ### Java ~~~ import java.io.*; import java.util.*; class Point { long x; long y; Point(long x, long y) { this.x = x; this.y = y; } } public class Main { public static void main(String[] args) { Scanner in = new Scanner(System.in); double xd = in.nextDouble(), yd = in.nextDouble(), rd = in.nextDouble(); // convert double to long(accurate) long x0 = (long)(xd * 1000), y0 = (long)(yd * 1000), r0 = (long)(rd * 1000); Point result = solve(x0, y0, r0); System.out.println(result.x + " " + result.y); } private static Point solve(long x0, long y0, long r) { Point res = new Point(Long.MIN_VALUE, Long.MIN_VALUE); long x_max = Long.MIN_VALUE, y_max = Long.MIN_VALUE; // p^2 + q^2 = 2R/divisor, p > q >= 1 // 1 <= q^2 < R/G < p^2 < 2R/G ==> p >= 2 List<Long> divisors = getDivisor(r << 1); for (long divisor : divisors) { long lower = Math.max(2, (long)Math.sqrt(r * 1.0/ divisor)); // long upper = (long)Math.sqrt(2.0 * r / divisor); for (long p = lower; p * p <= 2 * r / divisor; p++) { long q = (long)Math.sqrt(2.0 * r / divisor - p * p); // check if q is integer if (p * p + q * q != 2 * r / divisor) continue; // ensure p^2 and q^2 have no common divisor if (gcd(p * p, q * q) == 1) { long m = divisor * p * q; long n = r - p * p * divisor; List<Point> points = new ArrayList<Point>(); points.add(new Point(m + x0, n + y0)); points.add(new Point(m + x0, -1 * n + y0)); points.add(new Point(-1 * m + x0, n + y0)); points.add(new Point(-1 * m + x0, -1 * n + y0)); for (Point point : points) { updateAns(point, res); } } } } // axis point check List<Point> axis = new ArrayList<Point>(); axis.add(new Point(x0 + r, y0)); axis.add(new Point(x0 - r, y0)); axis.add(new Point(x0, y0 + r)); axis.add(new Point(x0, y0 - r)); for (Point point : axis) { updateAns(point, res); } // divide by 1000 res.x /= 1000; res.y /= 1000; return res; } public static void updateAns(Point p, Point res) { // point(x, y) in integer if ((p.x % 1000 == 0) && (p.y % 1000 == 0)) { if (p.x > res.x) { res.x = p.x; res.y = p.y; } else if (p.x == res.x && p.y > res.y) { res.y = p.y; } } } // enumerate all the divisor for n public static List<Long> getDivisor(long n) { List<Long> result = new ArrayList<Long>(); for (long i = 1; i * i <= n; i++) { if (n % i == 0) { result.add(i); // i * i <= n ==> i <= n / i if (i != n / i) result.add(n / i); } } Collections.sort(result); return result; } public static long gcd(long a, long b) { return (b == 0L) ? a : gcd(b, a % b); } } ~~~ ### 源碼分析由于更新結果的操作非常頻繁，單獨寫一個方法較好。 ### 復雜度分析求所有素數時間復雜度 O(n)O(\sqrt{n})O(√n), 判斷是否互質時間復雜度 O(logn)O(\log n)O(logn). 枚舉最大公約數時間復雜度約 (n)(\sqrt{n})(√n)，總的時間復雜度估算應該比 O(n)O(n)O(n) 小一些，但是小的不明顯。**所以說，這種方法費了老大勁，但是吃力不討好！筆試中這種方法極不可取！** ### 題解3 - 勾股數除了以上使用數論部分整數分解的方法外，還可以巧用勾股數的特性，這種方法需要熟知勾股數的特性。設正整數 m,n,rm, n, rm,n,r 滿足：m2+n2=r2m^2 + n^2 = r^2m2+n2=r2我們對上式兩邊進行平方可得：(m2?n2)2+(2mn)2=(m2+n2)2=(r2)2(m^2 - n^2)^2 + (2mn)^2 = (m^2 + n^2)^2 = (r^2)^2(m2?n2)2+(2mn)2=(m2+n2)2=(r2)2令 a=m2?n2a = m^2 - n^2a=m2?n2, b=2mnb = 2mnb=2mn, c=m2+n2c = m^2 + n^2c=m2+n2. 容易得到：a2+b2=c2a^2 + b^2 = c^2a2+b2=c2注意到上述推導可逆，那么也就是說只要我們找到正整數滿足`m > n`就能找到所有可能的勾股數。且根據素勾股數的特性，`m`, `n` 為一奇一偶，不妨設其為`2k-1`, `2k`. 代入`c`中可知`c`為`4K + 1`. 即`c % 4 = 1`. 根據 [Tree of primitive Pythagorean triples](https://en.wikipedia.org/wiki/Tree_of_primitive_Pythagorean_triples) 中提到的方法，我們只需找出小于給定的`r`的素勾股數即可，然后判斷是否能整除`r`. ### Java ~~~ import java.io.*; import java.util.*; import java.util.Queue; class Point { long x; long y; Point(long x, long y) { this.x = x; this.y = y; } } class Pythagorean { long x; long y; long z; Pythagorean(long x, long y, long z) { this.x = x; this.y = y; this.z = z; } } public class Main { public static void main(String[] args) { Scanner in = new Scanner(System.in); double xd = in.nextDouble(), yd = in.nextDouble(), rd = in.nextDouble(); // convert double to long(accurate) long x0 = (long)(xd * 1000), y0 = (long)(yd * 1000), r0 = (long)(rd * 1000); Point result = solve(x0, y0, r0); System.out.println(result.x + " " + result.y); } private static Point solve(long x0, long y0, long r) { Point res = new Point(Long.MIN_VALUE, Long.MIN_VALUE); long x_max = Long.MIN_VALUE, y_max = Long.MIN_VALUE; // init Pythagorean pyth0 = new Pythagorean(3, 4, 5); Queue<Pythagorean> q = new LinkedList<Pythagorean>(); q.offer(pyth0); boolean update = true; while (!q.isEmpty()) { int qSize = q.size(); for (int i = 0; i < qSize; i++) { Pythagorean pyth = q.poll(); if ((r % pyth.z) == 0) { // System.out.println("x = " + pyth.x + ", y = " + pyth.y + ", r = " + pyth.z); long k = r / pyth.z; long kx = k * pyth.x; long ky = k * pyth.y; List<Point> points = new ArrayList<Point>(); points.add(new Point(x0 + kx, y0 + ky)); points.add(new Point(x0 - kx, y0 + ky)); points.add(new Point(x0 + kx, y0 - ky)); points.add(new Point(x0 - kx, y0 - ky)); if (kx != ky) { points.add(new Point(y0 + ky, x0 + kx)); points.add(new Point(y0 - ky, x0 + kx)); points.add(new Point(y0 + ky, x0 - kx)); points.add(new Point(y0 - ky, x0 - kx)); } for (Point point : points) { updateAns(point, res); } } // add next level Pythagorean for (Pythagorean p : nextPyths(pyth)) { if (p.z > r) continue; q.offer(p); } } } // axis point check List<Point> axis = new ArrayList<Point>(); axis.add(new Point(x0 + r, y0)); axis.add(new Point(x0 - r, y0)); axis.add(new Point(x0, y0 + r)); axis.add(new Point(x0, y0 - r)); for (Point point : axis) { updateAns(point, res); } // divide by 1000 res.x /= 1000; res.y /= 1000; return res; } public static List<Pythagorean> nextPyths(Pythagorean pyth) { List<Pythagorean> pyths = new ArrayList<Pythagorean>(); // method 1 Pythagorean pyth1 = new Pythagorean(0, 0, 0); pyth1.x = pyth.x - 2 * pyth.y + 2 * pyth.z; pyth1.y = 2 * pyth.x - 1 * pyth.y + 2 * pyth.z; pyth1.z = 2 * pyth.x - 2 * pyth.y + 3 * pyth.z; pyths.add(pyth1); // method 2 Pythagorean pyth2 = new Pythagorean(0, 0, 0); pyth2.x = pyth.x + 2 * pyth.y + 2 * pyth.z; pyth2.y = 2 * pyth.x + 1 * pyth.y + 2 * pyth.z; pyth2.z = 2 * pyth.x + 2 * pyth.y + 3 * pyth.z; pyths.add(pyth2); // method 3 Pythagorean pyth3 = new Pythagorean(0, 0, 0); pyth3.x = -1 * pyth.x + 2 * pyth.y + 2 * pyth.z; pyth3.y = -2 * pyth.x + pyth.y + 2 * pyth.z; pyth3.z = -2 * pyth.x + 2 * pyth.y + 3 * pyth.z; pyths.add(pyth3); return pyths; } public static void updateAns(Point p, Point res) { // point(x, y) in integer if ((p.x % 1000 == 0) && (p.y % 1000 == 0)) { if (p.x > res.x) { res.x = p.x; res.y = p.y; } else if (p.x == res.x && p.y > res.y) { res.y = p.y; } } } } ~~~ ### 源碼分析根據鏈接中提到的數據結構，使用隊列按層次遍歷較好，但是空間消耗較大，所以在入隊時一定要剪枝。 ### 復雜度分析時間復雜度最壞情況下需要遍歷所有可能素勾股數。空間復雜度消耗也比較客觀... ### Reference - [BZOJ 1041 [HAOI2008] 圓上的整點題解與分析 - 初學者 - 博客頻道 - CSDN.NET](http://blog.csdn.net/csyzcyj/article/details/10044629) - [[BZOJ1041 [HAOI2008]圓上的整點]數論、勾股數相關定理 | edward_mj](http://edward-mj.com/archives/166) - [勾股數 - 維基百科，自由的百科全書](https://zh.wikipedia.org/wiki/%E5%8B%BE%E8%82%A1%E6%95%B0) - [hihoCoder](http://hihocoder.com/discuss/question/2619)