restart;
with( plots ):
Miscellaneous Definitions
ELL11 := (x-1)^2 + y^2 = 1:
ELL12 := (x-1)^2 + y^2/4 = 1:
ELL1a := (x-1)^2 + y^2/a^2 = 1:
ParabolaH := a*y^2=x:
ParabolaV := a*x^2=y:
Line := y = a*x:
PtInFirstQuad := proc (pt, P)
local q;
q := evalf(eval(pt, P));
return evalb((eval(x, q))::nonnegative and (eval(y, q))::nonnegative)
end proc:
PCr := r -> implicitplot( Cr(r), x=-2..2, y=-2..2, grid=[25,25] ):
Ppt := pts -> plot( pts, style= [point, point, point],
symbol=[box, diamond, diamond],
color= [blue, green, pink] ):
Pline := pts -> plot( pts, style=line, color=magenta ):
Pframe := rr -> display( [PFixedCurve,
PCr(rr),
Ppt( eval([P,Q,R], r=rr)),
Pline(eval([P,R], r=rr))], scaling=constrained):
ApplyLHRule := proc( LIMIT )
local f,x,u,L, dir, f1, f2;
if not member( op(0,LIMIT), {limit, Limit} ) then break end if;
f := normal(op(1,LIMIT));
x,u := [lhs,rhs](op(2,LIMIT))[];
if nops(LIMIT)>2 then dir := op(3,LIMIT) else dir := NULL end if;
f1 := numer(f);
f2 := denom(f);
if limit(f1,x=u,dir)=0 and limit(f2,x=u,dir)=0 then
return op(0,LIMIT)( simplify(diff(f1,x)/diff(f2,x)), x=u, dir );
else
WARNING("this limit is not an indeterminate form of type 0/0")
end if;
end proc:
MultNumerDenom := proc( f, term )
local expr;
expr := normal( f );
return numer(expr)*term / expand( denom(expr)*term );
end proc:
Case 1: An Explicit Example (Circle)
Symbolic Solution
SOLNs := [allvalues( solve( {FixedCurve,Cr(r)}, {x,y} ) )];
# map( PtInFirstQuad, SOLNs, [r=1,a=2] );
SOLN := select( s -> PtInFirstQuad(s,[r=1,a=2]), SOLNs )[];
P := [0,r];
Q := eval( [x,y], SOLN );
rayPQ := (y-P[2])/(x-P[1]) = (Q[2]-P[2])/(Q[1]-P[1]);
R := [ solve( eval(rayPQ,y=0), x ), 0 ];
Case 2: An Explicit Example (Specific Ellipse)
Symbolic Solution
SOLNs := [allvalues( solve( {FixedCurve,Cr(r)}, {x,y} ) )];
# map( PtInFirstQuad, SOLNs, [r=1,a=2] );
SOLN := select( s -> PtInFirstQuad(s,[r=1,a=2]), SOLNs )[];
P := [0,r];
Q := eval( [x,y], SOLN );
rayPQ := (y-P[2])/(x-P[1]) = (Q[2]-P[2])/(Q[1]-P[1]);
R := [ solve( eval(rayPQ,y=0), x ), 0 ];
Symbolic Evaluation of Limit
Rlimit := Limit( R[1], r=0, right ):
Rlimit = value( Rlimit );
Limit( R[1], r=0, right );
ApplyLHRule( % );
series( numer(R[1]), r=0 )/series( denom(R[1]), r=0 );
limit( %, r=0, right );
R[1];
MultNumerDenom( R[1], denom(R[1])-6*r );
MultNumerDenom( %, denom(%)+64-6*r^2 );
Case 3: The General Ellipse
Symbolic Solution
SOLNs := [allvalues( solve( {FixedCurve,Cr(r)}, {x,y} ) )];
# map( PtInFirstQuad, SOLNs, [r=1,a=2] );
SOLN := select( s -> PtInFirstQuad(s,[r=1,a=2]), SOLNs )[];
P := [0,r];
Q := eval( [x,y], SOLN );
rayPQ := (y-P[2])/(x-P[1]) = (Q[2]-P[2])/(Q[1]-P[1]);
R := [ solve( eval(rayPQ,y=0), x ), 0 ];
Symbolic Evaluation of Limit
Rlimit := Limit( R[1], r=0, right ):
Rlimit = value( Rlimit );# assuming a>1;
Rlimit = value( Rlimit ) assuming a>1;
Rlimit;
ApplyLHRule( % ) assuming a>1;
ApplyLHRule( % ) assuming a>1;
ApplyLHRule( % ) assuming a>1;
series( numer(R[1]), r=0 ) / series( denom(R[1]), r=0 );
collect(convert(series( convert(series( numer(R[1]), r=0, 4 ),polynom),a=1),polynom),r);
collect(convert(series( convert(series( denom(R[1]), r=0 ),polynom),a=1),polynom),r) assuming a::real;
%%/%;
simplify(%);
plot(%,a=0..2);
limit( R[1], r=0, right ) assuming a>1;
Case 4: The General Parabola
Symbolic Solution
SOLNs := [allvalues( solve( {FixedCurve,Cr(r)}, {x,y} ) )];
# map( PtInFirstQuad, SOLNs, [r=1,a=2] );
SOLN := select( s -> PtInFirstQuad(s,[r=1,a=2]), SOLNs )[];
P := [0,r];
Q := eval( [x,y], SOLN );
rayPQ := (y-P[2])/(x-P[1]) = (Q[2]-P[2])/(Q[1]-P[1]);
R := [ solve( eval(rayPQ,y=0), x ), 0 ];
Case 5: The Other General Parabola
Symbolic Solution
SOLNs := [allvalues( solve( {FixedCurve,Cr(r)}, {x,y} ) )];
# map( PtInFirstQuad, SOLNs, [r=1,a=2] );
SOLN := select( s -> PtInFirstQuad(s,[r=1,a=2]), SOLNs )[];
P := [0,r];
Q := eval( [x,y], SOLN );
rayPQ := (y-P[2])/(x-P[1]) = (Q[2]-P[2])/(Q[1]-P[1]);
R := [ solve( eval(rayPQ,y=0), x ), 0 ];
Case 6: The Straight Line
Symbolic Solution
SOLNs := [allvalues( solve( {FixedCurve,Cr(r)}, {x,y} ) )];
# map( PtInFirstQuad, SOLNs, [r=1,a=2] );
SOLN := select( s -> PtInFirstQuad(s,[r=1,a=2]), SOLNs )[];
P := [0,r];
Q := eval( [x,y], SOLN );
rayPQ := (y-P[2])/(x-P[1]) = (Q[2]-P[2])/(Q[1]-P[1]);
R := [ solve( eval(rayPQ,y=0), x ), 0 ];