## 1152 Reputation

12 years, 126 days

## A proof...

I got curious. I have converted away from the geometry package and used a rational parameterisation for cos and sin

[a cos(t), b sin(t)]  to  [a* (1t^2)/(1+t^2), b*2*t/(1+t^2)]   . In this case t is not the angle.

This proves the angle is constant upto sign changes for signum . But that just changes the reported angle from sfor example 50deg to 130deg.

I would like to know the name of this theorem.

 > restart; with(plots):with(LinearAlgebra):
 > x0 := 100; y0 := 40; a := 7; b := 5; c := sqrt(a^2 - b^2); e1:= x^2/a^2 + y^2/b^2 - 1; F1:=[ -c, 0]; F2:=[ c, 0]; eq := simplify((a^2 - x0^2)*(y - y0)^2 + (b^2 - y0^2)*(x - x0)^2 + 2*x0*y0*(x - x0)*(y - y0)) = 0; sol := solve({eq}, {y}); tang1:=(lhs-rhs)~(op(sol[1])); tang2:= (lhs-rhs)~(op(sol[2])); sol2 := op(solve({op(sol[1]), x^2/a^2 + y^2/b^2 - 1 }, {x, y})); xM2 := rhs(sol2[1]); yM2 := rhs(sol2[2]); A:=[ xM2, yM2]; sol3 := op(solve({op(sol[2]), x^2/a^2 + y^2/b^2 - 1 }, {x, y})); xM3 := rhs(sol3[1]); yM3 := rhs(sol3[2]); B:=[ xM3, yM3]; slpvecAB:=<(B[2]-A[2]),B[1]-A[1]>;  #slpvec Pol:=-slpvecAB[1]*(x-A[1])+slpvecAB[2]*(y-A[2]); simplify(Pol); isolate(%, y); TANG := proc(t) local xM, yM; xM :=a*(1-t^2)/(1+t^2);# a*cos(t);                yM :=b*2*t/(1+t^2);# b*sin(t);                return(expand(1/49*x*xM + 1/25*y*yM - 1 ));         end proc; #t := -0.25; TT:= TANG(t); #Equation(TT); M:=[a*(1-t^2)/(1+t^2),b*2*t/(1+t^2)];#[ a*cos(t), b*sin(t)];
 (1)
 >
 > P:=rhs~(solve( [tang1, TT],[x,y])[]);
 (2)
 > Q:=rhs~(solve( [tang2, TT],[x,y])[]);
 (3)
 > slpvecPF2:=simplify(<(P[2]-F2[2]),P[1]-F2[1]>); PF2:=simplify(-slpvecPF2[1]*(x-P[1])+slpvecPF2[2]*(y-P[2]));
 >
 (4)
 > slpvecQF2:=simplify(<(Q[2]-F2[2]),Q[1]-F2[1]>); QF2:=simplify(-slpvecQF2[1]*(x-Q[1])+slpvecQF2[2]*(y-Q[2]));
 (5)
 >
 >
 > alpha :=(simplify(VectorAngle(slpvecQF2,slpvecPF2),assume=real))
 (6)
 > Setof:=frontend(indets,[alpha])
 (7)
 > (sin(Setof[2]))^2
 (8)
 > simpascn:=arcsin(sqrt(factor( (8) ))); # the t's eliminate
 (9)
 > alpha1:=(subs(Setof[2]=simpascn,alpha))*180/Pi
 (10)
 > plot(alpha1,t=-10..10,view=[-10..10,0..180])
 >
 >
 >
 >

## One possible way...

Hope this helps. There are probably more efficient ways.

 (1)

 (2)

## How about using the Backup...

Click in Restore Backup and list the files by date saved,

## Use fsolve...

The set of equations; as far as I can see, cannot be solved explicitly do to the trig. functions.
I rearranged your equations so you only need to solve eq5 and eq6 for y and z.

I randomly made up values for  t, T, a,c,e,f just to show the use of fsolve.

Hope this helps you.

Edit I re-attached the worksheet because I had a mistake in it.

 > restart
 > t:=Pi/2; T:=t;
 (1)
 > a:=3;c:=2;e:=4;f:=2
 (2)
 > #z:=1
 > w := sqrt(a^2 - (a - y)^2); v := cot(t + arcsin(w/a)) ;  u := (sqrt(c^2 - (c - z)^2) - v*x - w)/x^2;  x := y + z + e; eq5 := v + cot(T + arcsin((u*x^2 + v*x + w)/c)) + 2*u*x = 0; eq6 :=f = Pi/30*(6*u^2*x^5 + 15*u*v*x^4 + (20*u*w + 10*v^2)*x^3 + 30*v*w*x^2 + 30*x*w^2) - Pi/3*(-z^3 + 3*z^2*c - y^3 + 3*y^2*a);
 (3)
 > fsolve({ eq5,eq6},complex);
 (4)
 >

## One way...

```f:=x-> (1-k*x)/(1+x^2)
f := proc (x) options operator, arrow; (1-k*x)/(1+x^2) end proc

slpf:=diff(f(x),x)
k      2 (-k x + 1) x
slp := - ------ - --------------
2                 2
x  + 1     / 2    \
\x  + 1/

m:=eval(slpf,x=3);# slope of tangent line at x=3

2      3
m := -- k - --
25     50

f(3)
3      1
- -- k + --
10     10

line:=unapply(m*x+c,x)
line := proc (x) options operator, arrow; ((2/25)*k-3/50)*x+c

end proc

line(3)=f(3);# x=3 is tangent point
6      9          3      1
-- k - -- + c = - -- k + --
25     50         10     10

c:=solve(line(3)=f(3),c)
27     7
c := - -- k + --
50     25

line(x)
/2      3 \     27     7
|-- k - --| x - -- k + --
\25     50/     50     25

```

## Check under Options Display...

Check under Tools Options Display Typesetting level to see if it is set to extended.

## It looks like I have it working correctl...

@dharr  I appear to have the setup working in principal. I created seperate worksheets. One for the Main Package and one for the Sub Package. Both attached.  I have hilighted the differences in both in red. Then I copied both help files across to the correct location. They now list correctly in Help.

## One way...

Interesting problem.  I found this method in the RealDomain help page.

```restart;

n:=3;m:=2;
eqx:=x^(n/m)=a;
use RealDomain in (maple_sol:=solve(eqx,[x]))[] end use;  #also tried solve()
F:=map(X->eval(eqx,X),maple_sol);
map(X->evalb(X),F);
n := 3

m := 2

(3/2)
eqx := x      = a

[     (2/3)]
[x = a     ]

F := [a = a]

[true]```

## A way....

This is rather similar to the 1st question I asked on MaplePrimes. I just happened to look at it today

 > restart

 (1)
 >

## Would this do?...

This is a way using seq

```M := Matrix(4, 4, [[m[1, 1], m[1, 2], m[1, 3], m[1, 4]],
[m[2, 1], m[2, 2], m[2, 3], m[2, 4]],
[m[3, 1], m[3, 2], m[3, 3], m[3, 4]],
[m[4, 1], m[4, 2], m[4, 3], m[4, 4]]]);

L:=[seq(M[1 .. 4, i], i = 1 .. 4)];

whattype(L[3]);

```

I am sure there are more direct ways.1st, I changed you numbers just to see if this works more generally. Works with you numbers too

```a:=evalf(1119.8*21/101)
a := 232.8297030

b:=trunc(a)
b := 232

c:=evalf(a-b)
c := 0.8297030

d:=evalf(c,2)
d := 0.83

b+d
232.83

a:=evalf(3*21/100,3);
a := 0.630

b:=trunc(a)
b := 0

c:=evalf(a-b)
c := 0.630

d:=evalf(c,2)
d := 0.63

b+d
0.63
```

## Is this what you expect?...

Your code has so many errors. Are you typing it in a text editor? Hard to see how you could enter that much code in maple and not notice or find errors along the way.

This runs but I cannot verify you logic.

I used the VS code editor to find the errors. and maplemint That might help you.

```Cong:=proc(n)
local  a,b,An,Bn,Cn,Dn:
if n mod 2 = 1  then
An:=0:
Bn:=0:
for a from (round(-sqrt(n/(2)))) to round(sqrt(n/(2)) ) do
for b  from (round(-sqrt(n)) )to round(sqrt(n) )do
if is( (sqrt(n-2*a^(2)-b^(2)) )/(32),integer  ) then
An:=An+1 ;
elif is((sqrt(n-2*a^(2)-b^(2)) )/(8) ,integer )  then
Bn:=Bn+1 ;
end if:
end do:
end do:
end if;
if 2*An=Bn  then
print(True);
else
print(False);
end if;
if n mod 2 = 0 then
Cn:=0:
Dn:=0:
end if;
for a  from (round(-sqrt(n/(8)))) to round(sqrt(n/(8)) ) do
for b  from (round(-sqrt(n/(2)))) to round(sqrt(n/(2))) do
if is((sqrt(n/(2)-4*a^(2)-b^(2)) )/(32), integer   ) then
Cn:=Cn+1 ;
elif is((sqrt(n/(2)-4*a^(2)-b^(2)) )/(8) ,integer) then
Dn:=Dn+1 ;
end if:
end do:
end do:
if 2*Cn=Dn then
print("True");
else
print("False");
end if;
end proc;

maplemint(Cong);
Cong(108);```

## Solve(diff.... for z...

If I am intrepetering things correctly.  Your

`diff((z - zA)*(z - zB)*(z - zC), z)*I);   # I removed semicolin inside brackets`

gives a 2nd order polynomial (quadratic)   .say it is Q(Z)    Why have you multiplied the whole diff eqn. by I

solve(Q(z) , z ,explicit)   #gives two roots which are supposed to be the focii of the ellipse.

Look on Marden's theorem - Wikipedia

## Factrix a really useful proc...

This is a really useful procedure I found here on Primes  years ago. Called factrix. though in this case it pulls out (-1/2) instead of (1/2)
Edit:- the link to the reference to factrix  How do I extract common factors in a matrix? - MaplePrimes

[moderator: Copyrighted code removed. Source is factrix by Robert Israel, part of his Maple Advisor Database.]

## Use expand...

You have to expand the equation.

`if  expand((x+1)^2=x^2+2*x+1) then print(1) else print(0) fi`

or

`if  expand((x+1)^2)=x^2+2*x+1 then print(1) else print(0) fi`

 1 2 3 4 5 6 7 Page 1 of 7
﻿