数学建模社区-数学中国

标题: Lu基于系统内置对象创建扩展数据类型，小矩阵乘效率测试 [打印本页]
作者: forcal    时间: 2011-10-23 09:37
标题: Lu基于系统内置对象创建扩展数据类型，小矩阵乘效率测试
    本例中，我们将基于系统内置实数数组创建矩阵（matrix）类型，即：基本类型为luDynData_realarray（标识实数数组），扩展类型为matrix（标识矩阵）。为了简单，我们仅处理二维实数数组即矩阵类型。+ j* t2 b! u% P0 U; e
9 X; N/ z) X) r: q
    基本要点：* c( |7 Y) F- G# t$ Z
5 S% ]& y& j6 S% a6 M3 v- E8 J' ~$ ~
    （1）为扩展类型matrix编写运算符重载函数OpMatrix。
0 R" {7 e$ J7 H. ?! X; K5 n" w1 \8 B/ F' V( J2 i- N
    （2）用函数LockKey将重载函数OpMatrix注册到Lu，锁定的键的类型即为matrix，要注册为常量，以便于使用。
2 h/ ^9 F8 L  [6 V1 V! _' f5 u; }: o7 l, ~: c
    （3）为扩展类型matrix编写其他操作函数（本例未提供）。8 J0 C, d" N; N) q5 v3 P
0 |% H6 G4 ?! H
    （4）用函数LockKey解锁键matrix（本例中，程序退出时会自动解锁，故可以不用）。
  1. #include <windows.h>
    $ V: l. i& t" m
  2. #include <iostream>! Q, d( d4 ~' f7 m6 a+ x& i8 B
  3. #include <math.h>
    , \4 o, a. A3 p
  4. #include "lu32.h"
    ( N6 M( U& S; O

  5. - @: {1 y6 d5 o9 g) y, e
  6. #pragma comment( lib, "lu32.lib" )  V8 E( Q1 E0 ^6 D4 |# t

  7. 2 O: N0 Y7 D0 X9 _/ W& r& `. i3 v
  8. using namespace std;; Y) J3 L: }# n  Y( N
  9. ' j5 [( s5 {. l7 Z2 E6 L
  10. luKEY Matrix=-1000;        //标识矩阵类型，最终的Matrix由LockKey决定
    ; c' y0 [+ b+ u. L* J

  11. ( ?- [6 J$ `: Z' s! h& f
  12. void _stdcall LuMessage(wchar_t *pch)//输出动态库信息，该函数注册到Lu，由Lu二级函数调用 1 b- M( t" t$ B2 I" g3 ?( b5 E
  13. {6 v  s7 w; F& N+ ~
  14.         wcout<<pch;
    # a. Q5 }. g; o  j
  15. }
    6 n( Y. F2 z) y2 t0 c. ?/ X( F
  16. void _stdcall DelMatrix(void *me)        //用于LockKey函数，因为是基于系统内置实数数组创建矩阵，故该函数什么也不做, G/ i4 K8 m9 E7 S4 C
  17. {9 A" i; t6 E) L) @+ ?
  18. }6 S. N. f1 t: u8 |+ m( z
  19. LuData _stdcall OpMatrix(luINT mm,LuData *xx,void *hFor,int theOperator)        //运算符重载函数，用于LockKey函数- B3 _4 ?, Y/ Y, i) Y
  20. {" j8 i1 }+ J" Q' D( G
  21.         LuData a;, B5 s  c6 h* d8 E( w7 z2 `
  22.         luRealArray *pRealArray1,*pRealArray2,*pRealArray3;2 P3 r! z0 r3 K8 I' Z9 }2 A/ \
  23.         luVOID i,j,k,m,n,u,v;
    ' e! o' f, F: W) T# i& ^6 I: u
  24.         double *pa,*pb,*pc;0 K, m  q2 l3 g6 v  z" V
  25.         luMessage pMessage;
    & ?0 K3 U2 g0 |5 B. x1 Q
  26.         wchar_t wchNum[32];- }( d6 [& t* O3 _+ o8 j
  27.         char chNum[32];
    / c6 i7 S" I6 J! S
  28.         a.BType=luStaData_nil; a.VType=luStaData_nil; a.x=0;
    4 j) x/ U6 v) G* y# y% _
  29.         switch(theOperator)
    + T- U! V7 Y7 q# n
  30.         {
    " ~: b# x4 J- _" |, e  V6 k
  31.         case 2:        //重载运算符*! W" b, J: u: V! U8 c
  32.                 pRealArray1=(luRealArray *)SearchKey((char *)&(xx->x),sizeof(luVOID),luDynData_realarray);8 l" ~: K5 M/ C4 a1 n; |3 b9 Y! g4 e
  33.                 pRealArray2=(luRealArray *)SearchKey((char *)&((xx+1)->x),sizeof(luVOID),luDynData_realarray);9 V6 ]/ ]' a( J
  34.                 if(!pRealArray1 || !pRealArray2) break;        //对象句柄无效，不是实数数组: Q. V; p0 A6 h" X
  35.                 if(pRealArray1->DimLen!=2 || pRealArray2->DimLen!=2) break;        //不是二维实数数组（矩阵）
      v8 [0 k  J% f( z4 R
  36.                 if(pRealArray1->Dim[1]!=pRealArray2->Dim[0]) break;        //维数不匹配' v2 d+ M8 t. O6 s# F1 I
  37.                 pRealArray3=(luRealArray *)NewSysObj(luDynData_realarray,pRealArray1->Dim[0]*pRealArray2->Dim[1],2);        //创建矩阵对象5 z. e! S6 [+ y) k. b
  38.                 if(!pRealArray3) break;
    5 p0 v7 p# a% b: ~  }
  39.                 pRealArray3->Dim[0]=pRealArray1->Dim[0]; pRealArray3->Dim[1]=pRealArray2->Dim[1];        //设置矩阵维数大小  J9 s: U4 J0 j+ N3 a/ m% k
  40.                 pa=pRealArray1->Array; pb=pRealArray2->Array; pc=pRealArray3->Array;
    ( t- O' u: O* D& [& H2 X5 k) p
  41.                 m=pRealArray1->Dim[0]; n=pRealArray1->Dim[1]; k=pRealArray2->Dim[1];: P3 p6 d2 U3 {; o! ~, j; s8 A
  42.                 for(i=0; i<m; i++)        //矩阵乘' C  y  ]# I1 |9 J1 D( o# _
  43.                 {
    7 J+ f; I' R8 M; L$ B$ o
  44.                         for(j=0; j<k; j++)% X6 J, J3 q: |3 c
  45.                         {
    9 e' Y% f* X0 }9 u6 y' S- x% L/ h4 y
  46.                                 u=i*k+j; pc[u]=0.0;( r" E8 L8 H; m% w" Y& ^/ k
  47.                                 for (v=0; v<n; v++)
    / v0 X- @* S' g+ E2 H1 p
  48.                                 {
    1 H- ?0 r2 l, ^1 A+ s
  49.                                         pc[u]=pc[u]+pa[i*n+v]*pb[v*k+j];+ B) `  u2 v, o1 {; U/ p. I+ S
  50.                                 }2 }( ^! r: u: _, p9 T3 g/ E
  51.                         }
    ! E4 c& H3 N8 s  t0 r' K% l5 h0 e+ g
  52.                 }& i$ |; _" h  v: X; s" s. L- l
  53.                 FunReObj(hFor);        //告诉Lu，返回一个动态对象+ c) a! o+ V3 P  g8 b
  54.                 a.BType=luDynData_realarray; a.VType=Matrix; a.x=0; *(luVOID *)&(a.x)=(luVOID)pRealArray3;# N5 T& W" S" G& j" o% l. M( ?  W
  55.                 break;
    1 t; A$ ?/ Y, f* ]9 i0 [
  56.         case 25:        //重载运算符.*/ Z5 J+ q8 ~9 Z6 ?1 ~# F; M3 F
  57.                 pRealArray1=(luRealArray *)SearchKey((char *)&(xx->x),sizeof(luVOID),luDynData_realarray);
    ( R1 v0 Q- o( b( N7 H* ]6 A
  58.                 pRealArray2=(luRealArray *)SearchKey((char *)&((xx+1)->x),sizeof(luVOID),luDynData_realarray);
    9 a' M( f* ~8 p. B
  59.                 if(!pRealArray1 || !pRealArray2) break;        //对象句柄无效，不是实数数组5 W2 v9 D  B% z
  60.                 if(pRealArray1->DimLen!=2 || pRealArray2->DimLen!=2) break;        //不是二维实数数组（矩阵）' w, q7 o2 b4 |  N9 D
  61.                 if(pRealArray1->Dim[0]!=pRealArray2->Dim[0] || pRealArray1->Dim[1]!=pRealArray2->Dim[1]) break;        //维数不相同6 G4 h/ K$ |- d
  62.                 pRealArray3=(luRealArray *)NewSysObj(luDynData_realarray,pRealArray1->ArrayLen,2);        //创建矩阵对象5 A3 u+ E$ ^7 w- [, p5 K* J
  63.                 if(!pRealArray3) break;
    # _( D6 k, g1 r7 O
  64.                 pRealArray3->Dim[0]=pRealArray1->Dim[0]; pRealArray3->Dim[1]=pRealArray1->Dim[1];        //设置矩阵维数大小7 |# k4 ?: a( z
  65.                 for(i=0;i<pRealArray1->ArrayLen;i++) pRealArray3->Array[i]=pRealArray1->Array[i]*pRealArray2->Array[i];//矩阵点乘
    " S1 o. s; U& I, `4 g4 ^; }* H
  66.                 FunReObj(hFor);        //告诉Lu，返回一个动态对象9 |& ^$ O4 q! ?( g
  67.                 a.BType=luDynData_realarray; a.VType=Matrix; a.x=0; *(luVOID *)&(a.x)=(luVOID)pRealArray3;+ }- O  W5 N! \- c$ J
  68.                 break;
    6 e2 r" ^) t% n* @5 r) w% H8 ^; v
  69.         case 46:        //重载函数new& ]8 o: r, n. F' T1 k' J
  70.                 a=ExeOperator(mm,xx,hFor,theOperator,luDynData_realarray);        //直接调用基本类型luDynData_realarray的new函数/ v" b' g# Z% Z2 j  N& ]1 o- U
  71.                 if(a.VType==luDynData_realarray) a.VType=Matrix;        //设置扩展类型为自定义的Matrix类型' |( k, R% p* P& P  G9 v& g
  72.                 break;4 q0 T/ _! h& F0 H2 ~, ]; h2 i
  73.         case 49:        //重载函数o
    ; w+ k3 Z* Q' i; r3 P
  74.                 pMessage=(luMessage)SearchKey("\0\0\0\0",sizeof(luVOID),luPubKey_User);- D  u) I1 u8 k$ X/ S: ~* }
  75.                 if(!pMessage) break;
    % `0 y, ^7 d5 s' _9 y4 r# r$ N
  76.                 pRealArray1=(luRealArray *)SearchKey((char *)&(xx->x),sizeof(luVOID),luDynData_realarray);$ [) O! `8 y. W9 r
  77.                 if(!pRealArray1) break;        //对象句柄无效，不是实数数组5 ^; Y( t2 _( ]0 `+ N/ v0 {5 u
  78.                 if(pRealArray1->DimLen!=2) break;        //不是二维实数数组（矩阵）# |  e" \0 ?, l1 Q2 N& _
  79.                 pa=pRealArray1->Array;
    " j0 {/ \* n/ Q2 l
  80.                 m=pRealArray1->Dim[0]; n=pRealArray1->Dim[1]; k=0;
    + ^6 I' s0 Q- Z
  81.                 for(i=0; i<m; i++)        //输出矩阵4 i8 m5 K" d2 z4 U- q6 E4 G
  82.                 {
      Y2 i2 N6 m, z
  83.                         pMessage(L"\r\n"); k+=2;
    0 O: T5 d- M) j- m9 P2 F
  84.                         for(j=0; j<n; j++)! Y8 V' _* D& j" I2 c( G- w
  85.                         {
    7 A% @+ H$ m: i, ~! J
  86.                                 _gcvt_s(chNum,pa[i*n+j],16);
    * D# a' H: q; g- w1 @
  87.                                 for(u=0;chNum[u];u++) {wchNum[u]=chNum[u]; k++;}4 r( C2 R( y8 d2 r+ I, ?. b( k
  88.                                 wchNum[u]='\0';
    0 a7 ^0 i5 C$ ?) \
  89.                                 pMessage(wchNum); pMessage(L"  "); k+=2;
    ; Z9 [. D0 g2 X) K2 b5 p7 I
  90.                         }; ^4 [1 N- e# K& k
  91.                 }% ~0 h, C% Y$ |2 F$ `* @) g
  92.                 pMessage(L"\r\n"); k+=2;  V8 C6 e* b# e& W2 U7 z
  93.                 a.BType=luStaData_int64; a.VType=luStaData_int64; a.x=k;        //按函数o的要求，返回输出的字符总数. _% D; K+ V6 f3 Z7 a, _
  94.                 break;
    / x1 |4 m6 j3 ?/ T: u1 N
  95.         default:# L5 ?# ]8 }7 y5 u2 c. q6 o( j
  96.                 break;
    5 H( Y( ]7 w- B7 m, E* `% T
  97.         }+ u: Y$ X# a4 S2 _/ k
  98.         return a;
    7 B  [. E+ K0 [8 g' x7 s
  99. }
    ) O: V' p3 W) o: d" R9 G& l
  100. void main(void)* Z7 j; ^; g6 g, ]2 i
  101. {  F" m- M- ~1 g2 u  I! u5 U: W* }
  102.         void *hFor;                //表达式句柄- y# U) L& b4 m! `# M
  103.         luINT nPara;                //存放表达式的自变量个数
      Y4 @6 _; D7 T) Z& n. q
  104.         LuData *pPara;                //存放输入自变量的数组指针. Y' m6 B9 D5 f3 F* l
  105.         luINT ErrBegin,ErrEnd;        //表达式编译出错的初始位置和结束位置2 [8 u: E8 y0 ~
  106.         int ErrCode;                //错误代码
    * q' o8 @! @; d
  107.         void *v;
    7 H) A3 Z8 n7 x( _1 D& `5 q
  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.]}";//字符串表达式，矩阵乘2 J, t% A* u1 i& U
  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.]}";//字符串表达式，矩阵点乘
    1 A; L1 c; m/ G6 l
  110.         LuData Val;
    ( M7 k2 O5 T) |# l
  111.         if(!InitLu()) return;        //初始化Lu* ^& H$ l) j' z+ i! C9 w  Q
  112.         while(LockKey(Matrix,DelMatrix,OpMatrix)){Matrix++;}        //锁定一个键，用于存储矩阵扩展类型! H: F7 l7 w8 ?4 K# f

  113. ( W* `0 S! }+ a3 q- h
  114.         Val.BType=luStaData_int64; Val.VType=luStaData_int64; Val.x=Matrix;        //定义整数常量
    4 x( l) ?- X! D4 U( c
  115.         SetConst(L"matrix",&Val);        //设置整数常量
    & J2 _1 |7 q& K. A! k2 u$ _  q7 _! g& N
  116.         InsertKey("\0\0\0\0",4,luPubKey_User,LuMessage,NULL,NULL,1,v); //使Lu运行时可输出函数信息' q, b; g4 T; q( N! P
  117.         wcout.imbue(locale("chs"));        //设置输出的locale为中文
    7 w+ L4 J" Z! Q- R. i4 ?8 b3 ]  i0 s: y0 z
  118.   
    7 \8 F: l* d$ X9 j0 W4 o& i2 z' j
  119.         ErrCode=LuCom(ForStr,0,0,0,hFor,nPara,pPara,ErrBegin,ErrEnd); //编译表达式1 m1 Z7 ~  |/ J  k6 J( M/ t- D
  120.         if(ErrCode). n# t8 G. ~( t( P- |) e$ i
  121.         {# @4 A% B% ~2 r7 ?
  122.                 wcout<<L"表达式有错误！错误代码："<<ErrCode<<endl;
    4 f' Z# F& `) a# _% I8 H
  123.         }& h: p* w( X4 E* g, M& ]4 I' ]' R
  124.         else
    5 X/ s% y! S2 Y% B+ C# C
  125.         {/ u% w* [% E+ O2 S
  126.                 LuCal(hFor,pPara);        //计算表达式的值7 P. @( X1 U6 H$ w
  127.         }7 ~# ~; c2 H3 e$ g7 B9 s
  128.         LockKey(Matrix,NULL,OpMatrix);//解锁键Matrix，本例中，该函数可以不用
    ) h0 N/ k+ @5 o. }, B
  129.         FreeLu();                        //释放Lu+ m- |& V% b; T
  130. }
复制代码
习题：
9 @5 _' A2 G3 n- F* @: H
0 k+ T- E$ u# T, N3 W/ s    （1）自定义矩阵的加、减、左除、右除、点左除等运算，自编测试字符串代码，重新编译运行程序，观察计算结果。
: N# p  }, q' Y( z0 p4 Z7 ~' Q! g7 L" I" v
    （2）小矩阵乘效率测试。编译运行以下Lu字符串代码：
  1. main(:a,b,c,d,t,i)=
    6 i% _  l; r8 E. O0 K: r
  2.     a=new[matrix,2,2,data:1.,2.,2.,1.],
    2 }$ g5 V  t6 ~5 w) B3 A. ]
  3.     b=new[matrix,2,2,data:2.,1.,1.,2.],2 M  m  |" o' H  t+ g8 @9 A8 F
  4.     c=new[matrix,2,2,data:2/3.,-1/3.,-1/3.,2/3.],9 }% D& t: Y0 R* X: q; l8 Q
  5.     t=clock(),/ g6 h$ S8 K5 V5 J% j; r& E
  6.     d=a*b, i=0, while{i<1000000, d=d*c*b, i++},' g4 q) W/ ?! F: _  K, f0 \$ m
  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.
    , k7 [( r1 B8 ?/ u/ n8 Q
  2. 5. 4.
    * h. J- t( Y; Y
  3. time=0.875 seconds.; C6 s9 E# \' n% @9 H3 U7 e
  4. 请按任意键继续. . .
复制代码
Matlab 2009a 代码：
  1. a=[1.,2.;2.,1.];! I9 [0 {0 ^/ r, s
  2. b=[2.,1.;1.,2.];
    * N$ w0 p# G! Z2 k% p
  3. c=[2/3.,-1/3.;-1/3.,2/3.];7 l. h% u( S2 T; g1 P
  4. tic,
    * @2 V5 H/ C' ?# Z7 ?$ _9 L7 S! g
  5. d=a*b;
    & r4 Q* Q! i0 w2 `% c5 {$ A6 m/ O
  6. for i=1:10000009 ]. v* J, l0 t8 i& S  I& M9 d& L4 E
  7.     d=d*c*b;; ~! {3 ^/ H3 f6 H7 X
  8. end. \" a# ]1 @# ]  F
  9. d,
    3 m( O) [* v' ~7 b9 r1 M. g- w
  10. toc
复制代码
结果：
  1. d =
    9 ^4 s5 a9 c$ |$ Y3 m) H# |/ e
  2.      4     5
    / A7 U/ v4 S- `/ e  A5 v! n
  3.      5     48 K. ~: i; Z0 A3 o" p
  4. Elapsed time is 2.903034 seconds.
复制代码
本例矩阵乘效率测试，Lu的速度超过了Matlab，主要在于Lu有更高的动态对象管理效率。
% T1 t" Y2 \: j$ J



欢迎光临 数学建模社区-数学中国 (http://www.madio.net/) Powered by Discuz! X2.5