数学建模社区-数学中国

标题: 在Lu中创建自定义数据类型，小矩阵乘效率测试 [打印本页]

---

作者: forcal    时间: 2011-10-24 18:54
标题: 在Lu中创建自定义数据类型，小矩阵乘效率测试
    本例中，我们将自定义矩阵（matrix）类型，基本类型和扩展类型均为matrix（标识矩阵）。
5 G1 F. C; y/ |- X9 u8 Y
" n0 e) H! _* C0 j* m6 m$ t$ Z    基本要点：
/ \2 Y  ^7 Y8 H! H! ?7 T  u9 \
    （1）编写生成矩阵（matrix）的函数NewMatrix和销毁矩阵的函数DelMatrix。
0 `: W# @" t* T
    （2）为自定义类型matrix编写运算符重载函数OpMatrix。

: h$ {% x' ^) i, Z) ^9 w    （3）用函数LockKey将重载函数OpMatrix注册到Lu，锁定的键的类型即为matrix，要注册为常量，以便于使用。
3 z$ g1 \' A( l# b$ _& K; Q7 G
    （4）为自定义类型matrix编写其他操作函数（本例未提供）。

; U% F# F5 c5 V    （5）用函数LockKey解锁键matrix（本例中，程序退出时会自动解锁，故可以不用）。

1. #include <windows.h>
   1 a& E+ W/ I- w% O
2. #include <iostream>
   
3. #include <math.h>
   / ~' g( D/ f* l( G  C, [% r
4. #include "lu32.h"
   ; o7 U4 n. H% A/ I- k2 p
5. #pragma comment( lib, "lu32.lib" )
   3 f7 z4 f4 h7 J
6. using namespace std;
   
7. //自定义矩阵
   
8. class myMatrix
   - h# G- k' A/ }  P/ x' w) E
9. {
   * j0 c7 X5 @( s# T9 v9 e3 Z
10. public:
    
11.         double *Array;        //数据缓冲区
    
12.         luVOID ArrayLen;        //数据缓冲区长度
    9 |/ X/ `* h- z5 T& p2 @; K
13.         luVOID Dim[2];        //矩阵维数
    ( {: F4 C) q8 Q
14.         myMatrix(){Array=NULL; ArrayLen=0; Dim[0]=0; Dim[1]=0;}
    3 }/ W) W6 k/ [- Q, U0 o
15.         ~myMatrix()
    
16.         {
    1 u( `$ i1 G2 ^, c
17.                 if(Array) delete[] Array;
    + u/ Z! E  w1 l9 b4 ]! p; p
18.         }
    
19. };
    % i' N/ S( N& i  _0 |' R5 w
20. luKEY Matrix=-1000;        //标识矩阵类型，最终的Matrix由LockKey决定
    
21. void _stdcall LuMessage(wchar_t *pch)//输出动态库信息，该函数注册到Lu，由Lu二级函数调用
    
22. {
      l5 F5 b" z) Q
23.         wcout<<pch;
    " M# U0 |* s# a2 V" F6 G( [3 _
24. }
    / T* ]& E% z4 J  N) Z
25. void _stdcall DelMatrix(void *me)        //用于LockKey函数及InsertKey函数，使Lu能自动销毁myMatrix对象
    
26. {
    $ n- d3 x7 c" j5 P
27.         delete (myMatrix *)me;
    
28. }
    * V: R& s! Q. m2 K
29. myMatrix * _stdcall NewMatrix(luVOID m,luVOID n)        //生成一个myMatrix对象
    ; f- @- T' g" t
30. {
    
31.         myMatrix *pMatrix;
    
32.         luVOID k;
    
33.         double *pa;
    / J: S. `: g* l& a' }
34.         char keyname[sizeof(luVOID)];
    
35.         void *NowKey;
    
36.         k=m*n;
    
37.         pMatrix=(myMatrix *)GetBufObj(Matrix,keyname);//先尝试从缓冲区中获取一个矩阵对象
    0 l8 `2 v5 i: J7 X
38.         if(pMatrix)
    
39.         {
    
40.                 if(pMatrix->ArrayLen!=k)        //重置矩阵的大小
    + z% {7 Q: ^  s2 S$ u% b
41.                 {
    ; U1 P' B% Z! H
42.                         pa=new double[k];
    ! V2 b# `5 p+ n$ Z
43.                         if(!pa)
    
44.                         {
    , ?7 |% _& z$ F* K
45.                                 DeleteKey(keyname,sizeof(luVOID),Matrix,DelMatrix,1);        //将矩阵对象放回缓冲区
    " |+ Y6 D: w) E4 j9 c) H1 _
46.                                 return NULL;
    
47.                         }
      A2 A: D" @9 N! V& M3 n& K4 J/ F
48.                         delete[] pMatrix->Array;
    
49.                         pMatrix->Array=pa;
    ; H- t7 I0 T% i
50.                 }
    . f" E# \6 n9 J. y0 n
51.         }
    
52.         else
    
53.         {
    * ^, |" U' N% F3 r
54.                 pMatrix=new myMatrix;        //创建矩阵对象
    " H7 r5 n( {; s# b
55.                 if(!pMatrix) return NULL;
    
56.                 pMatrix->Array=new double[k];
    
57.                 if(!pMatrix->Array)
    ; D6 ]* V7 j8 [7 {
58.                 {
    
59.                         delete pMatrix;
    9 k6 E: M: W' M' a1 ]& m
60.                         return NULL;
    . l" D3 w+ i' u4 u$ q3 Z
61.                 }
    1 a- g3 B5 H1 o) H
62.                 if(InsertKey((char *)&pMatrix,-1,Matrix,pMatrix,DelMatrix,NULL,0,NowKey))        //将矩阵对象注册到Lu
    0 G6 _' Z' a0 D$ S6 r
63.                 {
    
64.                         delete pMatrix;
    
65.                         return NULL;
    / L8 b$ O7 ^; f6 X1 Z
66.                 }
    
67.         }
    
68.         pMatrix->ArrayLen=k; pMatrix->Dim[0]=m; pMatrix->Dim[1]=n;
    / a$ N8 K1 @. V& o
69.         return pMatrix;
    9 p- r4 Y- w8 w
70. }
    
71. LuData _stdcall OpMatrix(luINT mm,LuData *xx,void *hFor,int theOperator)        //运算符重载函数，用于LockKey函数
    
72. {
    
73.         LuData a;
    
74.         myMatrix *pMatrix1,*pMatrix2,*pMatrix3;
    
75.         luVOID i,j,k,m,n,u,v;
    
76.         double *pa,*pb,*pc;
    . X  e# R+ q4 _( p7 x* T
77.         luMessage pMessage;
    - ^! G. W8 |" i! o9 j
78.         wchar_t wchNum[32];
    
79.         char chNum[32];
    
80.         a.BType=luStaData_nil; a.VType=luStaData_nil; a.x=0;
    + A6 H0 B! y7 k
81.         switch(theOperator)
    ! }2 j2 `; ~$ X4 l. L2 S1 W
82.         {
    ( I* h3 B2 T& G+ D! z, s% `
83.         case 2:        //重载运算符*
    
84.                 pMatrix1=(myMatrix *)SearchKey((char *)&(xx->x),sizeof(luVOID),Matrix);
    
85.                 pMatrix2=(myMatrix *)SearchKey((char *)&((xx+1)->x),sizeof(luVOID),Matrix);
    1 }! y  S% Y, {% @& T5 D
86.                 if(!pMatrix1 || !pMatrix2) break;        //对象句柄无效，不是矩阵
    . A0 Z; p$ F( O6 }0 x
87.                 if(pMatrix1->Dim[1]!=pMatrix2->Dim[0]) break;        //维数不匹配
    
88.                 pMatrix3=NewMatrix(pMatrix1->Dim[0],pMatrix2->Dim[1]);        //生成新矩阵
      K( ?+ b; C/ A8 L* r
89.                 if(!pMatrix3) break;
    
90.                 pa=pMatrix1->Array; pb=pMatrix2->Array; pc=pMatrix3->Array;
    : f1 P& s  [; c& u: F
91.                 m=pMatrix1->Dim[0]; n=pMatrix1->Dim[1]; k=pMatrix2->Dim[1];
    7 \! u  }5 j, B; r
92.                 for(i=0; i<m; i++)        //矩阵乘
    
93.                 {
    
94.                         for(j=0; j<k; j++)
    . n# w# y! [9 ^& F; [9 Y
95.                         {
    7 H! Z: J) Z5 E
96.                                 u=i*k+j; pc[u]=0.0;
    ; Y, j# Z! u9 d
97.                                 for (v=0; v<n; v++)
    ' Z, j/ B% F! N
98.                                 {
    9 a) g' m' c, \; r6 X5 F
99.                                         pc[u]=pc[u]+pa[i*n+v]*pb[v*k+j];
    
100.                                 }
     
101.                         }
     
102.                 }
     ' a; H9 Y/ p, U: g9 L8 w/ ~/ O
103.                 FunReObj(hFor);        //告诉Lu，返回一个动态对象
     7 K: D; H6 S5 S& J( t. A
104.                 a.BType=Matrix; a.VType=Matrix; a.x=0; *(luVOID *)&(a.x)=(luVOID)pMatrix3;
     ! R$ q3 N: i# a
105.                 break;
     
106.         case 25:        //重载运算符.*
     0 Y8 H/ D4 n" a0 u" F3 h
107.                 pMatrix1=(myMatrix *)SearchKey((char *)&(xx->x),sizeof(luVOID),Matrix);
     
108.                 pMatrix2=(myMatrix *)SearchKey((char *)&((xx+1)->x),sizeof(luVOID),Matrix);
     * N2 F, u+ m% e1 Y1 z% V) t+ Z  M7 a
109.                 if(!pMatrix1 || !pMatrix2) break;        //对象句柄无效，不是矩阵
     
110.                 if(pMatrix1->Dim[0]!=pMatrix2->Dim[0] || pMatrix1->Dim[1]!=pMatrix2->Dim[1]) break;        //维数不相同
     7 c; T' _7 D2 |6 m
111.                 pMatrix3=NewMatrix(pMatrix1->Dim[0],pMatrix1->Dim[1]);        //生成新矩阵
     2 }7 A1 `! z0 D: z* r, p0 ~0 I
112.                 if(!pMatrix3) break;
     ) r& D3 B2 g* T2 _# [7 s
113.                 for(i=0;i<pMatrix1->ArrayLen;i++) pMatrix3->Array[i]=pMatrix1->Array[i]*pMatrix2->Array[i];        //矩阵点乘
     
114.                 FunReObj(hFor);        //告诉Lu，返回一个动态对象
     ) ?1 ]5 S# r5 {! |
115.                 a.BType=Matrix; a.VType=Matrix; a.x=0; *(luVOID *)&(a.x)=(luVOID)pMatrix3;
     
116.                 break;
     
117.         case 46:        //重载函数new
     8 c0 R8 {$ f, B; R* E
118.                 if(mm<2) break;
     # ]$ h  {8 a! X, l6 y2 D: d
119.                 if((xx+1)->x<1 || (xx+2)->x<1 || (xx+1)->BType!=luStaData_int64 || (xx+2)->BType!=luStaData_int64) break;
     - [# V( l3 {" b* y. h+ h! Z
120.                 pMatrix3=NewMatrix((luVOID)(xx+1)->x,(luVOID)(xx+2)->x);//生成新矩阵
     
121.                 if(!pMatrix3) break;
     : W/ N0 |: p$ [& M: k
122.                 for(j=0,i=3;i<=mm;i++,j++)        //赋初值
     + R* \+ Q* B+ e9 O& h. U; {0 m3 J
123.                 {
     
124.                         if(j>=pMatrix3->ArrayLen) break;
     
125.                         if((xx+i)->BType!=luStaData_double) break;        //只接受实数参数
     ; s7 w: F: I9 N# s! m- o
126.                         pMatrix3->Array[j]=*(double *)&((xx+i)->x);
     
127.                 }
     0 y6 d7 }7 w1 |( P
128.                 FunReObj(hFor);        //告诉Lu，返回一个动态对象
     1 B4 T! p2 O3 f' M- v
129.                 a.BType=Matrix; a.VType=Matrix; a.x=0; *(luVOID *)&(a.x)=(luVOID)pMatrix3;
     # }1 f$ R- N, e
130.                 break;
     
131.         case 49:        //重载函数o
     * I/ U/ B; o' V) [' P
132.                 pMessage=(luMessage)SearchKey("\0\0\0\0",sizeof(luVOID),luPubKey_User);
     
133.                 if(!pMessage) break;
     , @. e* I) V% E+ ?2 W8 t; `5 h
134.                 pMatrix1=(myMatrix *)SearchKey((char *)&(xx->x),sizeof(luVOID),Matrix);
     
135.                 if(!pMatrix1) break;        //对象句柄无效，不是矩阵
     % O+ |" p! \; d1 s' f5 E& x
136.                 pa=pMatrix1->Array;
     + E) c- r6 Y" g- N) r$ y" f
137.                 m=pMatrix1->Dim[0]; n=pMatrix1->Dim[1]; k=0;
     * [2 z, b2 {8 F) k& t4 Z
138.                 for(i=0; i<m; i++)        //输出矩阵
     % W8 F5 t( O9 R: C6 J: P3 }
139.                 {
     - N2 H7 s( p  v! B1 Y* q! O$ Q
140.                         pMessage(L"\r\n"); k+=2;
     2 v1 s; L& |& K' h
141.                         for(j=0; j<n; j++)
     
142.                         {
     , R; a5 i7 y9 N' S
143.                                 _gcvt_s(chNum,pa[i*n+j],16);
     + m. t- }+ X/ N. S, _4 \! \9 ?' d
144.                                 for(u=0;chNum[u];u++) {wchNum[u]=chNum[u]; k++;}
       ]! m; _4 f1 g0 \. f
145.                                 wchNum[u]='\0';
     
146.                                 pMessage(wchNum); pMessage(L"  "); k+=2;
     
147.                         }
     
148.                 }
     ' i4 m% k& a( S3 `/ s
149.                 pMessage(L"\r\n"); k+=2;
     
150.                 a.BType=luStaData_int64; a.VType=luStaData_int64; a.x=k;        //按函数o的要求，返回输出的字符总数
     + ]- ~) K  u: ~$ @" E
151.                 break;
     
152.         default:
     
153.                 break;
     
154.         }
     / `4 a% A+ Y( W, v$ R: a( q
155.         return a;
     4 `  y! H0 P$ e5 p" g! R- @9 b
156. }
     
157. void main(void)
     % L: m* b5 C' k' V  ~
158. {
     
159.         void *hFor;                //表达式句柄
     
160.         luINT nPara;                //存放表达式的自变量个数
     
161.         LuData *pPara;                //存放输入自变量的数组指针
     
162.         luINT ErrBegin,ErrEnd;        //表达式编译出错的初始位置和结束位置
     
163.         int ErrCode;                //错误代码
     * q* T& V* [( l. I5 a8 `5 G: ~
164.         void *v;
     
165.         wchar_t ForStr[]=L"o{new[matrix,2,3: 0.,1.,2.;3.,4.,5.]*new[matrix,3,2: 1.,2.;3.,4.;5.,6.]}";//字符串表达式，矩阵乘
     ) ?0 J1 Q$ i+ j8 f
166.         //wchar_t ForStr[]=L"o{new[matrix,2,3: 0.,1.,2.;3.,4.,5.].*new[matrix,2,3: 1.,2.,3.;4.,5.,6.]}";//字符串表达式，矩阵点乘
     5 |( m& _7 k- ?6 W; F2 G  A" v5 Q
167.         LuData Val;
     
168.         if(!InitLu()) return;        //初始化Lu
     
169.         while(LockKey(Matrix,DelMatrix,OpMatrix)){Matrix--;}        //锁定一个键，用于存储矩阵扩展类型
     ' t. R- A+ b) V1 Y: p; G
170. 
     " ~# D1 T7 `: a5 A4 M
171.         Val.BType=luStaData_int64; Val.VType=luStaData_int64; Val.x=Matrix;        //定义整数常量
     $ K1 Z" V6 \8 k! Z8 X
172.         SetConst(L"matrix",&Val);        //设置整数常量
     
173.         InsertKey("\0\0\0\0",4,luPubKey_User,LuMessage,NULL,NULL,1,v); //使Lu运行时可输出函数信息
     
174.         wcout.imbue(locale("chs"));        //设置输出的locale为中文
     
175.   
     2 U6 q7 {, K' y- ?9 ~$ m1 N
176.         ErrCode=LuCom(ForStr,0,0,0,hFor,nPara,pPara,ErrBegin,ErrEnd); //编译表达式
     
177.         if(ErrCode)
     
178.         {
     0 F- f) u9 z' \  b
179.                 wcout<<L"表达式有错误！错误代码："<<ErrCode<<endl;
     4 _3 K1 }' C  n2 h, ?1 C- R
180.         }
     ' S& `8 m) }& A' p; M
181.         else
     4 R' T) b! m6 t( K3 M
182.         {
       b: J( M+ `( i, b. d+ f. O
183.                 LuCal(hFor,pPara);        //计算表达式的值
       U0 @1 q; S- ~: h. m
184.         }
     
185.         LockKey(Matrix,NULL,OpMatrix);//解锁键Matrix，本例中，该函数可以不用
     8 ~" q, K  t! ?
186.         FreeLu();                        //释放Lu
     
187. }

复制代码

习题：

, a9 {0 M: B3 F3 W2 Z: R    （1）自定义矩阵的加、减、左除、右除、点左除等运算，自编测试字符串代码，重新编译运行程序，观察计算结果。
! v9 X" r" Y+ u+ ^6 g3 c
" v  C. o) j: t! O. w+ {! U    （2）小矩阵乘效率测试。编译运行以下Lu字符串代码：

1. main(:a,b,c,d,t,i)=
   
2.     a=new[matrix,2,2: 1.,2.,2.,1.],
   % ^/ L1 q0 j) Q- B1 o! r( Q
3.     b=new[matrix,2,2: 2.,1.,1.,2.],
   
4.     c=new[matrix,2,2: 2/3.,-1/3.,-1/3.,2/3.],
   7 ]( {0 b2 {8 z7 X4 Q4 q% z
5.     t=clock(),
   
6.     d=a*b, i=0, while{i<1000000, d=d*c*b, i++},
   
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: 1.,2.,2.,1.], b=new[matrix,2,2: 2.,1.,1.,2.], c=new[matrix,2,2: 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.
   ; D- D$ F+ C2 U  |, H  `7 _* {
2. 5. 4.
   6 B: g+ {4 `) _  |' `8 H
3. time=0.797 seconds.
   : K- K: ^/ r2 n& h! `! O5 d6 B
4. 请按任意键继续. . .

复制代码

Matlab 2009a 代码：

1. a=[1.,2.;2.,1.];
   
2. b=[2.,1.;1.,2.];
   
3. c=[2/3.,-1/3.;-1/3.,2/3.];
   ' q' u1 Z  @' G# O
4. tic,
   
5. d=a*b;
   
6. for i=1:1000000
   # l2 b: }( d* ]  q4 J
7.     d=d*c*b;
   
8. end
   
9. d,
   
10. toc

复制代码

结果：

1. d =
   6 l& ^( [: ?; e( p2 H2 I; N* K* ^+ \
2.      4     5
   
3.      5     4
   & ]4 X. f7 m- i2 j  t
4. Elapsed time is 2.903034 seconds.

复制代码

本例矩阵乘效率测试，Lu的速度超过了Matlab，主要在于Lu有更高的动态对象管理效率。
8 j4 d* E9 g, I, r9 j
    由以上可以看出，自定义数据类型和系统内置类型有近乎相同的效率。

---

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

Powered by Discuz! X2.5