neorobin

回复 28# plp626


    刚看到你回复的 30 楼, 突然不见了
测试结果
pOffs := (d1st * 10000.0) / 3652425 - p;

The pOffsMin is: -0.004045;  the pOffsMax is: 0.001971
The pMin is: 8303;  the pMax is: 8496
pOffsMin = -0.004045,  pOffsMax + 365*10000 / 3652425 = 1.001307

-0.004045 * 3652425 = -14774.059125
14780 比较恰好的调整误差范围到 [0.000000, 1.005354]:

pOffs := (d1st * 10000.0 + 14780) / 3652425 - p;

The pOffsMin is: 0.000000;  the pOffsMax is: 0.006018
The pMin is: 4;  the pMax is: 1296
pOffsMin = 0.000000,  pOffsMax + 365*10000 / 3652425 = 1.005354


pOffs := (d1st * 10000.0 -7300) / 3652425 - p;

The pOffsMin is: -0.006044;  the pOffsMax is: -0.000027
The pMin is: 8303;  the pMax is: 2096
pOffsMin = -0.006044,  pOffsMax + 365*10000 / 3652425 = 0.999309

Pascal 测试代码

var
  d1st, p, pMin, pMax: longint;
  pOffs, pOffsMin, pOffsMax: real;
begin
  pOffsMin := 0;
  pOffsMax := -1;
  for p := 1 to 20000 do
  begin
    d1st := 365 * p + p div 4 - p div 100 + p div 400;
    pOffs := (d1st * 10000.0) / 3652425 - p;
    if pOffs < pOffsMin then begin pOffsMin := pOffs; pMin := p; end
    else if pOffs > pOffsMax then begin pOffsMax := pOffs; pMax := p; end;
    WriteLn(p, '  ', d1st, '  ', pOffs: 0: 6);
  end;
  WriteLn('The pOffsMin is: ', pOffsMin: 0: 6, ';  the pOffsMax is: ', pOffsMax: 0: 6);
  WriteLn('The pMin is: ', pMin: 0, ';  the pMax is: ', pMax: 0);
  WriteLn('pOffsMin = ', pOffsMin: 0: 6, ',  pOffsMax + 365*10000 / 3652425 = ', pOffsMax + 365 * 10000 / 3652425: 0: 6);
  Readln;
end.
复制代码

neorobin

回复 31# neorobin
   
http://alcor.concordia.ca/~gpkatch/gdate-method.html

. In the intervening years of this 400 year cycle, the error maxima are 1.4775 days in year 303 of the cycle, and -0.72 days in year 96 of the cycle. Because of this error, finding the year given a number of days (the inverse function) is not exact, but we can find a very close approximation. Given d days, the year number can be approximated as
y = d / 365.2425     (3)

查看 400 年周期内, 平均数估算天数的误差峰值, 在命令行运行:

>nul (set /a "mi=0,ma=0,ii=1,ia=1"&(for /l %y in (1 1 400) do set /a "y=%y,d=y*365+y/4-y/100+y/400,d*=10000,d=y*3652425-d,ii^=(d-mi>>31)&(y^ii),ia^=(ma-d>>31)&(y^ia),mi^=(d-mi>>31)&(d^mi),ma^=(ma-d>>31)&(d^ma)"))&set m&set i
复制代码

ma=14775
mi=-7200
ia=303
ii=96

plp626

回复 30# neorobin
理论上的分析没经得起实践测试；漏掉了一个变量，分析有误，便删了；
---------------------------
date2index的反函数index2date；
应该可以精简优化到两行150个字节以内
很期待这个函数。。。

Batcher

回复 31# neorobin


    现在还有学校在教Pascal啊，好神奇。

terse

没完全测试 翻版一下网上资料及万年历

@echo off&setlocal enabledelayedexpansion
FOR /l %%i in (0 1 7305155) DO (
    set JD=%%i
    IF !JD! GEQ 2299161 set /a "JD+=1+(JD*100-186721625)/3652425-(JD*100-186721625)/3652425/4"
    set /a "B=JD+1524,Y=(B*100-12210)/36525,D=36525*Y/100"
    set /a "M=(B-D)*10000/306001,D=B-D-306001*M/10000,M=(M-2)%%12+1,Y-=4715+^!(2/M)"
    set "test=%%i:   !y!.!m!.!d!"
    set /a "M=(M+9)%%12+3,Y-=M/13,JD=36525*(Y+4716)/100+306001*(M+1)/10000+D-1524"
    set /a "JD+=(2-Y/100+Y/400)*^!(2299161/JD)"
    set "test=!test!     !JD!"
    if "!JD!" neq "%%i" set "test=!test!     error"
    echo !test!
    if "!test:error=!" neq "!test!" pause
)
PASUE
复制代码

neorobin

回复 33# plp626

把 25 楼的代码改了一点点, 没有实质性的变化, 凑成 2 行 155 字节:

set /a "y=(i*99+145)/36159,y+=i-365*y-y/4+y/100-y/400>>31,dd=i-365*y-y/4+y/100-y/400"
set /a "mi=(5*dd+2)/153,m=(mi+2)%%12+1,y-=m-3>>31,d=1+dd-(153*mi+2)/5"
复制代码

neorobin

回复 35# terse


    不错, 我一直测试到 7305155 没有发现错误

neorobin

回复 34# Batcher

事实上我没有上过电脑课, Pascal: 端庄典雅秀丽的淑女; C: 自由奔放不羁的野马

plp626

回复 36# neorobin

其实你很容易就能精简到140字节以内。。。

neorobin

回复 39# plp626


    弄在一行了, 不必要的变量都省了, 137 字节

set /a "y=(i*99+145)/36159,y+=i-365*y-y/4+y/100-y/400>>31,d=i-365*y-y/4+y/100-y/400,m=(5*d+2)/153,d+=1-(153*m+2)/5,y+=m/10,m=(m+2)%%12+1"
复制代码

plp626

有一个问题
如何测试可以断定date2index一定正确？
如何测试可以断定index2date一定正确？

在date2index都不能断定是正确的情况下，
那么，测试index2date正确性的时候，去调用date2index得出的结论应该是不确定的。。

neorobin

回复 41# plp626

@echo off & setlocal enabledelayedexpansion
REM nextIndex := date2index(1.1.1);
for /l %%y in (1 1 20000) do (
  set "y=%%y"
  set /a "leap = (^!(%%y %% 4) & ^!^!(%%y %% 100)) | ^!(%%y %% 400)"
  if "!leap!"=="1" (set "eFeb=29") else set "eFeb=28"
  for %%t in ("1 31" "2 !eFeb!" "3 31" "4 30" "5 31" "6 30" "7 31" "8 31" "9 30" "10 31" "11 30" "12 31") do (
    for /f "tokens=1,2" %%m in ("%%~t") do (
      set "m=%%m"
      for /l %%d in (1 1 %%n) do (
        set "d=%%d"

        REM testing...
        REM if date2index(y.m.d) <> nextIndex (
          REM error
        REM )
        REM if index2date(date2index(y.m.d)) <> y.m.d (
          REM error
        REM )
        REM nextIndex += 1;
        echo !y!.!m!.!d!
      )
    )
    pause
  )
)
复制代码

和 13 楼说的差不多, 我没想到可以证明正确的方式:

设想一个并不聪明的一定范围内穷举式的测试方式:

公历的历法是明确的:
每年 12 个月: m∈[1..12]
4,6,9,11 月份固定 30 天, d ∈[1..30]
2 月根据年份平闰取 28 天或 29 天,  d∈[1..29]
其余月份固定 31 天, d∈[1..31]

闰年:
(被 4 整除 且 不被 100 整除) 或 被 400 整除 的年份
闰年外的年份都是平年

从这个公历历法的原型定义出发设计一个 尽可能直接符合 且 足够简单 的日期生成算法(假定它是正确的算法):
从一个初始日期开始, 依次生成每一个明天,
把这个算法生成的日期作为 date2index 的参数, 得到相应的 序号:

i_0: 初始日期的序号(值由算法自定);   (_0 表示 0 为下标)

i_p: 初始日期 或者 初始日期之后的 某个日期的序号;
i_q: i_p 对应的那个日期的下一日的序号;

那么, 如果对于任意 i_p 与 i_q, 都得到:

i_q = i_p + 1

就可以把 date2index 看作是测试范围内正确的算法.

在上面的测试中可以同时 作 index2date 测试:
如果始终都有:
index2date( date2index(y.m.d) ) = y.m.d

在 date2index 测试正确的情况下, 就可以认为 index2date 也是正确的.


date2index 和 index2date 在任意参数的情况下的输出都是唯一的, 且是确定的(对于一个确定的输入 I, 不会得到 O1, O2两种或更多可能的输出),

因为算法是基于 数字计算机 的, 输入参数 不具备不确定性(每次输入的值的集合都是确定的),
且算法中不涉及随机函数的调用, 所以 输出是 具有 唯一性和确定性的.

plp626

bat测试太慢；
set/a 表达式等同C的算术表达式；
这里给个超级精简版的C编译器，方便大家测试用；

terse

回复 42# neorobin
以此代码推算至100.2.28后就是100.3.1按实历算的话应该是100.2.29后才是100.3.1
这里是否考虑实历1582前后问题 即1582前每4年一闰  1582后才多了100和400余
也就是说如果这里仅仅代表一种序号的话 那这个序号的合理性应该是计算1582年10月15号后的
譬如 1582年10月4号 儒略日 (2299160)  后一日 便直接跳至1582年10月15日 儒略日 (2299161)
也可能个人理解的不同  故一点浅见而已

neorobin

回复 44# terse

惭愧, 对儒略日, 历法变化等问题完全没有了解或者了解甚少, 前面都完全是按一个固定的历法规则处理的, 如要适应相关的 规则约定 或者 实际情况, 肯定要作出进一步的修改.