Kitonum

21083 Reputation

26 Badges

16 years, 197 days

MaplePrimes Activity


These are Posts that have been published by Kitonum

Yesterday, I accidentally discovered a nasty bug in a fairly simple example:

restart;
Expr:=a*sin(x)+b*cos(x);
maximize(Expr, x=0..2*Pi);
minimize(Expr, x=0..2*Pi);
                                    

I am sure the correct answers are  sqrt(a^2+b^2)  and  -sqrt(a^2+b^2)  for any real values  a  and  b .  It is easy to prove in many ways. The simplest method does not require any calculations and can be done in the mind. We will consider  Expr  as the scalar product (or the dot product) of two vectors  <a, b>  and  <sin(x), cos(x)>, one of which is a unit vector. Then it is obvious that the maximum of this scalar product is reached if the vectors are codirectional and equals to the length of the first vector, that is, sqrt(a^2+b^2).

Bugs in these commands were noted by users and earlier (see search by keywords bug, maximize, minimize) but unfortunately are still not fixed. 

In this post an interesting geometric problem is solved: for an arbitrary convex polygon, find a straight line that divides both the area and the perimeter in half. The following results on this problem are well known:
1. For any convex polygon there is such a straight line.
2. For any convex polygon into which a circle can be inscribed, in particular for any triangle, the desired line must pass through the center of the inscribed circle.
3. For a triangle, the number of solutions can be 1, 2, or 3.
4. If the polygon has symmetry with respect to a point, then any straight line passing through this point is a solution.

The procedure called  InHalf  (the code below) symbolically solves this problem. The formal parameter of the procedure is the list of coordinates of the vertices of a convex polygon (vertices must be passed opposite or clockwise). The procedure returns all solutions in the form of a list of pairs of points, lying on the perimeter of the polygon, that are the ends of segments that implement the desired dividing.


 

restart;

InHalf:=proc(V::listlist)
local L, n, a, b, M, N, i, j, P, Q, L1, L2, Area, Area1, Area2, Perimeter, Perimeter1, Perimeter2, sol, m, k, Sol;
uses LinearAlgebra, ListTools;
L:=map(convert,[V[],V[1]],rational); n:=nops(L)-1;
a:=<(V[2]-V[1])[1],(V[2]-V[1])[2],0>; b:=<(V[n]-V[1])[1],(V[n]-V[1])[2],0>;
if is(CrossProduct(a,b)[3]<0) then L:=Reverse(L) fi;
M:=[seq([L[i],L[i+1]], i=1..n)]:
N:=0;
for i from 1 to n-1 do
for j from i+1 to n do
P:=map(t->t*(1-s),M[i,1])+map(t->t*s,M[i,2]); Q:=map(s->s*(1-t),M[j,1])+map(s->s*t,M[j,2]);
L1:=[P,L[i+1..j][],Q,P];
L2:=[Q,L[j+1..-1][],L[1..i][],P,Q];
Area:=L->(1/2)*add(L[k, 1]*L[k+1, 2]-L[k, 2]*L[k+1, 1], k = 1 .. nops(L)-1);
Area1:=Area(L1);
Area2:=Area(L2);
Perimeter:=L->add(sqrt((L[k,1]-L[k+1,1])^2+(L[k,2]-L[k+1,2])^2), k=1..nops(L)-2);
Perimeter1:=Perimeter(L1);
Perimeter2:=Perimeter(L2);
sol:=[solve({Area1=Area2,Perimeter1=Perimeter2,s>=0,s<1,t>=0,t<1}, {s,t}, explicit)] assuming real;
if sol<>[] then m:=nops(sol);
for k from 1 to m do
N:=N+1; if nops(sol[k])=2 then Sol[N]:=simplify(eval([P,Q],sol[k])) else Sol[N]:=simplify(eval([P,Q],s=t)) fi;
od; fi;
od; od;
Sol:=convert(Sol, list);
`if`(indets(Sol)={},Sol,op([Sol,t>=0 and t<1]));
end proc:  


Examples of use

# For the Pythagorean triangle with sides 3, 4, 5, we have a unique solution

L:=[[4,3],[4,0],[0,0]]:
P:=InHalf(L);
plots:-display(plot([L[],L[1]], color=green, thickness=3), plot(P,  color=red), scaling=constrained);

[[[8/5-(2/5)*6^(1/2), 6/5-(3/10)*6^(1/2)], [4, (1/2)*6^(1/2)]]]

 

 

# For an isosceles right triangle, there are 3 solutions. We see that all the cuts pass through the center of the inscribed circle

L:=[[0,0],[4,0],[4,4]]:
InHalf(L);
P:=InHalf(L);
r:=(4+4-4*sqrt(2))/2: a:=4-r: b:=r:
plots:-display(plot([L[],L[1]], color=green, thickness=3), plot(P,  color=red), plot([r*cos(t)+a,r*sin(t)+b, t=0..2*Pi], color=blue), scaling=constrained);

[[[2*2^(1/2), 0], [2*2^(1/2), 2*2^(1/2)]], [[4, 0], [2, 2]], [[4, -2*2^(1/2)+4], [-2*2^(1/2)+4, -2*2^(1/2)+4]]]

 

[[[2*2^(1/2), 0], [2*2^(1/2), 2*2^(1/2)]], [[4, 0], [2, 2]], [[4, -2*2^(1/2)+4], [-2*2^(1/2)+4, -2*2^(1/2)+4]]]

 

 

# There are 3 solutions for the quadrilateral below

L:=[[0,0],[4.5,0],[4,3],[0,2]]:
P:=InHalf(L);
plots:-display(plot([L[],L[1]], color=green, thickness=3), plot(P,  color=red), scaling=constrained);

[[[(1/44844)*6^(1/2)*(17^(1/2)-13/2)*37^(3/4)*(2*17^(1/2)+13)^(1/2)*((1836*37^(1/2)-6956)*17^(1/2)+7995*37^(1/2)-56425)^(1/2)-(1/4)*17^(1/2)-(1/8)*37^(1/2)+23/8, 0], [(6^(1/2)*37^(1/4)*((1836*37^(1/2)-6956)*17^(1/2)+7995*37^(1/2)-56425)^(1/2)*(2*17^(1/2)+13)^(1/2)+(-156*37^(1/2)+7770)*17^(1/2)-711*37^(1/2)+50505)/(1776*17^(1/2)+11544), (-6^(1/2)*37^(1/4)*((1836*37^(1/2)-6956)*17^(1/2)+7995*37^(1/2)-56425)^(1/2)*(2*17^(1/2)+13)^(1/2)+(156*37^(1/2)+222)*17^(1/2)+711*37^(1/2)+1443)/(296*17^(1/2)+1924)]], [[(1/90576)*17^(3/4)*(37^(1/2)+37)^(1/2)*(37^(1/2)-37)*((1461*37^(1/2)+29415)*17^(1/2)-986*37^(1/2)-149702)^(1/2)+(5/4)*17^(1/2)+(1/8)*37^(1/2)-27/8, 0], [(17^(1/4)*((1461*37^(1/2)+29415)*17^(1/2)-986*37^(1/2)-149702)^(1/2)*(37^(1/2)+37)^(1/2)+(37*17^(1/2)-51)*37^(1/2)+703*17^(1/2)-1887)/(17*37^(1/2)+629), (17^(1/4)*((1461*37^(1/2)+29415)*17^(1/2)-986*37^(1/2)-149702)^(1/2)*(37^(1/2)+37)^(1/2)+(37*17^(1/2)+85)*37^(1/2)+703*17^(1/2)+3145)/(68*37^(1/2)+2516)]], [[-(1/90576)*17^(3/4)*(37^(1/2)+37)^(1/2)*(37^(1/2)-37)*((1461*37^(1/2)+29415)*17^(1/2)-986*37^(1/2)-149702)^(1/2)+(5/4)*17^(1/2)+(1/8)*37^(1/2)-27/8, 0], [(-17^(1/4)*((1461*37^(1/2)+29415)*17^(1/2)-986*37^(1/2)-149702)^(1/2)*(37^(1/2)+37)^(1/2)+(37*17^(1/2)-51)*37^(1/2)+703*17^(1/2)-1887)/(17*37^(1/2)+629), (-17^(1/4)*((1461*37^(1/2)+29415)*17^(1/2)-986*37^(1/2)-149702)^(1/2)*(37^(1/2)+37)^(1/2)+(37*17^(1/2)+85)*37^(1/2)+703*17^(1/2)+3145)/(68*37^(1/2)+2516)]]]

 

 

# There are infinitely many solutions for a polygon with a center of symmetry. Any cut through the center solves the problem. The picture shows 2 solutions.

L:=[[1,0],[1+2*sqrt(3),2],[2*sqrt(3),sqrt(3)+2],[0,sqrt(3)]]:
P:=InHalf(L);
plots:-display(plot([L[],L[1]], color=green, thickness=3), plot(eval(P[1],t=1/3),  color=red), scaling=constrained);

[[[2*3^(1/2)*t+1, 2*t], [-2*3^(1/2)*(t-1), 3^(1/2)-2*t+2]], [[2*3^(1/2)-t+1, 3^(1/2)*t+2], [t, -3^(1/2)*(t-1)]]], 0 <= t and t < 1

 

 

 


 

Download In_Half.mw

Recently I looked through an interesting book D. Wells "The Penquin book of Curious and Interesting Geometry" and came across this result, which I did not know about before: starting with a given quadrilateral , construct a square on each side. Van Aubel's theorem states that the two line segments between the centers of opposite squares are of equal lengths and are at right angles to one another. See the picture below:

                                  

It is interesting that this is true not only for a convex quadrilateral, but for arbitrary one and even self-intersecting one. This post is devoted to proving this result in Maple. The proof was very short and simple. Note that the coordinates of points are given not by lists, but by vectors. This is convenient when using  LinearAlgebra:-DotProduct  and  LinearAlgebra:-Norm  commands.

The code of the proof (the notation of the points on the picture coincide with their names in the code):

restart;
with(LinearAlgebra):
assign(seq(P||i=<x[i],y[i]>, i=1..4)):
P||5:=P||1:
assign(seq(Q||i=(P||i+P||(i+1))/2+<0,1; -1,0>.(P||(i+1)-P||i)/2, i=1..4)):
expand(DotProduct((Q||3-Q||1),(Q||4-Q||2), conjugate=false));
is(Norm(Q||3-Q||1, 2)=Norm(Q||4-Q||2, 2));

The output:

                                                      0
                                                    true

 

VA.mw

This post is devoted to the rigorous proof of Miquel's five circles theorem, which I learned about from this question. The proof is essentially very simple and takes only 15 lines of code. The figure below, in which all the labels coincide with the corresponding names in the code, illustrates the basic ideas of the code. First, we symbolically define common points of intersection of blue circles with a red unit circle  (these parameters  s1 .. s5  are the polar coordinates of these points). All other parameters of this configuration can be expressed through them. Then we find the centers  M  and  N  of two circles. Then we find the coordinates of the point  K  from the condition that  CK  is perpendicular to  MN . Then we find the point  and using the result obtained, we easily find the coordinates  of all the points  A1 .. A5. Then we find the coordinates of the point   P  as the point of intersection of the lines  A1A2  and  A3A4 . Finally, we verify that the point  P  lies on a circle with center at the point  N , which completes the proof.

                      

 

Below - the code of the proof. Note that the code does not use any special (in particular geometric) packages, only commands from the Maple kernel. I usually try any geometric problems to solve in this style, it is more reliable,  and often shorter.

restart;
t1:=s1/2+s2/2: t2:=s2/2+s3/2:
M:=[cos(t1),sin(t1)]: N:=[cos(t2),sin(t2)]:
C:=[cos(s2),sin(s2)]: K:=(1-t)*~M+t*~N:
CK:=K-C: MN:=M-N:
t0:=simplify(solve(CK[1]*MN[1]+CK[2]*MN[2]=0, t)):
E:=combine(simplify(C+2*eval(CK,t=t0))):
s0:=s5-2*Pi: s6:=s1+2*Pi:
assign(seq(A||i=eval(E,[s2=s||i,s1=s||(i-1),s3=s||(i+1)]), i=1..5)):
Dist:=(p,q)->sqrt((p[1]-q[1])^2+(p[2]-q[2])^2):
LineEq:=(P,Q)->(y-P[2])*(Q[1]-P[1])=(x-P[1])*(Q[2]-P[2]):
Line1:=LineEq(A1,A2):
Line2:=LineEq(A3,A4):
P:=combine(simplify(solve({Line1,Line2},[x,y])))[]:
Circle:=(x-N[1])^2+(y-N[2])^2-Dist(N,C)^2:
is(eval(Circle, P)=0);  
# The final result

                                                                    true


It may seem that this proof is easy to repeat manually. But this is not so. Maple brilliantly coped with very cumbersome trigonometric transformations. Look at the coordinates of point  , expressed through the initial parameters  s1 .. s5 :

simplify(eval([x,y], P));  # The coordinates of the point  P

  

  

 

ProofMiquel.mw

A few days ago, I drew attention to the question in which OP talked about the generation of triangles in a plane, for which the lengths of all sides, the area and radius of the inscribed circle are integers. In addition, all vertices must have different integer coordinates (6 different integers), the lengths of all sides are different and the triangles should not be rectangular. I prepared the answer to this question, but the question disappeared somewhere, so I designed my answer as a separate post.

The triangles in the plane, for which the lengths of all sides and the area  are integers, are called as Heronian triangles. See this very interesting article in the wiki about such triangles
https://en.wikipedia.org/wiki/Integer_triangle#Heronian_triangles

The procedure finds all triangles (with the fulfillment of all conditions above), for which the lengths of the two sides are in the range  N1 .. N2 . The left side of the range is an optional parameter (by default  N1=5). It is not recommended to take the length of the range more than 100, otherwise the operating time of the procedure will greatly increase. The procedure returns the list in which each triangle is represented by a list of  [list of coordinates of the vertices, area, radius of the inscribed circle, list of lengths of the sides]. Without loss of generality, one vertex coincides with the origin (obviously, by a shift it is easy to place it at any point). 

The procedure works as follows: one vertex at the origin, then the other two must lie on circles with integer and different radii  x^2+y^2=r^2. Using  isolve  command, we find all integer points on these circles, and then in the for loops we select the necessary triangles.


 

 

restart;
HeronianTriangles:=proc(N2::posint,N1::posint:=5)
local k, r, S, L, Ch, Dist, IsOnline, c, P, p, A, B, C, a, b, s, ABC, cc, s1, T ;
uses combinat, geometry;
if N2<N1 then error "Should be N2>=N1" fi;
if N2<34 then return [] fi;
k:=0:
for r from max(N1,5) to N2 do
S:=[isolve(x^2+y^2=r^2)];
if nops(S)>4 then k:=k+1; L[k]:=select(s->s[1]<>0 and s[2]<>0,map(t->rhs~(convert(t,list)), S)); fi;
od:
L:=convert(L, list):
if type(L[1],symbol) then return [] fi;

Ch:=combinat:-choose([$1..nops(L)], 2):
Dist:=(A::list,B::list)->simplify(sqrt((A[1]-B[1])^2+(A[2]-B[2])^2));
IsOnline:=(A::list,B::list)->`if`(A[1]*B[2]-A[2]*B[1]=0, true, false);
k:=0:
for c in Ch do
for A in L[c[1]] do
for B in L[c[2]] do
if not IsOnline(A,B) and nops({A[],B[]})=4 then if type(Dist(A,B),posint) then
 k:=k+1; P[k]:=[A,B] fi; fi;
od: od: od:
P:=convert(P, list):
if type(P[1],symbol) then return [] fi;

k:=0:
for p in P do
point('A',0,0), point('B',p[1]), point('C',p[2]);
a:=simplify(distance('A','B')); b:=simplify(distance('A','C')); c:=simplify(distance('B','C'));
s:=sort([a,b,c]); s1:={a,b,c};
triangle(ABC,['A','B','C']);
incircle(cc,ABC);
r:=radius(cc);
if type(r,integer) and s[3]^2<>s[1]^2+s[2]^2 and nops(s1)=3 then k:=k+1; T[k]:=[[[0,0],p[]],area(ABC),r, [a,b,c]] fi;
od:
T:=convert(T,list);
if type(T[1],symbol) then return [] fi;
T;
end proc:

Examples of use of the procedure  HeronianTriangles

T:=HeronianTriangles(100): # All the Geronian triangles, whose lengths of two sides do not exceed 100
nops(T);

256

(1)

Tp:=select(p->p[1,2,1]>0 and p[1,2,2]>0 and p[1,3,1]>0 and p[1,3,2]>0, T);

[[[[0, 0], [16, 30], [28, 21]], 252, 6, [34, 35, 15]], [[[0, 0], [30, 16], [21, 28]], 252, 6, [34, 35, 15]], [[[0, 0], [21, 28], [15, 36]], 168, 4, [35, 39, 10]], [[[0, 0], [28, 21], [36, 15]], 168, 4, [35, 39, 10]], [[[0, 0], [27, 36], [13, 84]], 900, 10, [45, 85, 50]], [[[0, 0], [36, 27], [84, 13]], 900, 10, [45, 85, 50]], [[[0, 0], [33, 44], [48, 36]], 462, 7, [55, 60, 17]], [[[0, 0], [44, 33], [36, 48]], 462, 7, [55, 60, 17]], [[[0, 0], [33, 44], [96, 28]], 1650, 15, [55, 100, 65]], [[[0, 0], [44, 33], [28, 96]], 1650, 15, [55, 100, 65]], [[[0, 0], [16, 63], [72, 21]], 2100, 20, [65, 75, 70]], [[[0, 0], [63, 16], [21, 72]], 2100, 20, [65, 75, 70]], [[[0, 0], [39, 52], [18, 80]], 1092, 12, [65, 82, 35]], [[[0, 0], [52, 39], [80, 18]], 1092, 12, [65, 82, 35]], [[[0, 0], [32, 60], [56, 42]], 1008, 12, [68, 70, 30]], [[[0, 0], [60, 32], [42, 56]], 1008, 12, [68, 70, 30]], [[[0, 0], [42, 56], [30, 72]], 672, 8, [70, 78, 20]], [[[0, 0], [56, 42], [72, 30]], 672, 8, [70, 78, 20]]]

(2)

Tr:=map(p->p+[2,1],Tp[1,1]);
with(geometry):
point(A,Tr[1]), point(B,Tr[2]), point(C,Tr[3]):
triangle(ABC,[A,B,C]):
simplify(distance(A,B)), simplify(distance(A,C)), simplify(distance(B,C));
local O:
incircle(c,ABC, centername=O):
draw([A,B,C, ABC, c(color=blue)], color=red, thickness=2, symbol=solidcircle, tickmarks = [spacing(1)$2], gridlines, scaling=constrained, view=[0..31,0..33], size=[800,550], printtext=true, font=[times, 18], axesfont=[times, 10]);

[[2, 1], [18, 31], [30, 22]]

 

34, 35, 15

 

 



Examples of triangles with longer sides

T:=HeronianTriangles(1000,980):  # All the Geronian triangles, whose lengths of two sides lie in the range  980..1000
nops(T);

56

(3)

Tp:=select(p->p[1,2,1]>0 and p[1,2,2]>0 and p[1,3,1]>0 and p[1,3,2]>0, T);  # Triangles lying in the first quarter x>0, y>0
nops(%);

[[[[0, 0], [540, 819], [680, 714]], 85680, 80, [981, 986, 175]], [[[0, 0], [819, 540], [714, 680]], 85680, 80, [981, 986, 175]], [[[0, 0], [216, 960], [600, 800]], 201600, 168, [984, 1000, 416]], [[[0, 0], [960, 216], [800, 600]], 201600, 168, [984, 1000, 416]], [[[0, 0], [380, 912], [324, 945]], 31806, 31, [988, 999, 65]], [[[0, 0], [912, 380], [945, 324]], 31806, 31, [988, 999, 65]], [[[0, 0], [594, 792], [945, 324]], 277992, 216, [990, 999, 585]], [[[0, 0], [792, 594], [324, 945]], 277992, 216, [990, 999, 585]]]

 

8

(4)

 


 

Download Integer_Triangle1.mw

Edit.

1 2 3 4 5 6 7 Last Page 3 of 12