-- 作者：tgs28
-- 发布时间：2006/10/2 9:14:00
-- [原创]计算某日星期、干支的统一方法及公历转换成农历的公式
计算某日星期、干支、农历日期的统一方法及公历转换成农历的公式 作者：谭笑风 为方便叙述，首先介绍两个函数
高斯函数：f(x)＝[x]，[x]表示不超过数x的最大整数，
那么对数x进行取整有三种方法:用高斯函数分别表示为
去尾法取整f(x)＝[x];进一法取整f(x)=-[-x];四舍五入法取整f(x)=[x+0.5].
求余数函数：f(x)=xmodt，表示x除以t后所得的余数（即mod表示余数运算）
且被除数x=除数t*商+余数，所以有x=t*[x/t]+xmodt 如果以某年的1月1日为这一年第1天，那么这一年2月1日就是第32天……
如果设月份数为m，日期数为d，求m月d日是这一年的第几天(即这一天的序数)，
这也容易算出，但能否给出个公式，把m,d代入公式直接求出序数g(m,d)呢？
下面来推导这个公式：
要求出序数g(m,d)只需把日期d加上这个月之前各月天数的总和g(m)即可，
比如说求1月份某日序数，把日期加0即可，
求2月份某日序数，把日期加31即可……
容易得出平年和闰年的1-12月份的月前天数总和g(m)分别是：
月份``m` ，1，`2， `3，`4， ``5， ``6， ``7， ``8， ``9， `10，`11， `12。
平年g(m)，0，31，59，90，120，151，181，212，243，273，304，334。
闰年g(m)，0，31，60，91，121，152，182，213，244，274，305，335。
下面来求这个月前天数总和g(m)的表达式,
假设(1)如果一年中的每个月都是30天，很容易得出月前天数总和g(m)为30m-30，
假设(2)如果一年中的单月都是31天，双月都是30天，
月前天数总和比假设(1)每个双月多一天，在假设(1)基础加[m/2],(即对m/2取整)
即可以得出这种情况下g(m)为30m+[m/2]-30，
假设(3)如果一年中8月以前的单月都是31天，双月都是30天，
并且8月及其以后的单月都是30天，双月都是31天，
月前天数总和比假设(2)中9月份及其以后每个月多一天，在假设(2)基础加[m/9]，
即可以得出这种情况下g(m)为30m+[m/2]+[m/9]-30，
假设(3)跟实际情况已经很接近了，只差2月份的天数不同，
先考虑平年时的情况，2月份是28天，比假设(3)中少两天，
所以当m=1，2时g(m)为30m+[m/2]+[m/9]-30，
当m≥3时，g(m)为30m+[m/2]+[m/9]-30-2=30m+[m/2]+[m/9]-32，
可以用式子[2/m]和[1/m]把m=1，2时的情况化到上式当中，即得
平g(m)=30m+[m/2]+[m/9]-32+([2/m]-[1/m])*2，由此还可得闰年时
闰g(m)=30m+[m/2]+[m/9]-31+([2/m]-[1/m])，
设年份数为y，所以这个g(m)的统一式还和y有关。
前面已定义f(x)=xmodt表示x除以t后所得的余数，则闰年的条件可以这样来表述：
已知当ymod4=0且ymod100≠0或ymod400=0时，y年才为闰年，
由开篇介绍可知对(ymodt)/t用进一法取整表达式为-[-(ymodt)/t]，
其含义为当ymodt=0时，-[-(ymodt)/t]=0；当ymodt≠0时，-[-(ymodt)/t]=1，
再由开篇x=t*[x/t]+xmodt可得-[-(ymodt)/t]=-[[y/t]-y/t]，
令m(y)=-[[y/4]-y/4]+[[y/100]-y/100]-[[y/400]-y/400]，
那么可以得出当y年是闰年时，m(y)=0，当y年不是闰年时，m(y)=1，
以上几式都可以用讨论的方法加以证明，证明从略。
将m(y)应用到上面的g(m)式中，可以得到g(m)的统一式：
g(m)=30m+[m/2]+[m/9]-31-m(y)+([2/m]-[1/m])*(1+m(y))，
所以g(m,d)=g(m)+d=30m+[m/2]+[m/9]-31-m(y)+([2/m]-[1/m])*(1+m(y))+d。 假设从公元前一年到现在的历法都是现行历法(格里高利历法)，
即公历置闰方法未变，都是4年一闰且100年不闰400年又闰，
实际没有公元0年，现在假设有公元0年，并把公元0年1月1日设为第1天（序数为1），
那么公元1年1月1日就是第367天（序数为367）……
但能否给出个公式g(y,m,d)，把y,m,d代入公式直接求y年m月d日这一天的序数呢？
下面来推导这个公式：
要求出序数g(y,m,d)只需把g(m,d)加上这一年之前的天数总和g(y)，
假设(1)如果每年都是365天，很容易得出某一年之前的天数总和g(y)为365y，
假设(2)如果每4年一闰，即每4k+1的年份的序数都比假设(1)多一天，
由此可以得出g(y)为365y+[(y-1)/4]+1=365y+[y/4]-[[y/4]-y/4]
[(y-1)/4]+1=[y/4]-[[y/4]-y/4]可以用讨论的方法加以证明
由此可以得到g(y)的表达式为
g(y)=365y+[y/4]-[[y/4]-y/4]-[y/100]+[[y/100]-y/100]+[y/400]-[[y/400]-y/400]，
上面有m(y)=-[[y/4]-y/4]+[[y/100]-y/100]-[[y/400]-y/400]，
所以g(y)=365y+[y/4]-[y/100]+[y/400]+m(y)，
g(y,m,d)=365y+[y/4]-[y/100]+[y/400]+m(y)+g(m)+d，
g(y,m,d)=365y+[y/4]-[y/100]+[y/400]+30m+[m/2]+[m/9]-31+([2/m]-[1/m])*(1+m(y))+d。
根据这个式子很容易求出任意两天之间相差的天数是多少，
即分别求出两天的序数之后作差即可。 为方便计算，下面来化简这个式子
设公元y年m月d日（y≥1582且y+m/10≥1583且y+m+d/15≥1593，1≤m≤12）
设公元y年年份的末两位数为b（b=y-[y/100]），
年份去掉末两位数后剩下的数为a（a=[y/100]），即有y=100a+b。
设f(y)=365y+[y/4]-[y/100]+[y/400]，将y=100a+b代入式中得
f(y)=36500a+365b+[25a+b/4]-[a+b/100]+[a/4+b/400]
b为年份的末两位数，所以0≤b＜100，所以[b/100]=0，b/400＜1/4，
用讨论的方法(设a=4k+0，1，2，3)可以证明[a/4+b/400]=[a/4],所以
f(y)=36500a+365b+25a+[b/4]-a+[a/4]，
f(y)=36524a+[a/4]+365b+[b/4]，所以
g(y,m,d)=36524a+[a/4]+365b+[b/4]+30m+[m/2]+[m/9]-31+([2/m]-[1/m])*(1+m(y))+d。 设f(m)=30m+[m/2]+[m/9]+([2/m]-[1/m])*(1+m(y))+ε,(ε为一待定常数项)
设f(y,m,d)=f(y)+f(m)+d，则
f(y,m,d)=36524a+[a/4]+365b+[b/4]+30m+[m/2]+[m/9]+([2/m]-[1/m])*(1+m(y))+ε+d，
则f(y,m,d)就是公元y年m月d日一个相对的序数，
如果以某个历法周期p为模(除数)，对进行求余数运算就可以求出某一天的是这个周期的哪一天了，
(1)如果以7为模，就可以求某一天星期几，
f(y,m,d)≡36524a+[a/4]+365b+[b/4]+f(m)+d (mod7)，
≡5a+[a/4]+b+[b/4]+2m+[m/2]+[m/9]+([2/m]-[1/m])*(1+m(y))+ε+d (mod7)，
然后对照历书就可以得出ε的值，
在这里f(m)是以公式的形式给出的，如果能直接记住每个月f(m)(mod7)的值，
会使计算更加快捷的，总之记的数据越多，算起来越快。
(2)如果以60为模，就可以求公历某年某月某日的干支序数，
f(y,m,d)≡36524a+[a/4]+365b+[b/4]+f(m)+d (mod60)，
f(y,m,d)≡44a+[a/4]+5b+[b/4]+f(m)+d (mod60)，
f(y,m,d)≡4*((11*(amod15))mod15)+[a/4]+5*bmod12+[b/4]+f(m)+d (mod60)，
1-12月对应的f(m)的余数分别对应10，41，9，40，10，41，11，42，13，43，14，44；
闰年的1，2月份要减去1，最后求出的结果1-59分别对应干支甲子-壬戌，0对应癸亥，
即求出的结果用10除求余数，余数为0-9时分别对应天干：癸甲乙丙丁戊己庚辛壬，
求出的结果用12除求余数，余数为0-11时分别对应地支：亥子丑寅卯辰巳午未申酉戌亥，
当a=19时，44a+[a/4]≡0(mod60)；当a=20时，44a+[a/4]≡45(mod60)；
记住世纪的余数0和45就能很快算出上下100年某天的干支了，
如1949年10月1日干支序数为0+5*49mod12+[49/4]+43+1≡5+12+43+1≡61≡1(mod60)，
所以新中国成立的这一天是甲子日！
再如2006年10月1日干支序数为45+5*6+[6/4]+43+1≡15+1+43+1≡0(mod60)，
所以这一天是癸亥日，
也可以不用60作模，而分别用10和12作模对f(y,m,d)求余数，需要计算两次才能求出干支，
不知道这两种算法哪一种快，我还没有比较过。
(3)如果以28为模，就可以求某日对应的廿八宿，
当然廿八宿也可以根据星期和干支用中国剩余定理算出来。
(4)如果以29.5306为模，则可以比较精确的估计出公历某年某月某日的农历日期，
f(y,m,d)≡36524a+[a/4]+365b+[b/4]+f(m)+d (mod29.5306)，
≡[a/4]-5.3522a+10.6328b+[b/4]+0.4694m+[m/2]+[m/9]+([2/m]-[1/m])*(1+m(y))+ε+d(mod29.5306)
据说目前还没有把公历转换成农历的公式，我想上面这个公式应该算是把公历转换成农历的公式了，
如果把ε取一恰当的值，再把所得的结果四舍五入，也许会得到农历日期的精确值，
这需要计算来验证，粗略算了一下，2000年的ε值在6.0134和6.0746之间，
如果得不到农历日期的精确值，我猜想ε对于每年都是不同的值，并且ε是年份y的函数，
如果也不成立的话，就需要换成更精确的周期29.53058867+0.0000000019(y-1900)来支持了，
计算起来更为麻烦，不过用计算机算起来也许很快，只是我不会编程，暂时没法验证。 另外，对于儒略历y年m月d日，有8≤y≤1582且y+m/10≤1583且y+m+d/4≤1593，
此时f(m)=30m+[m/2]+[m/9]+ε+([2/m]-[1/m])*(1-[[y/4]-y/4])，（ε为一待定常数项）
f(y)=36525a+365b+[b/4]，设f(y,m,d)=f(y)+f(m)+d，则
f(y,m,d)=36525a+365b+[b/4]+30m+[m/2]+[m/9]+([2/m]-[1/m])*(1-[[y/4]-y/4])+ε+d，
则f(y,m,d)就是公元y年m月d日一个相对的序数，
应用这个式子可以计算儒略历时的星期和干支等数据，不再赘述。
|