|
Warning, the protected name MathML has been redefined and
unprotected
| |
最初の例は、どのようにしてパッケージがベクトル値関数を作成するかを示します。
>
|
f := 1-exp(-t)*x*y:
g := x-y*exp(-t)*x^2:
v := vector([f,g]);
|
| (2.1) |
>
|
fg := makeproc(v,parameters=[t,x,y]);
|
| (2.2) |
subroutine fg(t,x,y,crea_par)
doubleprecision t
doubleprecision x
doubleprecision y
doubleprecision crea_par(2)
doubleprecision t1
doubleprecision t5
doubleprecision v(2)
t1 = exp(-t)
v(1) = 1-t1*x*y
t5 = x**2
v(2) = x-y*t1*t5
crea_par(1) = v(1)
crea_par(2) = v(2)
return
return
end
| |
下に示す関数 f は、点 (x, y) を (C, D) だけ 平行移動した後、原点を中心に theta ラジアン回転移動させた点から、中心 (h, k)、半径 r の円への距離の2乗を計算します。関数 f の C, D, theta に関する勾配を計算し、結果として C言語のソースコードを生成するために自動微分を用います。そして、最適化されたコードの計算量を出力します。
>
|
f := proc(C,D,theta,x,y,h,k,r)
local s,c,xbar,ybar,d1,d2,d;
s := sin(theta);
c := cos(theta);
xbar := (x+C)*c + (y+D)*s;
ybar := (y+D)*c - (x+C)*s;
d1 := (h-xbar)^2;
d2 := (k-ybar)^2;
d := sqrt( d1+d2 ) - r;
d^2
end proc:
G := GRADIENT(f,[C,D,theta],result_type=array);
|
| (2.3) |
>
|
C(G,optimized, precision=single);
|
#include <math.h>
void G(C,D,theta,x,y,h,k,r,grd)
float C;
float D;
float theta;
float x;
float y;
float h;
float k;
float r;
float grd[3];
{
float c;
float d1;
float d2;
float dfr0[7];
float s;
float t1;
float t10;
float t19;
float t21;
float t3;
float t7;
float t8;
{
s = sin(theta);
c = cos(theta);
t1 = x+C;
t3 = y+D;
t7 = h-t1*c-t3*s;
d1 = t7*t7;
t8 = k-t3*c+t1*s;
d2 = t8*t8;
t10 = sqrt(d1+d2);
dfr0[6] = 2.0*t10-2.0*r;
dfr0[5] = dfr0[6]/t10/2.0;
dfr0[4] = dfr0[5];
dfr0[3] = -2.0*dfr0[4]*t8;
dfr0[2] = -2.0*dfr0[4]*t7;
t19 = dfr0[3];
t21 = dfr0[2];
dfr0[1] = t19*t3+t21*t1;
dfr0[0] = -t19*t1+t21*t3;
grd[0] = 0.0;
grd[1] = t19*c+t21*s;
grd[2] = -dfr0[1]*s+dfr0[0]*c;
return;
}
}
| |
| (2.4) |
この例では、1+x+x^2/2+x^3/6+...+x^n/n! を計算します。まず、記号的な和を用います。C 言語のソースコードを得るために、明示的に記号的な和を for ループに変換します。
>
|
f := Sum(x^i/i!,i=0..n);
|
| (2.5) |
>
|
F := makeproc(f,parameters=[n,x],locals=[i]);
|
| (2.6) |
| (2.7) |
>
|
F := declare(n::integer,F);
|
| (2.8) |
double F(n,x)
int n;
double x;
{
double i;
int i1;
double s1;
double t1;
{
if( n < 0.0 )
s1 = 0.0;
else
{
t1 = 1.0;
s1 = 1.0;
for(i1 = 1;i1 <= n;i1++)
{
t1 = x/i1*t1;
s1 += t1;
}
}
return(s1);
}
}
| |
2 番目のアプローチでは、プログラムは有限和を計算するために繰り返しを用います。codegen パッケージにおいて、提供されるプログラムのための中間表現を用います。この表現は式木です。これは intrep2maple コマンドを用いて、望むなら Fortran または C 言語のコードに変換できる Maple の手続きへ変換することができます。
>
|
f := Proc( Parameters(n::integer,x::float),
Locals(i,s,t),
StatSeq( Assign(s,1.0),
Assign(t,1.0),
For( i, 1, 1, n, true,
StatSeq( Assign(t,x*t/i), Assign(s,s+t) ) ),
s ));
|
| (2.9) |
| (2.10) |
この例において、単純な関数の勾配とヘシアンを計算し、2 つの手続きを操作し 1 つにまとめます。
| (2.11) |
>
|
F := makeproc(f,[x,y,t]);
|
| (2.12) |
>
|
G := GRADIENT(F, [x,y], result_type=array);
|
| (2.13) |
>
|
G := dontreturn(grd, makeparam(grd,G));
|
| (2.14) |
>
|
H := HESSIAN(F, [x,y], result_type=array):
H := renamevar(grd=hes, dontreturn(grd, makeparam(grd,H)) );
|
| (2.15) |
>
|
GH := optimize( joinprocs([G,H]) );
|
| (2.16) |
最後の例は maple2intrep および intrep2maple コマンドの機能を示します。maple2intrep コマンドは、Maple の手続きを操作に適した中間表現に変換します。コード中の関数が評価されるので、この表現がどう評価されているか注意しなければいけません。
>
|
f := proc(x,n) local A,i;
A := array(0..n);
A[0] := 1;
for i to n do A[i] := x*A[i-1] end do;
A
end proc:
IR := maple2intrep(f);
|
| (2.17) |
>
|
intrep2maple(IR); # array(0..n) is evaluated, giving an error
|
Error, non-integer ranges in array/table creation
| |
>
|
intrep2maple(eval(IR,1));
|
| (2.18) |
|