欢迎访问 Lu程序设计

Lu数值计算扩展动态库LuMath V1.0

目  录

函数及类别

说     明

1 什么是LuMath

2 LuMath数组(矩阵)对象基本用法
3 LuMath常量
math::real_s 用在函数new[math::real_s,... ...]中,申请数组或矩阵。该对象的基本类型为reals,扩展类型为real_s。
math::pi math::pi=3.141592653589793
math::zeros 用在函数math::subg中,表示返回全0矩阵。
math::ones 用在函数math::subg中,表示返回全1矩阵。
math::eye 用在函数math::subg中,表示返回单位矩阵。
math::rand 用在函数math::subg中,表示返回0~1之间均匀分布的随机数矩阵。
math::randn 用在函数math::subg中,表示返回均值为0,方差为1的标准正态分布的随机数矩阵。
4 LuMath运算符
加法 + 数组(矩阵)加法。运行错误:1:不可识别对象;2:数组不存在;3:数组(矩阵)不匹配;4:内存错误。
减法 - 数组(矩阵)减法。运行错误:1:不可识别对象;2:数组不存在;3:数组(矩阵)不匹配;4:内存错误。
乘法 * 矩阵乘法。运行错误:1:不可识别对象;2:矩阵不存在;3:矩阵不匹配;4:内存错误。
转置 ' 矩阵转置。运行错误:1:不可识别对象;2:矩阵不存在;3:内存错误。
点乘 .* 数组点乘。运行错误:1:不可识别对象;2:数组不存在;3:数组不匹配;4:内存错误。
点左除 ./ 数组点左除。运行错误:1:不可识别对象;2:数组不存在;3:数组不匹配;4:内存错误。
点乘方 .^ 数组点乘方。运行错误:1:不可识别对象;2:数组不存在;3:数组不匹配;4:内存错误。

5 LuMath函数

重载函数  
sin、cos、... ... 该库对以下函数进行了重载,支持real_s类型对象的相关运算:new、copy、len、oset、oget、o、sqrt、exp、ln、lg、sin、cos、tan、asin、acos、atan、sinh、cosh、tanh、abs、floor、ceil、atan2、fmod。
数(组)函数  
math::ra1[2,3.5,-5] 申请一维数组并初始化。 运行错误:1:不可识别参数,只接受实数或整数参数。
math::SetArray 对数组或矩阵进行初始化。
math::outa 输出数组或矩阵。
math::Sum 数组元素求和。
sum 该函数是重载函数,执行与math::Sum完全相同的运算。
math::ndgrid 生成n维函数的辅助格点阵列。
copy 该函数是重载函数。copy(a)将生成对象a的副本;copy(a,b)将数组b复制到数组a,两数组必须等长,并将多维数组视为一维维数;copy(a,abegin,b,bbegin,k)将数组b从地址bbegin开始的k个元素复制到数组a的自地址abegin开始的单元。运行错误:1:参数太多或者太少;2:数组不存在;3:数组不匹配 ;4:内存错误。
oset 该函数是重载函数。设a是一个二维数组对象,则执行形如a[i,j]=2.0的赋值语句时将调用该函数。运行错误:1:参数太少;2:不可识别对象;3:数组不匹配;4:参数不符合要求。
oget 该函数是重载函数。设a是一个二维数组对象,则执行形如a[i,j]的语句时将获得该对象元素的值。运行错误:1:参数太少;2:不可识别对象;3:数组不匹配;4:元素地址错误。
math::arrayfun 以实数数组元素为自变量进行函数调用。实际上,很容易实现以任意元素(lu元素)为自变量进行函数调用的函数lufun
math::reshape 改变数组(矩阵)形状。
矩阵函数  
math::matrix 申请m×n矩阵并初始化。
math::subg 返回一个矩阵的子矩阵。
math::subs 设置一个矩阵的子矩阵。
math::row(a,b,n)
math::row(a,b,x)
产生一个单行矩阵(行向量)。a、b、x为实数,n为整数。a和b是行向量的第一个和最后一个元素,元素总数是n。或者,a是第一个元素,后面的元素按x递增,最后一个元素不超过b。
math::col(a,b,n)
math::col(a,b,x)
产生一个单列矩阵(列向量)。a、b、x为实数,n为整数。a和b是列向量的第一个和最后一个元素,元素总数是n。或者,a是第一个元素,后面的元素按x递增,最后一个元素不超过b。
math::linspace(a,b,n)
math::linspace(a,b,x)
产生一个一维数组。a、b、x为实数,n为整数。a和b是数组的第一个和最后一个元素,元素总数是n。或者,a是第一个元素,后面的元素按x递增,最后一个元素不超过b。
math::zeros(a)
math::zeros(i,j)
将矩阵a重新初始化为全0矩阵,或者生成一个i×j全0矩阵。运行错误:1:数组不存在;2:不是二维数组(矩阵);3:内存错误;4:参数太多或者太少。
math::ones(a)
math::ones(i,j)
将矩阵a重新初始化为全1矩阵,或者生成一个i×j全1矩阵。运行错误:1:数组不存在;2:不是二维数组(矩阵);3:内存错误;4:参数太多或者太少。
math::eye(a)
math::eye(i,j)
将矩阵a重新初始化为单位矩阵,或者生成一个i×j单位矩阵。运行错误:1:数组不存在;2:不是二维数组(矩阵);3:内存错误;4:参数太多或者太少。
math::rand(a)
math::rand(i,j)
将矩阵a重新初始化为0~1之间均匀分布的随机数矩阵,或者生成一个i×j矩阵,并用0~1之间均匀分布的随机数初始化。 运行错误:1:数组不存在;2:不是二维数组(矩阵);3:内存错误;4:参数太多或者太少。
math::randn(a)
math::randn(i,j)
将矩阵a重新初始化为均值为0,方差为1的标准正态分布的随机数矩阵,或者生成一个i×j矩阵,并用均值为0,方差为1的标准正态分布的随机数初始化。 运行错误:1:数组不存在;2:不是二维数组(矩阵);3:内存错误;4:参数太多或者太少。
math::inv(a) 矩阵求逆。运行错误:1:不可识别对象;2:不是n×n矩阵;3:内存错误;4:矩阵奇异。
解方程函数  
math::netn 求非线性方程组一组实根的拟牛顿法。
math::btFindRoot 对分法求非线性方程的实根。

1 什么是LuMath [返回页首]

    LuMath32.dll是一个Lu数值计算扩展动态库,该库以线性代数特别是矩阵运算为基础。

    在LuMath中的函数是通过二级函数命名空间“math”输出的,所有函数均具有类似“math::matrix(...)”的格式。使用!!!using("math");可简化LuMath中的函数访问。

    LuMath库的数组是C格式的,元素序号是基于0的。

    LuMath库函数具有内存消耗低、执行效率高、代码简洁、实用性强的特点。

    LuMath的源代码是开放的,源代码下载:http://www.forcal.net/xiazai/lu1/lu1code.rar

    LuMath库中所用的算法或许不是最好的,如果您有好的算法,可以方便地进行替换,提升LuMath的性能。

    LuMath库可用于开发极致性能的应用程序,是熟悉C/C++、Fortran的数学爱好者的极佳的练手工具,同时也期望对一般的数值计算用户提供越来越多的方便。

    限于作者水平,期待与朋友们共同完善LuMath!如果您有什么好的算法,任何改进的意见或建议,请与作者联系。

    本库仅是抛砖引玉,任何喜欢LuMath的个人、团队或商业公司可基于此库开发商业程序。

2 LuMath数组(矩阵)对象基本用法 [返回页首]

    借助Lu编译器的支持,LuMath中的数组(矩阵)对象有以下简洁的用法。

    设A、B是LuMath数组(矩阵)对象:

      A.=B                 //将对象B赋值给对象A

      A.=A+B               //将对象A+B赋值给对象A

      A[2]=A[2]+1          //将一维数组A第2个单元的值加1

      A[2,3,5]=A[3,1,2]+1  //将三维数组A的第[2,3,5]个单元的值赋值为A[3,1,2]+1

      A.*2.0               //数组(矩阵)A的每一个元素乘以2.0

    若A是一个LuMath矩阵:

      A(type)              //返回与矩阵A同样大小的矩阵,并根据type的值进行初始化。type的可选值为zeros(全0矩阵)、ones(全1矩阵)、eye(单位矩阵)、rand(0~1之间均匀分布的随机数)、randn(均值为0,方差为1的标准正态分布的随机数)。

      A(all:j)             //取矩阵A第j列所有元素。all是一个常量,在这里表示取所有行

      A(i:all)             //取矩阵A第i行所有元素。all是一个常量,在这里表示取所有列

      A(i,i+m:all)         //取矩阵A第i~i+m行所有元素。all是一个常量,在这里表示取所有列

      A(all:j,j+m)         //取矩阵A第j~j+m列所有元素。all是一个常量,在这里表示取所有行

      A(i,i+m:j,j+n)       //取矩阵A第i~i+m行,并在第j~j+n列的所有元素

      A(all:j)=B           //设置矩阵A第j列所有元素,B是一个数组。all是一个常量,在这里表示取所有行

      A(i:all)=B           //设置矩阵A第i行所有元素,B是一个数组。all是一个常量,在这里表示取所有列

      A(i,i+m:all)=B       //设置矩阵A第i~i+m行所有元素,B是一个矩阵。all是一个常量,在这里表示取所有列

      A(all:j,j+m)=B       //设置矩阵A第j~j+m列所有元素 ,B是一个矩阵。all是一个常量,在这里表示取所有行

      A(i,i+m:j,j+n)=B     //设置矩阵A第i~i+m行,并在第j~j+n列的所有元素 ,B是一个矩阵。

    [例子]

!!!using["math"];
main(:t0,k,i,a,b)=
{
    t0=clock(),
    k=zeros(5,5),        
//生成5×5矩阵k,初始化为0
    i=0,(++i<=1000).while{
//循环计算1000次
        a=rand(5,7), b=rand(7,5),
//生成5×7矩阵a,7×5矩阵b,用0~1之间的随机数初始化
        k.=k+a*b+a(0,4:1,5)*b(1,5:0,4)+a(all:6)*b(3:all)
    },
    k.outa(),  
          //输出矩阵k
    [clock()-t0]/1000.
   //得到计算时间,秒
};

5 LuMath函数 [返回页首]

[返回页首] math::SetArray(a : x1,x2,... ...):对数组或矩阵进行初始化

    用法1:math::SetArray(a : x1,x2,... ...):对数组或矩阵a进行初始化,从地址0开始存放数据x1,x2,... ...,初始化数据可以少于或多于数组大小,多余的数据被忽略。

    用法2:math::SetArray(a : "1  -2 3 2.3 ... ..."): 对数组或矩阵a进行初始化,从地址0开始存放数据,数据存放在字符串"1  -2 3 2.3 ... ..."中,初始化数据可以少于或多于数组大小,多余的数据被忽略。

    运行错误:1:非法的数组对象;2:非法的Lu字符串及其他参数;3:非法的数字字符串。

[返回页首] math::outa(a, precision, rl, index):输出数组或矩阵

    a数组或矩阵,按列表方式输出 。
   
precision:有效数字位数(2~18),缺省是6位有效数字。
   
rl:左对齐或右对齐输出,非0表示右对齐。 缺省是右对齐。
    index:index非0时输出数组元素下标。缺省是不输出数组元素下标。

    运行错误:1:参数太多或者太少;2: 参数不符合要求;3:数组不存在;3:内存错误。

[返回页首] math::Sum(a,i):数组元素求和

    设有m维数组(含矩阵):a(n1,n2,... ...,nm),则Sum(a,i)对第i维求和,返回一个m-1维数组。若a是一维数组、1×k矩阵、k×1矩阵、i<1或者i>m,则Sum函数返回所有数组元素的和。Sum(a)相当于Sum(a,m)。

    运行错误:1:参数太多或者太少;2:数组不存在 ;3:数组维数太多;4:内存错误。

[返回页首] math::ndgrid(x1,x2,...,xn : &y1,&y2,...,yn):生成n维函数的辅助格点阵列

    x1,x2,...,xn:均是一维数组。
    &y1,&y2,...,yn:返回n维函数的辅助格点阵列。
    返回值:无意义。

    运行错误:1:参数太多或者太少;2:数组不存在;3: 不是一维数组;4:内存错误。

    例子:

!!!using["math"];
mvar:
ndgrid[linspace(0.0,1.0,0.0011),linspace(1.0,2.0,0.0009),&x,&y],
Sum[
cos(1.0-sin((1.2).*x.^(y./2.0)+cos(1.0-sin((1.2).*y.^(x./2.0))))),0];

    相当于:

f(x,y)=cos(1.0-sin(1.2*x^(y/2.0)+cos(1-sin(1.2*y^(x/2.0)))));
sum[@f:0.0,1.0,0.0011:1.0,2.0,0.0009];

[返回页首] math::arrayfun(hFor:x1,x2,...,xn : xx):以 实数数组元素为自变量进行函数调用。

    hFor为函数句柄;x1,x2,...,xn : xx为等长的 实数数组(忽略维数);xx用于存放结果,可缺省。若缺省xx,返回一个实数数组存放结果,否则该函数返回xx。该函数以数组x1,x2,...,xn的第i个元素为自变量,用函数hFor进行计算,结果存入xx的第i个单元。

    运行错误:1:至少需要2个参数;2:指定的表达式不存在;3:数组个数与自变量个数不相等;4:不是实数数组;5:数组长度不匹配;6:内存错误 ;7:函数应返回一个实数。

    例子:

!!!using["math"];
f(x,y)=x+y;
arrayfun[@f,linspace(1.,5.,5),linspace(1.,5.,5)].o[];

[返回页首] math::reshape(a,... ...):改变数组(矩阵)形状

    格式1:math::reshape(a):将矩阵重置为一维数组。
    格式2:math::reshape(a,i,j):将矩阵a(或一维数组)重置为i×j矩阵,保持总元素不变。
    格式3:math::reshape(a,i,j,k,...):改变多维数组的维数大小,保持维数和总元素不变。

    运行错误代码:1:参数太多或者太少;2:实数数组不存在;3:数组维数小于1;4:数组不匹配。

[返回页首] math::matrix(m,n:x11,x12,... ...):或者math::matrix(str):申请m×n矩阵并初始化。

    m,n:矩阵的行和列。初始化数据 x11,x12,...,x1n,x21,x22,...,x2n,... ...,xm1,xm2,...,xmn 可以缺省部分及全部,多余的数据被忽略。
    str:字符串,包含矩阵数据,数据之间可以用空格或逗号分隔。此时,将自动确定矩阵的行和列。例如以下是一个2×3矩阵字符串:

"1    1.2  3,
 2.2  3,   5"

    运行错误:1:参数不符合要求 ;2:字符串矩阵每行中数据个数不同;3:字符串中数据非法。

[返回页首] math::subg(a, ... ...):返回一个矩阵的子矩阵

    subg(a):返回矩阵a的副本。
   
subg(a:type):返回与矩阵a同样大小的矩阵,并根据type的值进行初始化。type的可选值为zeros(全0矩阵)、ones(全1矩阵)、eye(单位矩阵)、rand(0~1之间均匀分布的随机数)、randn(均值为0,方差为1的标准正态分布的随机数)。
   
subg(a:all:j): 取矩阵a第j列所有元素。all是一个常量,在这里表示取所有行。
   
subg(a:i:all):取矩阵a第i行所有元素。all是一个常量,在这里表示取所有列。
   
subg(a:i,i+m:all):取矩阵a第i~i+m行所有元素。all是一个常量,在这里表示取所有列。
   
subg(a:all:j,j+m):取矩阵a第j~j+m列所有元素。all是一个常量,在这里表示取所有行。
   
subg(a:i,i+m:j,j+n):取矩阵a第i~i+m行,并在第j~j+n列的所有元素。

    运行错误:1:参数太多或者太少;2:数组不存在;3:不是二维数组(矩阵);4:不可识别的类型标识type;5:无法定位矩阵的行或列 ;6:内存错误。

    [例子]

!!!using["math"];
main(:a)=
    a=matrix[6,5].rand(),
    a.outa(),
    a.subg(all:3).outa(),
    a.subg(3:all).outa(),
    a.subg(3,5:2,3).outa();

[返回页首] math::subs(a, ... ...,b):设置一个矩阵的子矩阵

    subs(a:b):将矩阵b赋值给矩阵a,a与b应具有相同的结构。
   
subs(a:all:j,b):对矩阵a第j列所有元素赋值,b是一个等长的一维数组。all是一个常量,在这里表示取所有行。
   
subs(a:i:all,b):对矩阵a第i行所有元素赋值,b是一个等长的一维数组。all是一个常量,在这里表示取所有列。
   
subs(a:i,i+m:all,b):对矩阵a第i~i+m行所有元素赋值,b是一个相同大小的矩阵。all是一个常量,在这里表示取所有列。
   
subs(a:all:j,j+m,b):对矩阵a第j~j+m列所有元素赋值,b是一个相同大小的矩阵。all是一个常量,在这里表示取所有行。
   
subs(a:i,i+m:j,j+n,b):对矩阵a第i~i+m行,并在第j~j+n列的所有元素赋值,b是一个相同大小的矩阵。
    返回值:
返回a。

    运行错误:1:参数太多或者太少;2:数组不存在;3:不是二维数组(矩阵);4:数组不匹配;5:无法定位矩阵的行或列。

[返回页首] math::netn(f,x,eps,t,h,k):求非线性方程组一组实根的拟牛顿法

f:自定义二元或2n(n>1)元函数句柄,但必须由二级函数HFor获得该句柄。由用户自编。

    如果f是二元函数,则两个参数为等长的一维数组:

f(x,y)=    //二元函数,x[i]为自变量,y[i]为方程右端点函数
{
    y[0]=f1(x[0],x[1],...,x[n-1]),
    y[1]=f2(x[0],x[1],...,x[n-1]),
    ... ...
    y[n-1]=fn(x[0],x[1],...,x[n-1])
};

    或者:

f(x1,x2,...,xn,y1,y2,...,yn)=    //2n元函数,xi为自变量,yi为方程右端点函数
{
    y1=f1(x1,x2,...,xn),
    y2=f2(x1,x2,...,xn),
    ... ...
    yn=fn(x1,x2,...,xn)
};

x:双精度实型一维数组,长度不小于n,存放一组初值,返回方程组的一组实根。
eps:控制精度(eps>0)。可缺省该参数,缺省值eps=1e-6。
t:控制h大小的变量,0<t<1。可缺省该参数,缺省值h=0.1。
h:增量初值,在本函数中将被破坏。可缺省该参数,缺省值h=0.1。
k:允许的最大迭代次数。可缺省该参数,缺省值k=100。

返回值:返回值为实际迭代次数。若返回值为-1表示代数方程奇异,若返回值为-2表示β=0,这两种情况可放宽精度要求、改变各个初值或改变各个方程顺序再试;若返回值等于0说明迭代了k次仍未满足精度要求,程序工作失败。

运行错误代码:1:参数太多或者太少;2:指定的表达式不存在;3:方程f参数个数不符合要求;4:实数数组不存在;5:数组元素太少;6:精度、控制变量、迭代次数等不符合要求;7:内存错误 ;8:用户删除了自定义二元函数中的工作数组。

    [例子]

f(x1,x2,x3,y1,y2,y3)=      //函数定义
{
    y1=x1*x1+x2*x2+x3*x3-1.0,
    y2=2.0*x1*x1+x2*x2-4.0*x3,
    y3=3.0*x1*x1-4.0*x2+x3*x3
};
!!!using["math"];
main(:x,i)={
    x=new{real_s,data : 1.0,1.0,1.0},
   //设置初值为1.0,1.0,1.0
    i=netn[@f,x],
          //拟牛顿法解方程
    x.outa(), 
            //输出结果
    i
                      //返回迭代次数
};

    或者:

f(x,y)=                   //函数定义
{
    y[0]=x[0]*x[0]+x[1]*x[1]+x[2]*x[2]-1.0,
    y[1]=2.0*x[0]*x[0]+x[1]*x[1]-4.0*x[2],
    y[2]=3.0*x[0]*x[0]-4.0*x[1]+x[2]*x[2]
};
!!!using["math"];
main(:x,i)={
    x=new{real_s,data : 1.0,1.0,1.0},
   //设置初值为1.0,1.0,1.0
    i=netn[@f,x],
          //拟牛顿法解方程
    x.outa(), 
            //输出结果
    i
                      //返回迭代次数
};

[返回页首] math::btFindRoot(f,a,b,h,k,eps):对分法求非线性方程的实根

f:自定义一元函数句柄。函数自变量不能重新赋值。
a,b:求根区间的左端点和右端点(b>a)。
h:搜索求根时采用的步长(h>0)。
k: 可能的实根个数。可缺省,缺省值为10。
eps:控制精度(eps>0)。可缺省,缺省值为1e-6。

返回值:返回一个数组(存放搜索到的实根),或者返回nil。若返回数组的长度为k,则有可能在求根区间[a,b]内的实根未搜索完。

运行错误代码:1:参数太多或者太少;2:指定的函数不存在;3:不是一元方程;4:自变量不能重新赋值;5:参数设置不符合要求;6:内存错误;7:指定的函数应返回一个实数。

例子:

f(x)=x^6-5*x^5+3*x^4+x^3-7*x^2+7*x-20;    //函数定义
math::btFindRoot[@f,-2.,5.,0.2].o[];


版权所有© Lu程序设计 2011-2012,保留所有权利
E-mail: forcal@sina.com
  QQ:630715621
最近更新: 2013年11月08日