数学建模社区-数学中国

标题: Lu基于系统内置对象创建扩展数据类型，小矩阵乘效率测试 [打印本页]

---

作者: forcal    时间: 2011-10-23 09:37
标题: Lu基于系统内置对象创建扩展数据类型，小矩阵乘效率测试
    本例中，我们将基于系统内置实数数组创建矩阵（matrix）类型，即：基本类型为luDynData_realarray（标识实数数组），扩展类型为matrix（标识矩阵）。为了简单，我们仅处理二维实数数组即矩阵类型。

    基本要点：

- c* H6 [- g4 L, W; @9 ^, g) f9 Z    （1）为扩展类型matrix编写运算符重载函数OpMatrix。

, M8 X( U  ^) i, B7 k& }    （2）用函数LockKey将重载函数OpMatrix注册到Lu，锁定的键的类型即为matrix，要注册为常量，以便于使用。

    （3）为扩展类型matrix编写其他操作函数（本例未提供）。

    （4）用函数LockKey解锁键matrix（本例中，程序退出时会自动解锁，故可以不用）。

1. #include <windows.h>
   0 u/ j! @/ z2 r6 X6 t4 W3 b
2. #include <iostream>
   
3. #include <math.h>
   + F2 `; T1 _  h! W# v& D
4. #include "lu32.h"
   - T* H/ y5 e7 G! r! }% h& q1 j: ]
5. 
   . `$ h2 z* Y$ m- T
6. #pragma comment( lib, "lu32.lib" )
   * b1 u) _" F  K8 V$ J/ ~
7. 
   ; p+ ~$ }# q, I
8. using namespace std;
   & y/ u) g8 X  O+ w
9. 
   
10. luKEY Matrix=-1000;        //标识矩阵类型，最终的Matrix由LockKey决定
    ! h0 v- O9 D- R- ]/ ~
11. 
    3 I5 d2 A/ e* i, \- Z) H' o: x
12. void _stdcall LuMessage(wchar_t *pch)//输出动态库信息，该函数注册到Lu，由Lu二级函数调用
    9 R1 V& O5 p% p
13. {
    2 A6 N7 B  Z" a* T0 H" z) K; ^
14.         wcout<<pch;
    
15. }
    & d( b& }$ w/ z0 ^
16. void _stdcall DelMatrix(void *me)        //用于LockKey函数，因为是基于系统内置实数数组创建矩阵，故该函数什么也不做
    
17. {
    
18. }
    
19. LuData _stdcall OpMatrix(luINT mm,LuData *xx,void *hFor,int theOperator)        //运算符重载函数，用于LockKey函数
    
20. {
    
21.         LuData a;
    : B+ o: W; R) V% H6 y2 d4 n
22.         luRealArray *pRealArray1,*pRealArray2,*pRealArray3;
    5 R% W2 S, h4 h4 `2 D+ }
23.         luVOID i,j,k,m,n,u,v;
    
24.         double *pa,*pb,*pc;
    $ x, `; N+ h, A" `+ O2 }7 ~
25.         luMessage pMessage;
    6 T- N( ]$ R& d: x/ \$ G5 E4 o
26.         wchar_t wchNum[32];
    
27.         char chNum[32];
    
28.         a.BType=luStaData_nil; a.VType=luStaData_nil; a.x=0;
    
29.         switch(theOperator)
    4 u# e1 h8 i9 F% O: `+ p' ?' u9 V
30.         {
    
31.         case 2:        //重载运算符*
    % V3 u' R2 j6 r/ R
32.                 pRealArray1=(luRealArray *)SearchKey((char *)&(xx->x),sizeof(luVOID),luDynData_realarray);
    
33.                 pRealArray2=(luRealArray *)SearchKey((char *)&((xx+1)->x),sizeof(luVOID),luDynData_realarray);
    9 b9 K* i. H# c
34.                 if(!pRealArray1 || !pRealArray2) break;        //对象句柄无效，不是实数数组
    5 _2 w3 w& @8 V0 n. ~! y9 N7 A- b
35.                 if(pRealArray1->DimLen!=2 || pRealArray2->DimLen!=2) break;        //不是二维实数数组（矩阵）
    1 C, e" _+ s! ^3 X+ B( b
36.                 if(pRealArray1->Dim[1]!=pRealArray2->Dim[0]) break;        //维数不匹配
    
37.                 pRealArray3=(luRealArray *)NewSysObj(luDynData_realarray,pRealArray1->Dim[0]*pRealArray2->Dim[1],2);        //创建矩阵对象
    
38.                 if(!pRealArray3) break;
    , A' J1 ^4 h4 f- a
39.                 pRealArray3->Dim[0]=pRealArray1->Dim[0]; pRealArray3->Dim[1]=pRealArray2->Dim[1];        //设置矩阵维数大小
    
40.                 pa=pRealArray1->Array; pb=pRealArray2->Array; pc=pRealArray3->Array;
    5 E, ~) t( |7 a$ P# i* F
41.                 m=pRealArray1->Dim[0]; n=pRealArray1->Dim[1]; k=pRealArray2->Dim[1];
    
42.                 for(i=0; i<m; i++)        //矩阵乘
    
43.                 {
    5 s4 M+ ^0 N$ q
44.                         for(j=0; j<k; j++)
    / U" W" \1 K) e* H7 e
45.                         {
    
46.                                 u=i*k+j; pc[u]=0.0;
    " W& M9 M: p7 c7 ~; g- S
47.                                 for (v=0; v<n; v++)
    
48.                                 {
    
49.                                         pc[u]=pc[u]+pa[i*n+v]*pb[v*k+j];
    
50.                                 }
    / x' J7 d! d( e( I, j& K3 _
51.                         }
    6 x# k9 p! M, J, O, v6 c- S8 _, g
52.                 }
    1 x: i& g! E* m# G1 F3 G
53.                 FunReObj(hFor);        //告诉Lu，返回一个动态对象
    7 F! k0 J2 Y* i1 l" K& W5 O3 B- B( e7 x
54.                 a.BType=luDynData_realarray; a.VType=Matrix; a.x=0; *(luVOID *)&(a.x)=(luVOID)pRealArray3;
    ' }% `; x0 T- M& j
55.                 break;
    
56.         case 25:        //重载运算符.*
    
57.                 pRealArray1=(luRealArray *)SearchKey((char *)&(xx->x),sizeof(luVOID),luDynData_realarray);
    
58.                 pRealArray2=(luRealArray *)SearchKey((char *)&((xx+1)->x),sizeof(luVOID),luDynData_realarray);
    
59.                 if(!pRealArray1 || !pRealArray2) break;        //对象句柄无效，不是实数数组
    
60.                 if(pRealArray1->DimLen!=2 || pRealArray2->DimLen!=2) break;        //不是二维实数数组（矩阵）
    
61.                 if(pRealArray1->Dim[0]!=pRealArray2->Dim[0] || pRealArray1->Dim[1]!=pRealArray2->Dim[1]) break;        //维数不相同
    
62.                 pRealArray3=(luRealArray *)NewSysObj(luDynData_realarray,pRealArray1->ArrayLen,2);        //创建矩阵对象
    " f& A: O; X5 Y' W, L+ S
63.                 if(!pRealArray3) break;
    % ^$ Z' X3 q% Q
64.                 pRealArray3->Dim[0]=pRealArray1->Dim[0]; pRealArray3->Dim[1]=pRealArray1->Dim[1];        //设置矩阵维数大小
    ' m" s" Y& j3 u
65.                 for(i=0;i<pRealArray1->ArrayLen;i++) pRealArray3->Array[i]=pRealArray1->Array[i]*pRealArray2->Array[i];//矩阵点乘
    ; k( w: j" W# x/ h2 o7 ?6 V
66.                 FunReObj(hFor);        //告诉Lu，返回一个动态对象
    
67.                 a.BType=luDynData_realarray; a.VType=Matrix; a.x=0; *(luVOID *)&(a.x)=(luVOID)pRealArray3;
    
68.                 break;
    
69.         case 46:        //重载函数new
    ! b+ R' T& _* p: }4 J$ X4 A# b
70.                 a=ExeOperator(mm,xx,hFor,theOperator,luDynData_realarray);        //直接调用基本类型luDynData_realarray的new函数
    
71.                 if(a.VType==luDynData_realarray) a.VType=Matrix;        //设置扩展类型为自定义的Matrix类型
    
72.                 break;
    ! Q$ g( Q4 e; S9 I' E
73.         case 49:        //重载函数o
    9 Z- A: R& }6 u3 \6 N
74.                 pMessage=(luMessage)SearchKey("\0\0\0\0",sizeof(luVOID),luPubKey_User);
    & }" v5 [4 v1 J% l' t* r3 |$ \
75.                 if(!pMessage) break;
    
76.                 pRealArray1=(luRealArray *)SearchKey((char *)&(xx->x),sizeof(luVOID),luDynData_realarray);
    
77.                 if(!pRealArray1) break;        //对象句柄无效，不是实数数组
    
78.                 if(pRealArray1->DimLen!=2) break;        //不是二维实数数组（矩阵）
    : t- c9 I9 {, E4 q. c
79.                 pa=pRealArray1->Array;
    
80.                 m=pRealArray1->Dim[0]; n=pRealArray1->Dim[1]; k=0;
    / R, P. a  [6 {# W# A( ], f
81.                 for(i=0; i<m; i++)        //输出矩阵
    / _4 w+ [4 U/ d/ D
82.                 {
    
83.                         pMessage(L"\r\n"); k+=2;
    
84.                         for(j=0; j<n; j++)
    9 M/ L+ D* N0 V2 G1 P
85.                         {
    
86.                                 _gcvt_s(chNum,pa[i*n+j],16);
    8 N% z0 `: O# a& H; G& ^) J
87.                                 for(u=0;chNum[u];u++) {wchNum[u]=chNum[u]; k++;}
    
88.                                 wchNum[u]='\0';
    
89.                                 pMessage(wchNum); pMessage(L"  "); k+=2;
    , f6 S( {/ `$ n) U
90.                         }
    
91.                 }
    % Z: L5 W/ [* h6 O+ H7 W0 D6 v5 a
92.                 pMessage(L"\r\n"); k+=2;
    
93.                 a.BType=luStaData_int64; a.VType=luStaData_int64; a.x=k;        //按函数o的要求，返回输出的字符总数
    * k- z& m1 K* T
94.                 break;
    
95.         default:
    7 R. {+ V% a: Q
96.                 break;
    
97.         }
      n0 W. V9 n% d! W' _3 t
98.         return a;
    
99. }
    ) M* S* ]/ S9 K. E, \$ m7 o: x
100. void main(void)
     
101. {
     & N# a! `+ _" E( V2 Q5 ?
102.         void *hFor;                //表达式句柄
     , ]% F( r& u& _- z" F1 g
103.         luINT nPara;                //存放表达式的自变量个数
     
104.         LuData *pPara;                //存放输入自变量的数组指针
     " H2 \- T1 Q) ~$ t
105.         luINT ErrBegin,ErrEnd;        //表达式编译出错的初始位置和结束位置
     
106.         int ErrCode;                //错误代码
     * I2 I, j- Y7 U
107.         void *v;
     
108.         wchar_t ForStr[]=L"o{new[matrix,2,3,data: 0.,1.,2.;3.,4.,5.]*new[matrix,3,2,data: 1.,2.;3.,4.;5.,6.]}";//字符串表达式，矩阵乘
     & t. s! ^4 T. b; A" X3 K+ ^" k
109.         //wchar_t ForStr[]=L"o{new[matrix,2,3,data: 0.,1.,2.;3.,4.,5.].*new[matrix,2,3,data: 1.,2.,3.;4.,5.,6.]}";//字符串表达式，矩阵点乘
     
110.         LuData Val;
     4 M* @) }+ s8 O1 Q6 @6 l; ~
111.         if(!InitLu()) return;        //初始化Lu
     ) y3 d& z1 @7 r. B' t# l, s) d3 @
112.         while(LockKey(Matrix,DelMatrix,OpMatrix)){Matrix++;}        //锁定一个键，用于存储矩阵扩展类型
     
113. 
       P' P# U9 U* [: L4 e9 i4 M
114.         Val.BType=luStaData_int64; Val.VType=luStaData_int64; Val.x=Matrix;        //定义整数常量
     
115.         SetConst(L"matrix",&Val);        //设置整数常量
     / N5 O6 T' t" Y2 k5 o
116.         InsertKey("\0\0\0\0",4,luPubKey_User,LuMessage,NULL,NULL,1,v); //使Lu运行时可输出函数信息
     & w1 l, |! l! T: k  n- h
117.         wcout.imbue(locale("chs"));        //设置输出的locale为中文
     1 Z; p! U$ R( v
118.   
     
119.         ErrCode=LuCom(ForStr,0,0,0,hFor,nPara,pPara,ErrBegin,ErrEnd); //编译表达式
     
120.         if(ErrCode)
     
121.         {
     3 S- m8 k6 L# N" S: l
122.                 wcout<<L"表达式有错误！错误代码："<<ErrCode<<endl;
     
123.         }
     
124.         else
     $ X0 [$ d" G% K- Z% Z# {6 {8 V6 m
125.         {
     
126.                 LuCal(hFor,pPara);        //计算表达式的值
     
127.         }
     1 Z! I. X( c; {
128.         LockKey(Matrix,NULL,OpMatrix);//解锁键Matrix，本例中，该函数可以不用
     
129.         FreeLu();                        //释放Lu
     5 m% e: }' {; v& U* T5 q: R
130. }

复制代码

习题：

    （1）自定义矩阵的加、减、左除、右除、点左除等运算，自编测试字符串代码，重新编译运行程序，观察计算结果。
5 \4 T/ I7 Q" Z& d7 d" P( ^( Y
7 C4 W* u: c$ y* ?    （2）小矩阵乘效率测试。编译运行以下Lu字符串代码：

1. main(:a,b,c,d,t,i)=
   7 c1 f5 g4 ^, m* F
2.     a=new[matrix,2,2,data:1.,2.,2.,1.],
   ) }' o4 ?, ^, {, f6 t: ?) A' v
3.     b=new[matrix,2,2,data:2.,1.,1.,2.],
   . c$ X4 J' V' N9 ]
4.     c=new[matrix,2,2,data:2/3.,-1/3.,-1/3.,2/3.],
   - T: x  u5 R6 k/ V
5.     t=clock(),
   
6.     d=a*b, i=0, while{i<1000000, d=d*c*b, i++},
   1 n) E/ l, ]3 C: c/ c" V
7.     o{d, "time=",[clock()-t]/1000.," seconds.\r\n"}

复制代码

C/C++中的字符串定义为：

1. wchar_t ForStr[]=L"main(:a,b,c,d,t,i)= a=new[matrix,2,2,data:1.1,2.,2.,1.], b=new[matrix,2,2,data:2.,1.,1.,2.], c=new[matrix,2,2,data:2/3.,-1/3.,-1/3.,2/3.], t=clock(), d=a*b, i=0, while{i<1000000, d=d*c*b, i++}, o{d, \"time=\",[clock()-t]/1000.,\" seconds.\r\n\"}";//字符串表达式

复制代码

结果：

1. 4. 5.
   / r8 e8 |* R5 G2 e
2. 5. 4.
   3 F- K2 K9 V! K6 c, N
3. time=0.875 seconds.
   
4. 请按任意键继续. . .

复制代码

Matlab 2009a 代码：

1. a=[1.,2.;2.,1.];
   
2. b=[2.,1.;1.,2.];
   8 w5 ~* m  f  [6 Q% j6 N
3. c=[2/3.,-1/3.;-1/3.,2/3.];
   / L" ?$ K! y2 s/ z
4. tic,
   $ R( f- N0 I9 r! ~$ Z4 ?
5. d=a*b;
   
6. for i=1:1000000
   
7.     d=d*c*b;
   0 @3 v; C, w4 ]  R
8. end
   5 @$ N$ o2 N3 D0 C
9. d,
   
10. toc

复制代码

结果：

1. d =
   5 ]6 g: P2 i4 a5 t4 @% w* v
2.      4     5
   
3.      5     4
   
4. Elapsed time is 2.903034 seconds.

复制代码

本例矩阵乘效率测试，Lu的速度超过了Matlab，主要在于Lu有更高的动态对象管理效率。

---

欢迎光临 数学建模社区-数学中国 (http://www.madio.net/)

Powered by Discuz! X2.5