aboutsummaryrefslogtreecommitdiff
path: root/src/gr
diff options
context:
space:
mode:
authorjthorn <jthorn@f88db872-0e4f-0410-b76b-b9085cfa78c5>2002-04-22 15:10:32 +0000
committerjthorn <jthorn@f88db872-0e4f-0410-b76b-b9085cfa78c5>2002-04-22 15:10:32 +0000
commit74ef815aab12c0f9602b22b4f0542ca7c3a962be (patch)
treed6ab428f9b7bf9f14c19f981f929a08087bad310 /src/gr
parentca337a6176d82455b04e3552079f389852f49813 (diff)
change names alpha --> xcos, beta --> ycos, gamma --> zcos
because maple predefines gamma = Euler's constant git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinAnalysis/AHFinderDirect/trunk@575 f88db872-0e4f-0410-b76b-b9085cfa78c5
Diffstat (limited to 'src/gr')
-rw-r--r--src/gr/ellipsoid.c20
-rw-r--r--src/gr/ellipsoid.maple7
-rw-r--r--src/gr/ellipsoid.out196
3 files changed, 111 insertions, 112 deletions
diff --git a/src/gr/ellipsoid.c b/src/gr/ellipsoid.c
index b600b76..49a6457 100644
--- a/src/gr/ellipsoid.c
+++ b/src/gr/ellipsoid.c
@@ -4,25 +4,25 @@ fp t55;
t1 = a*a;
t2 = b*b;
t3 = t1*t2;
- t5 = t3*0.5772156649015329*CW;
+ t5 = t3*zcos*CW;
t6 = c*c;
t7 = t1*t6;
- t9 = t7*beta*BV;
+ t9 = t7*ycos*BV;
t10 = t2*t6;
- t12 = t10*alpha*AU;
- t28 = alpha*alpha;
+ t12 = t10*xcos*AU;
+ t28 = xcos*xcos;
t30 = CW*CW;
t33 = BV*BV;
t35 = t10*t28;
- t36 = beta*beta;
+ t36 = ycos*ycos;
t40 = AU*AU;
t42 = t7*t36;
- t43 = 0.5772156649015329*0.5772156649015329;
+ t43 = zcos*zcos;
t48 = t3*t43;
- t49 = -2.0*t1*0.5772156649015329*CW*beta*BV-2.0*t2*0.5772156649015329*CW*
-alpha*AU-2.0*t6*beta*BV*alpha*AU+t2*t28*t30+t6*t28*t33-t35+t1*t36*t30+t6*t36*
-t40-t42+t1*t43*t33+t2*t43*t40-t48;
+ t49 = -2.0*t1*zcos*CW*ycos*BV-2.0*t2*zcos*CW*xcos*AU-2.0*t6*ycos*BV*xcos*
+AU+t2*t28*t30+t6*t28*t33-t35+t1*t36*t30+t6*t36*t40-t42+t1*t43*t33+t2*t43*t40-
+t48;
t52 = sqrt(-t3*t6*t49);
- t55 = 1/(t48+t42+t35);
+ t55 = 1/(t35+t42+t48);
r_plus = (t5+t9+t12+t52)*t55;
r_minus = (t5+t9+t12-t52)*t55;
diff --git a/src/gr/ellipsoid.maple b/src/gr/ellipsoid.maple
index c9b0ae4..0c485fa 100644
--- a/src/gr/ellipsoid.maple
+++ b/src/gr/ellipsoid.maple
@@ -6,10 +6,11 @@
# angular coordinate system has center (U,V,W)
#
# direction cosines wrt angular coordinate center are (alpha,beta,gamma)
-# i.e. a point has coordinates (U+alpha*r, V+beta*r, W+gamma*r)
+# but Maple predefines gamma = Euler's constant, so we use (xcos,ycos,zcos)
+# instead, i.e. a point has coordinates (U+xcos*r, V+ycos*r, W+zcos*r)
#
# then the equation of the ellipsoid is
-# (U+alpha*r - A)^2 (V+beta*r - B)^2 (W+gamma*r - C)^2
+# (U+xcos*r - A)^2 (V+ycos*r - B)^2 (W+zcos*r - C)^2
# ----------------- + ---------------- + ----------------- = 1
# a^2 b^2 c^2
#
@@ -18,7 +19,7 @@
# BV = B - V
# CW = C - W
#
-eqn := (alpha*r - AU)^2/a^2 + (beta*r - BV)^2/b^2 + (gamma*r - CW)^2/c^2 = 1;
+eqn := (xcos*r - AU)^2/a^2 + (ycos*r - BV)^2/b^2 + (zcos*r - CW)^2/c^2 = 1;
read "../maple/util.mm";
read "../maple/codegen2.mm";
diff --git a/src/gr/ellipsoid.out b/src/gr/ellipsoid.out
index 27b5f21..2e77483 100644
--- a/src/gr/ellipsoid.out
+++ b/src/gr/ellipsoid.out
@@ -4,17 +4,18 @@
<____ ____> Waterloo Maple Inc.
| Type ? for help.
# ellipsoid.maple -- compute equations for offset ellipsoid setup
-# $Header$
+# $Header: /numrelcvs/ThornburgCVS/AHFinderDirect/src/gr/ellipsoid.maple,v 1.1 2002/04/22 11:02:27 jthorn Exp $
>
#
# ellipsoid has center (A,B,C), radius (a,b,c)
# angular coordinate system has center (U,V,W)
#
# direction cosines wrt angular coordinate center are (alpha,beta,gamma)
-# i.e. a point has coordinates (U+alpha*r, V+beta*r, W+gamma*r)
+# but Maple predefines gamma = Euler's constant, so we use (xcos,ycos,zcos)
+# instead, i.e. a point has coordinates (U+xcos*r, V+ycos*r, W+zcos*r)
#
# then the equation of the ellipsoid is
-# (U+alpha*r - A)^2 (V+beta*r - B)^2 (W+gamma*r - C)^2
+# (U+xcos*r - A)^2 (V+ycos*r - B)^2 (W+zcos*r - C)^2
# ----------------- + ---------------- + ----------------- = 1
# a^2 b^2 c^2
#
@@ -23,12 +24,12 @@
# BV = B - V
# CW = C - W
#
-> eqn := (alpha*r - AU)^2/a^2 + (beta*r - BV)^2/b^2 + (gamma*r - CW)^2/c^2 = 1;
- 2 2 2
- (alpha r - AU) (beta r - BV) (gamma r - CW)
- eqn := --------------- + -------------- + --------------- = 1
- 2 2 2
- a b c
+> eqn := (xcos*r - AU)^2/a^2 + (ycos*r - BV)^2/b^2 + (zcos*r - CW)^2/c^2 = 1;
+ 2 2 2
+ (xcos r - AU) (ycos r - BV) (zcos r - CW)
+ eqn := -------------- + -------------- + -------------- = 1
+ 2 2 2
+ a b c
>
> read "../maple/util.mm";
@@ -384,136 +385,132 @@ end proc
>
> [solve(eqn, r)];
-bytes used=1000108, alloc=917336, time=0.09
-bytes used=2000516, alloc=1376004, time=0.15
- 2 2 2 2 2 2
-[1/2 (2 a b gamma CW + 2 a c beta BV + 2 b c alpha AU + 2 (
+bytes used=1000572, alloc=851812, time=0.07
+ 2 2 2 2 2 2
+[1/2 (2 a b zcos CW + 2 a c ycos BV + 2 b c xcos AU + 2 (
- 4 2 2 2 4 2
- 2 a b gamma CW c beta BV + 2 a b gamma CW c alpha AU
+ 4 2 2 2 4 2
+ 2 a b zcos CW c ycos BV + 2 a b zcos CW c xcos AU
- 2 4 2 4 2 2 2 2 2 4 2 2 2
- + 2 a c beta BV b alpha AU - b c alpha a CW - b c alpha a BV
+ 2 4 2 4 2 2 2 2 2 4 2 2 2
+ + 2 a c ycos BV b xcos AU - b c xcos a CW - b c xcos a BV
- 4 4 2 2 4 2 2 2 2 2 4 2 2 2
- + b c alpha a - a c beta b CW - a c beta b AU
+ 4 4 2 2 4 2 2 2 2 2 4 2 2 2
+ + b c xcos a - a c ycos b CW - a c ycos b AU
- 4 4 2 2 4 2 2 2 2 2 4 2 2 2
- + a c beta b - a b gamma c BV - a b gamma c AU
+ 4 4 2 2 4 2 2 2 2 2 4 2 2 2
+ + a c ycos b - a b zcos c BV - a b zcos c AU
- 4 4 2 2 1/2 / 2 2 2 2 2 2 2 2 2
- + a b gamma c ) ) / (a b gamma + a c beta + b c alpha ),
- /
+ 4 4 2 2 1/2 / 2 2 2 2 2 2 2 2 2
+ + a b zcos c ) ) / (b c xcos + a c ycos + a b zcos ), 1/2 (
+ /
- 2 2 2 2 2 2
- 1/2 (2 a b gamma CW + 2 a c beta BV + 2 b c alpha AU - 2 (
+ 2 2 2 2 2 2
+ 2 a b zcos CW + 2 a c ycos BV + 2 b c xcos AU - 2 (
- 4 2 2 2 4 2
- 2 a b gamma CW c beta BV + 2 a b gamma CW c alpha AU
+ 4 2 2 2 4 2
+ 2 a b zcos CW c ycos BV + 2 a b zcos CW c xcos AU
- 2 4 2 4 2 2 2 2 2 4 2 2 2
- + 2 a c beta BV b alpha AU - b c alpha a CW - b c alpha a BV
+ 2 4 2 4 2 2 2 2 2 4 2 2 2
+ + 2 a c ycos BV b xcos AU - b c xcos a CW - b c xcos a BV
- 4 4 2 2 4 2 2 2 2 2 4 2 2 2
- + b c alpha a - a c beta b CW - a c beta b AU
+ 4 4 2 2 4 2 2 2 2 2 4 2 2 2
+ + b c xcos a - a c ycos b CW - a c ycos b AU
- 4 4 2 2 4 2 2 2 2 2 4 2 2 2
- + a c beta b - a b gamma c BV - a b gamma c AU
+ 4 4 2 2 4 2 2 2 2 2 4 2 2 2
+ + a c ycos b - a b zcos c BV - a b zcos c AU
- 4 4 2 2 1/2 / 2 2 2 2 2 2 2 2 2
- + a b gamma c ) ) / (a b gamma + a c beta + b c alpha )]
- /
+ 4 4 2 2 1/2 / 2 2 2 2 2 2 2 2 2
+ + a b zcos c ) ) / (b c xcos + a c ycos + a b zcos )]
+ /
> map(simplify, %);
-bytes used=3007136, alloc=1769148, time=0.22
-bytes used=4007324, alloc=1834672, time=0.29
- 2 2 2 2 2 2 2 2 2
-[(a b gamma CW + a c beta BV + b c alpha AU + (-a b c (
+bytes used=2000744, alloc=1376004, time=0.12
+bytes used=3001976, alloc=1769148, time=0.19
+ 2 2 2 2 2 2 2 2 2
+[(a b zcos CW + a c ycos BV + b c xcos AU + (-a b c (
- 2 2 2
- -2 a gamma CW beta BV - 2 b gamma CW alpha AU - 2 c beta BV alpha AU
+ 2 2 2
+ -2 a zcos CW ycos BV - 2 b zcos CW xcos AU - 2 c ycos BV xcos AU
- 2 2 2 2 2 2 2 2 2 2 2 2
- + b alpha CW + c alpha BV - b c alpha + a beta CW
+ 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
+ + b xcos CW + c xcos BV - b c xcos + a ycos CW + c ycos AU
- 2 2 2 2 2 2 2 2 2 2 2 2
- + c beta AU - a c beta + a gamma BV + b gamma AU
+ 2 2 2 2 2 2 2 2 2 2 2 2 1/2 /
+ - a c ycos + a zcos BV + b zcos AU - a b zcos )) ) / (
+ /
- 2 2 2 1/2 / 2 2 2 2 2 2 2 2 2
- - a b gamma )) ) / (a b gamma + a c beta + b c alpha ), (
- /
+ 2 2 2 2 2 2 2 2 2 2 2 2 2
+ b c xcos + a c ycos + a b zcos ), (a b zcos CW + a c ycos BV
- 2 2 2 2 2 2 2 2 2
- a b gamma CW + a c beta BV + b c alpha AU - (-a b c (
+ 2 2 2 2 2 2 2
+ + b c xcos AU - (-a b c (-2 a zcos CW ycos BV - 2 b zcos CW xcos AU
- 2 2 2
- -2 a gamma CW beta BV - 2 b gamma CW alpha AU - 2 c beta BV alpha AU
+ 2 2 2 2 2 2 2 2 2 2
+ - 2 c ycos BV xcos AU + b xcos CW + c xcos BV - b c xcos
- 2 2 2 2 2 2 2 2 2 2 2 2
- + b alpha CW + c alpha BV - b c alpha + a beta CW
+ 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
+ + a ycos CW + c ycos AU - a c ycos + a zcos BV + b zcos AU
- 2 2 2 2 2 2 2 2 2 2 2 2
- + c beta AU - a c beta + a gamma BV + b gamma AU
-
- 2 2 2 1/2 / 2 2 2 2 2 2 2 2 2
- - a b gamma )) ) / (a b gamma + a c beta + b c alpha )]
- /
+ 2 2 2 1/2 / 2 2 2 2 2 2 2 2 2
+ - a b zcos )) ) / (b c xcos + a c ycos + a b zcos )]
+ /
> [r_plus = %[1], r_minus = %[2]];
- 2 2 2 2 2 2 2 2 2
-[r_plus = (a b gamma CW + a c beta BV + b c alpha AU + (-a b c (
+ 2 2 2 2 2 2 2 2 2
+[r_plus = (a b zcos CW + a c ycos BV + b c xcos AU + (-a b c (
- 2 2 2
- -2 a gamma CW beta BV - 2 b gamma CW alpha AU - 2 c beta BV alpha AU
+ 2 2 2
+ -2 a zcos CW ycos BV - 2 b zcos CW xcos AU - 2 c ycos BV xcos AU
- 2 2 2 2 2 2 2 2 2 2 2 2
- + b alpha CW + c alpha BV - b c alpha + a beta CW
+ 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
+ + b xcos CW + c xcos BV - b c xcos + a ycos CW + c ycos AU
- 2 2 2 2 2 2 2 2 2 2 2 2
- + c beta AU - a c beta + a gamma BV + b gamma AU
+ 2 2 2 2 2 2 2 2 2 2 2 2 1/2 /
+ - a c ycos + a zcos BV + b zcos AU - a b zcos )) ) / (
+ /
- 2 2 2 1/2 / 2 2 2 2 2 2 2 2 2
- - a b gamma )) ) / (a b gamma + a c beta + b c alpha ),
- /
+ 2 2 2 2 2 2 2 2 2 2 2
+ b c xcos + a c ycos + a b zcos ), r_minus = (a b zcos CW
- 2 2 2 2 2 2 2 2 2
- r_minus = (a b gamma CW + a c beta BV + b c alpha AU - (-a b c (
+ 2 2 2 2 2 2 2 2
+ + a c ycos BV + b c xcos AU - (-a b c (-2 a zcos CW ycos BV
- 2 2 2
- -2 a gamma CW beta BV - 2 b gamma CW alpha AU - 2 c beta BV alpha AU
+ 2 2 2 2 2
+ - 2 b zcos CW xcos AU - 2 c ycos BV xcos AU + b xcos CW
- 2 2 2 2 2 2 2 2 2 2 2 2
- + b alpha CW + c alpha BV - b c alpha + a beta CW
+ 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
+ + c xcos BV - b c xcos + a ycos CW + c ycos AU - a c ycos
- 2 2 2 2 2 2 2 2 2 2 2 2
- + c beta AU - a c beta + a gamma BV + b gamma AU
+ 2 2 2 2 2 2 2 2 2 1/2 /
+ + a zcos BV + b zcos AU - a b zcos )) ) / (
+ /
- 2 2 2 1/2 / 2 2 2 2 2 2 2 2 2
- - a b gamma )) ) / (a b gamma + a c beta + b c alpha )]
- /
+ 2 2 2 2 2 2 2 2 2
+ b c xcos + a c ycos + a b zcos )]
> solnlist := [codegen[optimize](%)];
- 2 2 2
-solnlist := [t1 = a , t2 = b , t3 = t1 t2, t5 = t3 gamma CW, t6 = c ,
+ 2 2 2
+solnlist := [t1 = a , t2 = b , t3 = t1 t2, t5 = t3 zcos CW, t6 = c , t7 = t1 t6,
- 2
- t7 = t1 t6, t9 = t7 beta BV, t10 = t2 t6, t12 = t10 alpha AU, t28 = alpha ,
+ 2 2
+ t9 = t7 ycos BV, t10 = t2 t6, t12 = t10 xcos AU, t28 = xcos , t30 = CW ,
- 2 2 2 2
- t30 = CW , t33 = BV , t35 = t10 t28, t36 = beta , t40 = AU , t42 = t7 t36,
+ 2 2 2 2
+ t33 = BV , t35 = t10 t28, t36 = ycos , t40 = AU , t42 = t7 t36, t43 = zcos ,
- 2
- t43 = gamma , t48 = t3 t43, t49 = -2 t1 gamma CW beta BV
+ t48 = t3 t43, t49 = -2 t1 zcos CW ycos BV - 2 t2 zcos CW xcos AU
- - 2 t2 gamma CW alpha AU - 2 t6 beta BV alpha AU + t2 t28 t30 + t6 t28 t33
+ - 2 t6 ycos BV xcos AU + t2 t28 t30 + t6 t28 t33 - t35 + t1 t36 t30
- - t35 + t1 t36 t30 + t6 t36 t40 - t42 + t1 t43 t33 + t2 t43 t40 - t48,
+ 1/2
+ + t6 t36 t40 - t42 + t1 t43 t33 + t2 t43 t40 - t48, t52 = (-t3 t6 t49) ,
- 1/2 1
- t52 = (-t3 t6 t49) , t55 = ---------------,
- t48 + t42 + t35
+ 1
+ t55 = ---------------, r_plus = (t5 + t9 + t12 + t52) t55,
+ t35 + t42 + t48
- r_plus = (t5 + t9 + t12 + t52) t55, r_minus = (t5 + t9 + t12 - t52) t55]
+ r_minus = (t5 + t9 + t12 - t52) t55]
>
> ftruncate("ellipsoid.c");
@@ -522,5 +519,6 @@ solnlist := [t1 = a , t2 = b , t3 = t1 t2, t5 = t3 gamma CW, t6 = c ,
> print_name_list_dcl(temps_in_eqnlist(solnlist, [r_plus,r_minus]),
> "fp", "ellipsoid.c");
> codegen[C](solnlist, filename="ellipsoid.c");
+bytes used=4002168, alloc=1769148, time=0.26
> quit
-bytes used=4749864, alloc=1834672, time=0.37
+bytes used=4027884, alloc=1769148, time=0.27
generated by cgit v1.2.3 (git 2.39.1) at 2024-07-04 16:04:50 +0000