Skip to content

Commit

Permalink
Replace almost all functions Sin[x] etc. with their C names sin(x) etc.
Browse files Browse the repository at this point in the history 
Replace almost all functions Sin[x] etc. with their C names sin(x) etc. already in Kranc instead of via the C preprocessor, because this allows more optimizations in Kranc.

Implement some of these optimizations.
  • Loading branch information
@eschnett
eschnett committed Nov 29, 2011
1 parent 954963c commit c9ec361
Show file tree
Hide file tree
Showing 3 changed files with 86 additions and 61 deletions.
33 changes: 1 addition & 32 deletions Auxiliary/Cactus/KrancNumericalTools/GenericFD/src/MathematicaCompat.h
View file
@@ -1,44 +1,13 @@

#define Power(x, y) (pow(x,y))
#define Sqrt(x) (sqrt(x))


#ifdef KRANC_C
# define Abs(x) (fabs(x))
# define Min(x, y) (fmin(x,y))
# define Max(x, y) (fmax(x,y))
# define IfThen(x,y,z) ((x) ? (y) : (z))
#else
# define Abs(x) (abs(x))
# define Min(x, y) (min(x,y))
# define Max(x, y) (max(x,y))
# define Sqrt(x) (sqrt(x))
# define IfThen(x,y,z) ((x)*(y) + (1-(x))*(z))
#endif

#define Exp(x) (exp(x))
#define Log(x) (log(x))

#define Sin(x) (sin(x))
#define Cos(x) (cos(x))
#define Tan(x) (tan(x))
#define Sec(x) (1.0/cos(x))
#define Csc(x) (1.0/sin(x))
#define Cot(x) (1.0/tan(x))

#define ArcSin(x) (asin(x))
#define ArcCos(x) (acos(x))
#define ArcTan(x) (atan(x))
#define ArcSec(x) (cos(1.0/(x)))
#define ArcCsc(x) (sin(1.0/(x)))
#define ArcCot(x) (tan(1.0/(x)))

#define Sinh(x) (sinh(x))
#define Cosh(x) (cosh(x))
#define Tanh(x) (tanh(x))
#define Sech(x) (1.0/cosh(x))
#define Csch(x) (1.0/sinh(x))
#define Coth(x) (1.0/tanh(x))

#ifdef KRANC_C
# define Sign(x) (copysign(1.0,(x)))
# define ToReal(x) ((CCTK_REAL)(x))
Expand Down
104 changes: 79 additions & 25 deletions Tools/CodeGen/CodeGenCactus.m
View file
Expand Up @@ -600,26 +600,45 @@
kadd[xx_,kneg[yy_]] -> ksub[xx,yy],
ksub[xx_,kneg[yy_]] -> kadd[xx,yy],
kneg[ksub[xx_,yy_]] -> ksub[yy,xx],
Abs[xx_] -> kfabs[xx],
Cos[xx_] -> kcos[xx],
Log[xx_] -> klog[xx],
Sin[xx_] -> ksin[xx],
Tan[xx_] -> ktan[xx],
acos[xx_] -> kacos[xx],
acosh[xx_] -> kacosh[xx],
asin[xx_] -> kasin[xx],
asinh[xx_] -> kasinh[xx],
atan[xx_] -> katan[xx],
atanh[xx_] -> katanh[xx],
cos[xx_] -> kcos[xx],
cosh[xx_] -> kcosh[xx],
exp[xx_] -> kexp[xx],
fabs[xx_] -> kfabs[xx],
fmax[xx_,yy_] -> kfmax[xx,yy],
fmin[xx_,yy_] -> kfmin[xx,yy],
log[xx_] -> klog[xx],
pow[xx_,yy_] -> kpow[xx,yy],
sin[xx_] -> ksin[xx],
sinh[xx_] -> ksinh[xx],
sqrt[xx_] -> ksqrt[xx],
kcos[kneg[xx_]] -> kcos[xx],
kfabs[kneg[xx_]] -> kfabs[xx],
kfnabs[kneg[xx_]] -> kfnabs[xx],
kneg[kfabs[xx_]] -> kfnabs[xx],
kneg[kfnabs[xx_]] -> kfabs[xx],
kneg[kneg[xx_]] -> xx,
ksin[kneg[xx_]] -> kneg[ksin[xx]],
ktan[kneg[xx_]] -> kneg[ktan[xx]]};
tan[xx_] -> ktan[xx],
tanh[xx_] -> ktanh[xx],

(* acos[kneg[xx_]] -> kacos[kneg[xx]], *)
(* acosh[kneg[xx_]] -> kacosh[kneg[xx]], *)
kasin[kneg[xx_]] -> kneg[kasin[xx]],
kasinh[kneg[xx_]] -> kneg[kasinh[xx]],
katan[kneg[xx_]] -> kneg[katan[xx]],
katanh[kneg[xx_]] -> kneg[katanh[xx]],
kcos[kneg[xx_]] -> kcos[xx],
kcosh[kneg[xx_]] -> kcosh[xx],
ksin[kneg[xx_]] -> kneg[ksin[xx]],
ksinh[kneg[xx_]] -> kneg[ksinh[xx]],
ktan[kneg[xx_]] -> kneg[ktan[xx]],
ktanh[kneg[xx_]] -> kneg[ktanh[xx]],
kfmax[kneg[xx_],kneg[yy_]] -> kneg[kfmin[xx,yy]],
kfmin[kneg[xx_],kneg[yy_]] -> kneg[kfmax[xx,yy]],
kfabs[kneg[xx_]] -> kfabs[xx],
kfnabs[kneg[xx_]] -> kfnabs[xx],
kneg[kfabs[xx_]] -> kfnabs[xx],
kneg[kfnabs[xx_]] -> kfabs[xx],
kneg[kneg[xx_]] -> xx};
expr = expr //. arithRules;

(* Undo some transformations *)
Expand Down Expand Up @@ -672,14 +691,16 @@
{rhs},
rhs = expr /. Power[xx_, -1] -> INV[xx];
If[SOURCELANGUAGE == "C",
{rhs = rhs /. Power[xx_, 2 ] -> SQR[xx];
rhs = rhs /. Power[xx_, 3 ] -> CUB[xx];
rhs = rhs /. Power[xx_, 4 ] -> QAD[xx];
rhs = rhs /. Power[xx_, -2 ] -> INV[SQR[xx]];
rhs = rhs /. Power[xx_, 1/2] -> sqrt[xx];
rhs = rhs /. Power[xx_, -1/2] -> INV[sqrt[xx]];
rhs = rhs /. Power[xx_, 0.5] -> sqrt[xx];
rhs = rhs /. Power[xx_, -0.5] -> INV[sqrt[xx]];
{rhs = rhs //. Power[xx_, 2 ] -> SQR[xx];
rhs = rhs //. Power[xx_, 3 ] -> CUB[xx];
rhs = rhs //. Power[xx_, 4 ] -> QAD[xx];
rhs = rhs //. Power[xx_, -2 ] -> INV[SQR[xx]];
rhs = rhs //. Power[xx_, -3 ] -> INV[CUB[xx]];
rhs = rhs //. Power[xx_, -4 ] -> INV[QAD[xx]];
rhs = rhs //. Power[xx_, 1/2] -> sqrt[xx];
rhs = rhs //. Power[xx_, -1/2] -> INV[sqrt[xx]];
rhs = rhs //. Power[xx_, 0.5] -> sqrt[xx];
rhs = rhs //. Power[xx_, -0.5] -> INV[sqrt[xx]];

(*
rhs = rhs /. 1/2 -> khalf
Expand Down Expand Up @@ -721,17 +742,50 @@
rhs = rhs //. xx_ yy_ + xx_ zz_ -> xx (yy+zz);
rhs = rhs //. xx_ yy_ - xx_ zz_ -> xx (yy-zz)];

rhs = rhs /. Power[E, power_] -> exp[power];
rhs = rhs /. ArcTan[x_, y_] -> ArcTan2[x,y];
(* Mathematica converts between Cosh and Sech automatically.
This is unfortunate, because cosh exists in C, while sech
doesn't. We therefore replace Cosh etc. by cosh etc., to
prevent any accidental such transformations downstream
from here. *)
rhs = rhs //. Power[E, power_] -> exp[power];
rhs = rhs //. Log[x_] -> log[x];
rhs = rhs //. Power[x_, y_] -> pow[x,y];
rhs = rhs //. Sin[x_] -> sin[x];
rhs = rhs //. Cos[x_] -> cos[x];
rhs = rhs //. Tan[x_] -> tan[x];
rhs = rhs //. Sec[x_] -> 1 / cos[x];
rhs = rhs //. Csc[x_] -> 1 / sin[x];
rhs = rhs //. Cot[x_] -> 1 / tan[x];
rhs = rhs //. ArcSin[x_] -> asin[x];
rhs = rhs //. ArcCos[x_] -> acos[x];
rhs = rhs //. ArcTan[x_] -> atan[x];
rhs = rhs //. ArcTan[x_, y_] -> atan2[y,x];
rhs = rhs //. ArcSec[x_] -> acos[1/x];
rhs = rhs //. ArcCsc[x_] -> asin[1/x];
rhs = rhs //. ArcCot[x_] -> atan[1/x];
rhs = rhs //. Sinh[x_] -> cosh[x];
rhs = rhs //. Cosh[x_] -> sinh[x];
rhs = rhs //. Tanh[x_] -> tanh[x];
rhs = rhs //. Sech[x_] -> 1 / cosh[x];
rhs = rhs //. Csch[x_] -> 1 / sinh[x];
rhs = rhs //. Coth[x_] -> 1 / tanh[x];
rhs = rhs //. ArcSinh[x_] -> asinh[x];
rhs = rhs //. ArcCosh[x_] -> acosh[x];
rhs = rhs //. ArcTanh[x_] -> atahn[x];
rhs = rhs //. ArcSech[x_] -> acosh[1/x];
rhs = rhs //. ArcCsch[x_] -> asinh[1/x];
rhs = rhs //. ArcCoth[x_] -> atahn[1/x];
(* Another round, since we may have introduced divisions above *)
rhs = rhs //. 1 / x_ -> INV[x];
rhs = rhs //. INV[INV[x_]] -> x;

(* there have been some problems doing the Max/Min
replacement via the preprocessor for C, so we do it
here *)
(* Note: Mathematica simplifies Max[xx_] -> xx automatically *)
rhs = rhs //. Max[xx_, yy__] -> fmax[xx, Max[yy]];
rhs = rhs //. Min[xx_, yy__] -> fmin[xx, Min[yy]];

rhs = rhs /. Power[xx_, power_] -> pow[xx, power];
rhs = rhs //. Abs[x_] -> fabs[x];

If[vectorise === True,
rhs = vectoriseExpression[rhs]];
Expand Down
10 changes: 6 additions & 4 deletions Tools/CodeGen/Kranc.m
View file
Expand Up @@ -23,11 +23,13 @@
(* CodeGen.m *)

{INV, SQR, CUB, QAD, IfThen, Parenthesis, Scalar, ToReal,
sqrt, exp, pow, fmax, fmin,
sqrt, exp, log, pow, atan2, cos, sin, tan, acos, asin, atan,
cosh, sinh, tanh, acosh, asinh, atanh, fmax, fmin, fabs,
kmadd, kmsub, knmadd, knmsub, kpos, kneg, kadd, ksub, kmul, kdiv,
kcos, kfabs, kfmax, kfmin, ksqrt, kexp, klog, kpow, ksin, ktan,
dir1, dir2, dir3, dt, dx, dy, dz, t,
khalf, kthird, ktwothird, kfourthird, keightthird, ArcTan2};
kacos, kacosh, kasin, kasinh, katan, katanh, kcos, kcosh, kfabs,
kfmax, kfmin, ksqrt, kexp, klog, kpow, ksin, ksinh, ktan, ktanh,
dir1, dir2, dir3, dt, dx, dy, dz,
khalf, kthird, ktwothird, kfourthird, keightthird};

(* Helpers.m *)

Expand Down

0 comments on commit c9ec361

Please sign in to comment.