数学研发论坛

 找回密码
 欢迎注册
楼主: KeyTo9_Fans

[原创] 谢尔宾斯基地毯与圆的交集

[复制链接]
发表于 2017-9-26 21:00:41 | 显示全部楼层
当考察整个小方格是否位于圆内/外或者被圆穿越时,有没有人担心两个对顶格所公共的对顶格点刚好在圆上,使得两格刚好一内一外?
或者像@lsr314 那样仅仅考察格子的中心时,有没有考虑过一个格子的中心刚好在圆上?

一个简单而有趣的事实是:无论分辨率多大,不可能有一个格子,它的中心或者一个顶点刚好在圆上。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2017-9-26 21:07:07 | 显示全部楼层
其实这里面有大量的重复计算。 当确定了某个正方形已经是 整体的在圆内,或者圆外时,下一次迭代的时,由此产生的8个小正方形跟圆的关系 没必要重复判断 。我们只需进行简单的计数就行。 这样程序只进行 圆上的方格个数的统计。这样能大大的降低空间复杂度和时间复杂度。

点评

^_^,以前是冲着画图去的,现在冲着计算点的个数了。见下:  发表于 2017-9-26 22:04
对对对,俺在8#不就是这么计划的么,只是俺不擅长编程实现  发表于 2017-9-26 21:21
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2017-9-26 22:01:42 | 显示全部楼层
wayne 发表于 2017-9-26 21:07
其实这里面有大量的重复计算。 当确定了某个正方形已经是 整体的在圆内,或者圆外时,下一次迭代的时,由此 ...


按次思路写了一个很笨的代码,迭代10次,得到:
下面的第一个三元组分别代表圆外,圆上,圆内的格子的个数。
  1. {{0,8,0},8,{0,1.0000000,0}}
  2. {{4,28,32},64,{0.062500000,0.43750000,0.50000000}}
  3. {{104,76,332},512,{0.20312500,0.14843750,0.64843750}}
  4. {{984,204,2908},4096,{0.24023438,0.049804688,0.70996094}}
  5. {{8288,580,23900},32768,{0.25292969,0.017700195,0.72937012}}
  6. {{67744,1556,192844},262144,{0.25842285,0.0059356689,0.73564148}}
  7. {{545904,4180,1547068},2097152,{0.26030731,0.0019931793,0.73769951}}
  8. {{4377944,11204,12388068},16777216,{0.26094580,0.00066781044,0.73838639}}
  9. {{35052688,29724,99135316},134217728,{0.26116288,0.00022146106,0.73861566}}
  10. {{280499768,79276,793162780},1073741824,{0.26123577,0.000073831528,0.73869040}}
  11. {{2244206376,212076,6345516140},8589934592,{0.26126001,0.000024688896,0.73871530}}
  12. {{17954209288,565692,50764701756},68719476736,{0.26126813,8.2319020*10^-6,0.73872364}}
复制代码
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2017-9-26 22:21:45 | 显示全部楼层
wayne 发表于 2017-9-26 22:01
按次思路写了一个很笨的代码,迭代10次,得到:
下面的第一个三元组分别代表圆外,圆上,圆内的格子的 ...

迭代$18$次,可得:
  1. {0,8,0},8,{0.0000000000000000,1.0000000000000000,0.0000000000000000}
  2. {4,28,32},64,{0.0625000000000000,0.4375000000000000,0.5000000000000000}
  3. {104,76,332},512,{0.2031250000000000,0.1484375000000000,0.6484375000000000}
  4. {984,204,2908},4096,{0.2402343750000000,0.0498046875000000,0.7099609375000000}
  5. {8288,580,23900},32768,{0.2529296875000000,0.0177001953125000,0.7293701171875000}
  6. {67744,1556,192844},262144,{0.2584228515625000,0.0059356689453125,0.7356414794921875}
  7. {545904,4180,1547068},2097152,{0.2603073120117188,0.0019931793212891,0.7376995086669922}
  8. {4377944,11204,12388068},16777216,{0.2609457969665527,0.0006678104400635,0.7383863925933838}
  9. {35052688,29724,99135316},134217728,{0.2611628770828247,0.0002214610576630,0.7386156618595123}
  10. {280499768,79276,793162780},1073741824,{0.2612357661128044,0.0000738315284252,0.7386904023587704}
  11. {2244206376,212076,6345516140},8589934592,{0.2612600075080991,0.0000246888957918,0.7387153035961092}
  12. {17954209288,565692,50764701756},68719476736,{0.2612681315513328,0.0000082319020294,0.7387236365466379}
  13. {143635171632,1509332,406119132924},549755813888,{0.2612708551750984,0.0000027454589144,0.7387263993659872}
  14. {1149085383304,4026028,3248957101772},4398046511104,{0.2612717669999256,0.0000009154127838,0.7387273175872906}
  15. {9192693786392,10740796,25991667561644},35184372088832,{0.2612720716795138,0.0000003052717830,0.7387276230487032}
  16. {73541578891856,28646804,207933369171996},281474976710656,{0.2612721732896830,0.0000001017738924,0.7387277249364246}
  17. {588332707461496,76396620,1663467029827132},2251799813685248,{0.2612722071855238,0.0000000339269146,0.7387277588875616}
  18. {4706661863240488,203728972,13307736442512524},18014398509481984,{0.2612722184847365,0.0000000113092298,0.7387277702060338}
复制代码

点评

19次:{37653295449008288,543283204,106461892083564380},{0.26127222225314139514,3.7697845123307871518*10^-9,0.73872777397707409253}  发表于 2017-9-28 13:32
18次,跑多长时间了  发表于 2017-9-26 22:27
用C语言写的吧,^_^  发表于 2017-9-26 22:27
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2017-9-26 22:25:31 | 显示全部楼层
KeyTo9_Fans 发表于 2017-9-26 22:21
迭代$18$次,可得:

, 刚补上了 2个结果,结果 fans一下子给出了 18个,不好玩

点评

使用$64bit$的整数计算,只能出$18$个,第$19$个超出范围算不出来……  发表于 2017-9-26 22:28
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2017-9-26 22:44:59 | 显示全部楼层
hujunhua 发表于 2017-9-26 21:00
当考察整个小方格是否位于圆内/外或者被圆穿越时,有没有人担心两个对顶格所公共的对顶格点刚好在圆上,使 ...


是因为1/2 无法被表达成三进制里的有理数吧。
===
额,这种表达 好像有点问题,从迭代过程看, 正方形的顶点应该 都是从 三进制里的有理数。

点评

是因为3不是 费马素数(可表为两个整数的平方和者)  发表于 2017-9-26 23:10
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2017-9-26 23:06:18 | 显示全部楼层
设上述圆内、圆上的非白格数分别为 a(n),b(n),记@KeyTo9Fans要的面积为A,则\[
A=\lim_{n\to\infty}\frac{a(n)}{8^n},0=\lim_{n\to\infty}\frac{b(n)}{8^n}
\Longrightarrow A=\lim_{n\to\infty}\frac{a(n)+\lambda b(n)}{8^n}
\]`\lambda`为 0~1 的一个常系数。如果取\(\lambda_1\)序列递增而\(\lambda_2\)序列递减,这就得到了一个上界和下界。
显然取 0 和 1 可以得到一个上、下界,不过这未免也太懒了。貌似可取\(\lambda_1=1/2\),能再大一点么?。

点评

$b(n) ~ 4.3815 *(8/3)^n$  发表于 2017-9-27 15:20
不能。$1/2$已是最大。  发表于 2017-9-26 23:32
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2017-9-26 23:28:35 | 显示全部楼层
取\(\lambda=1/2\)的“面积”序列
{{0, 8, 0}, 0.50000000000000000000},
{{4, 28, 32},  0.71875000000000000000},
{{104, 76, 332},  0.72265625000000000000},
{{984, 204, 2908},  0.73486328125000000000},
{{8288, 580, 23900},  0.73822021484375000000},
{{67744, 1556, 192844},  0.73860931396484375000},
{{545904, 4180, 1547068},  0.73869609832763671875},
{{4377944, 11204, 12388068},  0.73872029781341552734},
{{35052688, 29724, 99135316},  0.73872639238834381104},
{{280499768, 79276, 793162780},  0.73872731812298297882},
{{2244206376, 212076, 6345516140}, 0.73872764804400503635},
{{17954209288, 565692, 50764701756},  0.73872775249765254557},
{{143635171632, 1509332, 406119132924},  0.73872777209544437937},
{{1149085383304, 4026028, 3248957101772},  0.73872777529368249816},
{{9192693786392, 10740796, 25991667561644},  0.73872777568459468966},
{{73541578891856, 28646804, 207933369171996}, 0.73872777582337079139},
{{588332707461496, 76396620, 1663467029827132}, 0.73872777585101889741},
{{4706661863240488, 203728972, 13307736442512524}, 0.73872777586064863886}

点评

据此可得极限0.738727775865795  发表于 2017-9-26 23:33
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2017-9-28 13:29:20 | 显示全部楼层
此图还存在一个y=x的对称性。于是,只需计算1/8的面积就行了。迭代至19次, 计算了大约半个小时,占用内存超过6G,这个内存占用主要来自于保存1/8的圆上的格子坐标。终于突破了Fans的记录,第19次,20次,……的迭代结果为:
  1. 19: {37653295449008288,543283204,106461892083564380},{0.26127222225314139514,3.7697845123307871518*10^-9,0.73872777397707409253}
  2. 20: {301226365040331144,1448779164,851695138117736668},{0.26127222350931090383,1.2566156136484263328*10^-9,0.73872777523407348252}
  3. 21: {2409810924185402784,3863345612,6813561108806027412},{0.26127222392811149179,4.1886477055927406887*10^-10,0.73872777565302373765}
  4. 22: {19278487403784147304,10302538780,54508488880751520380},{0.26127222406771505802,1.3962543659240198024*10^-10,0.73872777579265950539}
  5. 23: {154227899257743649824,27473690092,436067911073488311796}
  6. 24: {1233823194135208730288,73263231116,3488543288661173252292}
  7. 25: {9870585553277031899992,195369181668,27908346309484760627908}
  8. 26: {78964684426737227956608,520985280228,223266770476399080439708}
  9. 27: {631717475415287101141792,1389296277316,1786134163812581951993244}
  10. 28: {5053739803326001562892000,3704793953044,14289073310504360438453772}
  11. 29: {40429918426617891895941808,9879454623900,114312586484044763011824820}
  12. 30: {323439347412969480294205216,26345219429236,914500691872384449385489772}
复制代码

面积比例: $0.73872777586246736918$

算出第20次迭代的结果,果断终止程序的运行,因为此时虽然只是过去了一个小时12分钟,我有的是时间等待,但是内存已经占用了18G 了。
修改Fans的代码,运行1分34秒完成22次迭代,10分22秒完成24次迭代,1小时16分完成26次迭代,8小时27分完成28次迭代,2天12小时27分完成30次迭代,预计在2017年10月23日完成32次迭代。

点评

看到了,哈哈,谢谢。不用保存圆上的坐标倒是让我没想到,回头我拿过来,看能不能再往前推进!  发表于 2017-10-1 23:02
为保持讨论的连贯性,我没有发新贴,直接把我的代码附在$31$楼你发的代码后面了。  发表于 2017-10-1 22:56
这么强大!,把代码贡献出来呗,让我学习学习!  发表于 2017-10-1 22:20
冒着超出$64bit$整数范围的风险,把$24$楼的迭代次数直接改成$20$,硬着头皮重新运行程序,得到了相同的结果,运行时间为$14.79$秒。  发表于 2017-10-1 22:18

评分

参与人数 1威望 +3 金币 +3 贡献 +3 经验 +3 鲜花 +3 收起 理由
KeyTo9_Fans + 3 + 3 + 3 + 3 + 3 冒险硬着头皮迭代20次,得到了相同的结果

查看全部评分

毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2017-9-28 16:17:47 来自手机 | 显示全部楼层
我们可以查看圆弧和一个3*3网格相交的一般情况得出的函数。一个确定圆弧需要三个参数,圆心坐标和半径x,y,r.
我们把3*3正方形网格中间部分挖掉,然后看余下部分的圆弧下方部分占总面积比例为f0(x,y,r)
然后我们定义f_{n+1}(x,y,r)为8个小正方形各自看成一个大正方形计算对应fn值后的平均值。
那么n趋向无穷时结果就是我们所要的极限函数
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

小黑屋|手机版|数学研发网 ( 苏ICP备07505100号 )

GMT+8, 2019-10-18 15:20 , Processed in 0.063484 second(s), 16 queries .

Powered by Discuz! X3.4

© 2001-2017 Comsenz Inc.

快速回复 返回顶部 返回列表