数学研发论坛

 找回密码
 欢迎注册
楼主: 云梦

[原创] 真正的高精度科学计算器

[复制链接]
发表于 2013-1-13 14:28:59 | 显示全部楼层
50# 云梦


能上传部分源代码看看吗?
我想看看源代码,并不是我想窃取你的劳动果实,而是我想学一学如何用visual basic
编程序
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2013-1-13 14:29:54 | 显示全部楼层
软件感觉不错,很想学一学你的程序!我是个vb菜鸟,只懂一点点,所以。。。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2013-1-14 08:29:51 | 显示全部楼层
计算伽马函数的VB原代码示例:
Function Gamma()
Dim xm, MMU As Long
Dim pRes As String
Dim L0 As String
Dim Gmz, srd As Single
srd = DEG
DEG = 1
fsgamma = "0"
Gmz = 0    'x 是否为指定值标志
    Zx = Ax
If Ax.zs > -1 Then   '判断 x 是否整数
  Bx = Ax
  Bx.St = Left(Ax.St, Ax.zs + 1)
  Bx.zs = Len(Bx.St) - 1
  hsfx = "-"
  Call Addition
  If Left(Ax.St, 1) = "0" Then  '如果 x 是整数
    xm = Val(Left(Zx.St, Zx.zs + 1))   '取整数部分进行判断
    If Zx.bz = "" And (xm = 1 Or xm = 2) Then Gmz = 1   'x=1 x=2 时
    If Zx.bz = "-" Or xm = 0 Then Gmz = 2   'x=0,-1,-2,-3.....时
  End If
End If
Ax = Zx
If Gmz = 0 Then
If Ax.bz = "-" Then  '负数的伽马
      fsgamma = 1
      Ax.bz = ""
      Bx.bz = ""
      Bx.St = "1"
      Bx.zs = 0
      hsfx = "+"
      Call Addition
      Zx = Ax
End If
    Px = Ax
    Wx = Ax
    Hg = Ax
    Select Case ychen
       Case Is = 64
          MMU = 90
       Case Is = 128
          MMU = 130
       Case Is = 256
          MMU = 200
       Case Is = 512
          MMU = 540
       Case Is = 1024
         MMU = 8400
    End Select
   
    If Zx.zs < 0 Then
              Ax = Px
              Wx = Px
              Bx.bz = ""
              Bx.St = "1"
              Bx.zs = 0
              hsfx = "+"
              Call Addition
              Wx = Ax
              L0 = Mid(Ax.St, Ax.zs + 2)
              Cx = Hg
              Call mult_
              Hg = Cx
              
        For xm = 2 To MMU 'Γ(10000+z)
            Ax.St = Trim(Str(xm))
            Wx.zs = Len(Ax.St) - 1
            Wx.St = Ax.St + L0
            
            Ax = Wx
            Cx = Hg
            Call mult_
            Hg = Cx
        Next xm
        Px = Wx
    Else
    L0 = Mid(Px.St, Px.zs + 2)
    If Val(Left(Zx.St, Zx.zs + 1)) < MMU Then
            For xm = 1 To MMU 'Γ(10000+z)
                Ax.St = Trim(Str(1 + Val(Left(Px.St, Px.zs + 1))))
                Px.zs = Len(Ax.St) - 1
                Px.St = Ax.St + L0
                Ax = Px
                Cx = Hg
                Call mult_
                Hg = Cx
            Next xm
        Else
           Hg = Zx
           Px = Zx
        End If
    End If
   
    Ux = Hg
    Ax = Px  'Px=1000+z
    Call Lnx 'lnz
    ix = Ax
    '
    Ax = Px
    Bx.bz = ""
    Bx.St = "5"
    Bx.zs = -1
    hsfx = "-"
    Call Addition 'z-1/2
    Cx = Ax
    Ax = ix
    Call mult_  '(z-1/2)*lnz
    Bx = Px
    hsfx = "-"
    Call Addition  '(z-1/2)Lnz-z
    Bx.bz = ""
    Bx.St = epx(15) 'ln(√2Pi)
    Bx.zs = -1
    hsfx = "+"
    Call Addition 'Sx=(z-1/2)Lnz-z+ln(√2Pi)
    Sx = Ax
    '
    Cx = Px   '1000+z
    Ax = Cx
    Call mult_ 'z2
    Tx = Ax    'Tx=z2
    Lx = Px   'z
    Cx.bz = ""
    Cx.St = "12"
    Cx.zs = 1
    Ax = Px
    Call mult_ 'z*12
    Ax = Cx
    Cx.bz = ""
    Cx.St = "1"
    Cx.zs = 0
    Call multx_ '1/(z*12)
    Ax.bz = ""
    Hg = Ax
    Hg.bz = ""
    '求和 ∑
    For xm = 2 To 200 '(-1)[n-1]B2xm/(2xm(2xm-1)z[2xm-1])
        Cx = Lx   'Lx=z(2xm-1)
        Ax = Tx   'Tx=z2
        Call mult_  'z3
        Lx = Cx    'Lx=z3
        Ax = Cx
        
        Cx.bz = ""
        Cx.St = CStr(2 * xm * (2 * xm - 1))
        Cx.zs = Len(Cx.St) - 1
        Call mult_    '2xm(2xm-1)z[2xm-1]
        Ax = Cx
        Cx = pm2(2 * xm)
        Call mult_    'B2xm(分母)*2xm(2xm-1)z[2xm-1]
        Ax = Cx
        Cx = pz2(2 * xm)
        Call multx_   'B2xm(分子)/{B2xm(分母)*2xm(2xm-1)z[2xm-1]}
        Wn = Cx
        Ax = Cx
        Bx = Hg
        hsfx = "+"
        Call Addition
        Hg = Ax
        If Wn.zs < -xChen Then Exit For
    Next xm
    Bx = Hg
    Ax = Sx
    hsfx = "+"
    Call Addition
    Sx = Ax  'Sum 求和结束
    Ax = Ux
    Call Lnx
    Bx = Ax
    Ax = Sx
    hsfx = "-"
    Call Addition
    Sx = Ax
    Ax = Px
    Call Lnx
    Bx = Sx
    hsfx = "+"
    Call Addition
If fsgamma = 1 Then
      Vx = Ax
      Ax.bz = ""
      Ax.St = Pi.Tag
      Ax.zs = 0
      Cx = Zx
      Call mult_
      Ax = Cx
      Call xSin
      Gammbz = Ax.bz
      Cx.bz = ""
      Cx.St = Pi.Tag
      Cx.zs = 0
      Call multx_
      Ax = Cx
      Call Lnx
      Bx = Vx
      hsfx = "-"
      Call Addition
End If
Else
   If Gmz = 1 Then  'x=1,2 时 Γ(x)=1
       Ax.bz = ""
       Ax.St = "0"
       Ax.zs = "0"
   Else
       Text12.Text = "-→∞ (无穷大)"   'x=0,-1,-2,-3....时 Γ(x)=∞
   End If
End If
DEG = srd '恢复三角函数制式(弧度,度,梯度)
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2013-1-14 08:32:21 | 显示全部楼层
计算实例:
Gamma[Pi]=2.2880377953400324179595889090602339228896881533562224411993807454704710066085042825007253044679284747968492456163619736990086693068610684720719929526723030053489815429191959102611884193668794490992029541491572661438362658155059636402680610560412995023428572826183864816514132279851720320360861870024470511009961650028651601111561894239076005577921680490380611694414733267466945500752573250190012765539633365280130189120138691770838670048888356517427638773037998537252598048473914927354467723366483939306412743620 E0

Gamma[-0.9]=-1.0570564109631924262547207974739335769500642917875974825611111967149183759247578901707959771407163566441550287690113503822012272204536656437724552070348384652010074010435712807443049327495784005882882331444625260150178173814720139546897073184335868849137176626529366498372865716422604903386218900984502412837245905669753669097096301130708558392298748589544459136929945703037527541817131388712011444761809695207975411846305421133341651770623266953056232548078614752536919951874359416895723263937984783317257170929 E1
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2013-1-14 08:40:08 | 显示全部楼层
这款计算器不是提供给数学研究用的,它是工程计算人员的计算工具。适合、中、高、大学及非数学专业的学生、技术人员、工程师、数学爱好者等等。因此没有要求很高的计算精度,1024位已经足够了,一般常用64位和128位,这已经大大超出现有计算器的计算能力,且目前已经包含近200个常用函数、分布函数、积分函数、反函数、统计、解一元4次以下根、常用图形的面积和体积及众多的数学常数。
该计算器最大的特点是可以随时把自己常用和喜欢的函数、表达式纳入计算器,为特殊用途提供更多方便。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2013-1-14 08:45:24 | 显示全部楼层
百度里可以找到该计算器前期vb源码,只包括常用函数。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2013-1-14 08:53:37 | 显示全部楼层
上面的是伽马的内部调用函数,下面才是主函数。
Private Sub Command39_Click()  '2012-08.09
Call StoF
Dim pRes As String

Select Case Text8.Text
Case Is = ""  'Gamma函数
    NameF.Text = "Γ(x) ="
    Call SJbz
    Call Gamma
    If Ax.zs > 15 Then
        Cx = Ax
        Ax.St = epx(11)
        Ax.bz = ""
        Ax.zs = -1
        Call mult_
        Ax = Cx
        Ax.St = Mid(Cx.St, Cx.zs + 2)
        For i = 1 To xChen
           If Mid(Ax.St, i, 1) <> "0" Then Ax.St = Mid(Ax.St, i): Exit For
        Next i
        Ax.zs = -i
        Sx.St = Left(Cx.St, Cx.zs + 1)
        Cx.St = epx(10)
        Cx.bz = ""
        Cx.zs = 0
        Call mult_
        Call Exp
        
        'Ln(gamma)
        If Val(Sx.St) < 10 ^ 16 Then
            Ax.zs = Val(Sx.St)
            Call sc(Ax.bz, Ax.St, Ax.zs)
        Else
            If Val(Sx.St) < 32 Then
                pRes = Left(Ax.St, 1) + "." + Mid(Ax.St, 2, 10) + " E" + Sx.St
            Else
                pRes = Left(Ax.St, 1) + "." + Mid(Ax.St, 2, 10) + " E" + Str(Val(Sx.St))
            End If
                Prints.Text = pRes
                pRes = Left(Ax.St, 1) + "." + Mid(Ax.St, 2, ychen - 1) + " E" + Sx.St
            Text11.Text = pRes
        End If
    Else
        Call Exp
        Ax.bz = IIf(fsgamma = 1, Gammbz, Ax.bz)
        Call sc(Ax.bz, Ax.St, Ax.zs)
    End If
可以计算到9x10^11次方的伽马。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2013-1-14 12:01:28 | 显示全部楼层
1. 作为一个同样致力于数学函数的高精度计算的工作者,首先,我对楼主的产品表示肯定和鼓励。
  
  2。这类软件,往大了说,叫做计算机代数系统(Computer algebra system), 大块头的商业软件Maple,Mathcad,Mathematica就属于此类。轻量级的OPen Source的PARI/GP也属于此类。见http://en.wikipedia.org/wiki/Computer_algebra_system。 往小了说,叫做任意精度计算(Arbitrary-precision arithmetic),见http://en.wikipedia.org/wiki/Arbitrary-precision_arithmetic。 以库形式发布的有GMP,CLN,MPFR,MIRACL,NTL。提供库和一体化程序的bc,PARI/GP,Maple,Maxima,Mathematica,sage等。

  3. 关于界面问题。对于计算器这种软件来说。命令行/脚本形式的接口是最好的。他可以完成批量计算,可以解决复杂的问题。但对于初学者,命令行语法不变记忆,较难学习,所以GUI形式的接口可能是最简单的。

  4. 我推荐PARI/GP 和bc(Linux下的命令行计算器),这两个软件都是free的,而且是轻量级的,易于学习,功能也够用。昨天晚上,我还给我儿子讲了45分钟的课,关于PARI用法的。


  顺便说一下,我自己也致力于这类软件的开发,23年前,我就用只有8位10数精度的计算器,用了一个寒假的时间,算出精度达26位数字的对数表,那时我没有见过电脑。考大学时,报考的是计算机专业,这完全是从我的兴趣出发,那时我以为计算机可以直接计算几百位10进制数,上了大学,就有了电脑,就可以将初等函数算到更高的精度了。我从未放弃此类研究。我的一个重要的理想就是开发出比GMP更快的大数运算库。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2013-1-14 14:14:49 | 显示全部楼层
1. 作为一个同样致力于数学函数的高精度计算的工作者,首先,我对楼主的产品表示肯定和鼓励。
  
  2。这类软件,往大了说,叫做计算机代数系统(Computer algebra system), 大块头的商业软件Maple,Mathcad,Mathe ...
liangbch 发表于 2013-1-14 12:01

你把mathematica软件破解了,你就能超过GMP了

评分

参与人数 1经验 -2 收起 理由
wayne -2 风马牛不相及

查看全部评分

毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2013-1-14 14:18:32 | 显示全部楼层
我也想搞一个类似的计算器,只不过是用来求解专业的问题,最近在看vb的书籍,学习是很重要的,
但是我却是一个菜鸟,对很多东西都不懂!
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

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

GMT+8, 2019-7-16 06:19 , Processed in 0.057575 second(s), 19 queries .

Powered by Discuz! X3.4

© 2001-2017 Comsenz Inc.

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