Kitonum

19511 Reputation

26 Badges

14 years, 324 days

MaplePrimes Activity


These are Posts that have been published by Kitonum

In this post, the Numbrix Puzzle is solved by the branch and bound method (see the details of this puzzle in  https://www.mapleprimes.com/posts/210643-Solving-A-Numbrix-Puzzle-With-Logic). The main difference from the solution using the  Logic  package is that here we get not one but all possible solutions. In the case of a unique solution, the  NumbrixPuzzle procedure is faster than the  Numbrix  one (for convenience, I inserted the code for Numbrix procedure into the worksheet below). In the case of many solutions, the  Numbrix  procedure is usually faster (see all the examples below).

 

restart;

NumbrixPuzzle:=proc(A::Matrix)
local A1, L, N, S, MS, OneStepLeft, OneStepRight, F1, F2, m, L1, p, q, a, b, T, k, s1, s, H, n, L2, i, j, i1, j1, R;
uses ListTools;
S:=upperbound(A); N:=nops(op(A)[3]); MS:=`*`(S);
A1:=convert(A, listlist);
for i from 1 to S[1] do
for j from 1 to S[2] do
for i1 from i to S[1] do
for j1 from 1 to S[2] do
if A1[i,j]<>0 and A1[i1,j1]<>0 and abs(A1[i,j]-A1[i1,j1])<abs(i-i1)+abs(j-j1) then return `no solutions` fi;
od; od; od; od;
L:=sort(select(e->e<>0, Flatten(A1)));
L1:=[`if`(L[1]>1,seq(L[1]-k, k=0..L[1]-2),NULL)];
L2:=[seq(seq(`if`(L[i+1]-L[i]>1,L[i]+k,NULL),k=0..L[i+1]-L[i]-2), i=1..nops(L)-1), `if`(L[-1]<MS,seq(L[-1]+k,k=0..MS-L[-1]-1),NULL)];
  

OneStepLeft:=proc(A1::listlist)
local s, M, m, k, T;
uses ListTools;
s:=Search(a, Matrix(A1));   
M:=[[s[1]-1,s[2]],[s[1]+1,s[2]],[s[1],s[2]-1],[s[1],s[2]+1]];
T:=table(); k:=0;
for m in M do
if m[1]>=1 and m[1]<=S[1] and m[2]>=1 and m[2]<=S[2] and A1[op(m)]=0 then k:=k+1; T[k]:=subsop(m=a-1,A1);
fi;
od;
convert(T, list);
end proc;

 
OneStepRight:=proc(A1::listlist)
local s, M, m, k, T, s1;
uses ListTools;
s:=Search(a, Matrix(A1));  s1:=Search(a+2, Matrix(A1));  
M:=[[s[1]-1,s[2]],[s[1]+1,s[2]],[s[1],s[2]-1],[s[1],s[2]+1]];
T:=table(); k:=0;
for m in M do
if m[1]>=1 and m[1]<=S[1] and m[2]>=1 and m[2]<=S[2] and A1[op(m)]=0 and `if`(a+2 in L, `if`(is(abs(s1[1]-m[1])+abs(s1[2]-m[2])>1),false,true),true) then k:=k+1; T[k]:=subsop(m=a+1,A1);
fi;
od;
convert(T, list);   
end proc;

F1:=LM->ListTools:-FlattenOnce(map(OneStepLeft, LM));
F2:=LM->ListTools:-FlattenOnce(map(OneStepRight, LM));

T:=[A1];
for a in L1 do
T:=F1(T);
od;

for a in L2 do
T:=F2(T);
od;

R:=map(t->convert(t,Matrix), T);
if nops(R)=0 then return `no solutions` else R[] fi;

end proc:

Numbrix := proc( M :: ~Matrix, { inline :: truefalse := false } )

local S, adjacent, eq, i, initial, j, k, kk, m, n, one, single, sol, unique, val, var, x;

    (m,n) := upperbound(M);

    initial := &and(seq(seq(ifelse(M[i,j] = 0
                                   , NULL
                                   , x[i,j,M[i,j]]
                                  )
                            , i = 1..m)
                        , j = 1..n));

    adjacent := &and(seq(seq(seq(x[i,j,k] &implies &or(NULL
                                                       , ifelse(i>1, x[i-1, j, k+1], NULL)
                                                       , ifelse(i<m, x[i+1, j, k+1], NULL)
                                                       , ifelse(j>1, x[i, j-1, k+1], NULL)
                                                       , ifelse(j<n, x[i, j+1, k+1], NULL)
                                                      )
                                 , i = 1..m)
                             , j = 1..n)
                         , k = 1 .. m*n-1));

    one := &or(seq(seq(x[i,j,1], i=1..m), j=1..n));   


    single := &not(&or(seq(seq(seq(seq(x[i,j,k] &and x[i,j,kk], kk = k+1..m*n), k = 1..m*n-1)
                                , i = 1..m), j = 1..n)));

    sol := Logic:-Satisfy(&and(initial, adjacent, one, single));
    
    if sol = NULL then
        error "no solution";
    end if;
if inline then
        S := M;
     else
        S := Matrix(m,n);
    end if;

    for eq in sol do
        (var, val) := op(eq);
        if val then
            S[op(1..2, var)] := op(3,var);
        end if;
    end do;
    S;
end proc:

           Two simple examples

A:=<0,0,5; 0,0,0; 0,0,9>;
# The unique solution
NumbrixPuzzle(A);

A:=<0,0,5; 0,0,0; 0,8,0>;
# 4 solutions
NumbrixPuzzle(A);

Matrix(3, 3, {(1, 1) = 0, (1, 2) = 0, (1, 3) = 5, (2, 1) = 0, (2, 2) = 0, (2, 3) = 0, (3, 1) = 0, (3, 2) = 0, (3, 3) = 9})

 

Matrix(3, 3, {(1, 1) = 3, (1, 2) = 4, (1, 3) = 5, (2, 1) = 2, (2, 2) = 7, (2, 3) = 6, (3, 1) = 1, (3, 2) = 8, (3, 3) = 9})

 

Matrix(3, 3, {(1, 1) = 0, (1, 2) = 0, (1, 3) = 5, (2, 1) = 0, (2, 2) = 0, (2, 3) = 0, (3, 1) = 0, (3, 2) = 8, (3, 3) = 0})

 

Matrix(%id = 18446746210121682686), Matrix(%id = 18446746210121682806), Matrix(%id = 18446746210121674750), Matrix(%id = 18446746210121674870)

(1)


Comparison with Numbrix procedure. The example is taken from
http://rosettacode.org/wiki/Solve_a_Numbrix_puzzle 

 A:=<0, 0, 0, 0, 0, 0, 0, 0, 0;
 0, 0, 46, 45, 0, 55, 74, 0, 0;
 0, 38, 0, 0, 43, 0, 0, 78, 0;
 0, 35, 0, 0, 0, 0, 0, 71, 0;
 0, 0, 33, 0, 0, 0, 59, 0, 0;
 0, 17, 0, 0, 0, 0, 0, 67, 0;
 0, 18, 0, 0, 11, 0, 0, 64, 0;
 0, 0, 24, 21, 0, 1, 2, 0, 0;
 0, 0, 0, 0, 0, 0, 0, 0, 0>;
CodeTools:-Usage(NumbrixPuzzle(A));
CodeTools:-Usage(Numbrix(A));

Matrix(9, 9, {(1, 1) = 0, (1, 2) = 0, (1, 3) = 0, (1, 4) = 0, (1, 5) = 0, (1, 6) = 0, (1, 7) = 0, (1, 8) = 0, (1, 9) = 0, (2, 1) = 0, (2, 2) = 0, (2, 3) = 46, (2, 4) = 45, (2, 5) = 0, (2, 6) = 55, (2, 7) = 74, (2, 8) = 0, (2, 9) = 0, (3, 1) = 0, (3, 2) = 38, (3, 3) = 0, (3, 4) = 0, (3, 5) = 43, (3, 6) = 0, (3, 7) = 0, (3, 8) = 78, (3, 9) = 0, (4, 1) = 0, (4, 2) = 35, (4, 3) = 0, (4, 4) = 0, (4, 5) = 0, (4, 6) = 0, (4, 7) = 0, (4, 8) = 71, (4, 9) = 0, (5, 1) = 0, (5, 2) = 0, (5, 3) = 33, (5, 4) = 0, (5, 5) = 0, (5, 6) = 0, (5, 7) = 59, (5, 8) = 0, (5, 9) = 0, (6, 1) = 0, (6, 2) = 17, (6, 3) = 0, (6, 4) = 0, (6, 5) = 0, (6, 6) = 0, (6, 7) = 0, (6, 8) = 67, (6, 9) = 0, (7, 1) = 0, (7, 2) = 18, (7, 3) = 0, (7, 4) = 0, (7, 5) = 11, (7, 6) = 0, (7, 7) = 0, (7, 8) = 64, (7, 9) = 0, (8, 1) = 0, (8, 2) = 0, (8, 3) = 24, (8, 4) = 21, (8, 5) = 0, (8, 6) = 1, (8, 7) = 2, (8, 8) = 0, (8, 9) = 0, (9, 1) = 0, (9, 2) = 0, (9, 3) = 0, (9, 4) = 0, (9, 5) = 0, (9, 6) = 0, (9, 7) = 0, (9, 8) = 0, (9, 9) = 0})

 

memory used=7.85MiB, alloc change=-3.01MiB, cpu time=172.00ms, real time=212.00ms, gc time=93.75ms

 

Matrix(9, 9, {(1, 1) = 49, (1, 2) = 50, (1, 3) = 51, (1, 4) = 52, (1, 5) = 53, (1, 6) = 54, (1, 7) = 75, (1, 8) = 76, (1, 9) = 81, (2, 1) = 48, (2, 2) = 47, (2, 3) = 46, (2, 4) = 45, (2, 5) = 44, (2, 6) = 55, (2, 7) = 74, (2, 8) = 77, (2, 9) = 80, (3, 1) = 37, (3, 2) = 38, (3, 3) = 39, (3, 4) = 40, (3, 5) = 43, (3, 6) = 56, (3, 7) = 73, (3, 8) = 78, (3, 9) = 79, (4, 1) = 36, (4, 2) = 35, (4, 3) = 34, (4, 4) = 41, (4, 5) = 42, (4, 6) = 57, (4, 7) = 72, (4, 8) = 71, (4, 9) = 70, (5, 1) = 31, (5, 2) = 32, (5, 3) = 33, (5, 4) = 14, (5, 5) = 13, (5, 6) = 58, (5, 7) = 59, (5, 8) = 68, (5, 9) = 69, (6, 1) = 30, (6, 2) = 17, (6, 3) = 16, (6, 4) = 15, (6, 5) = 12, (6, 6) = 61, (6, 7) = 60, (6, 8) = 67, (6, 9) = 66, (7, 1) = 29, (7, 2) = 18, (7, 3) = 19, (7, 4) = 20, (7, 5) = 11, (7, 6) = 62, (7, 7) = 63, (7, 8) = 64, (7, 9) = 65, (8, 1) = 28, (8, 2) = 25, (8, 3) = 24, (8, 4) = 21, (8, 5) = 10, (8, 6) = 1, (8, 7) = 2, (8, 8) = 3, (8, 9) = 4, (9, 1) = 27, (9, 2) = 26, (9, 3) = 23, (9, 4) = 22, (9, 5) = 9, (9, 6) = 8, (9, 7) = 7, (9, 8) = 6, (9, 9) = 5})

 

memory used=1.21GiB, alloc change=307.02MiB, cpu time=37.00s, real time=31.88s, gc time=9.30s

 

Matrix(%id = 18446746210094669942)

(2)


In the example below, which has 104 solutions, the  Numbrix  procedure is faster.

C:=Matrix(5,{(1,1)=1,(5,5)=25});
CodeTools:-Usage(NumbrixPuzzle(C)):
nops([%]);
CodeTools:-Usage(Numbrix(C)):

Matrix(5, 5, {(1, 1) = 1, (1, 2) = 0, (1, 3) = 0, (1, 4) = 0, (1, 5) = 0, (2, 1) = 0, (2, 2) = 0, (2, 3) = 0, (2, 4) = 0, (2, 5) = 0, (3, 1) = 0, (3, 2) = 0, (3, 3) = 0, (3, 4) = 0, (3, 5) = 0, (4, 1) = 0, (4, 2) = 0, (4, 3) = 0, (4, 4) = 0, (4, 5) = 0, (5, 1) = 0, (5, 2) = 0, (5, 3) = 0, (5, 4) = 0, (5, 5) = 25})

 

memory used=0.94GiB, alloc change=-22.96MiB, cpu time=12.72s, real time=11.42s, gc time=2.28s

 

104

 

memory used=34.74MiB, alloc change=0 bytes, cpu time=781.00ms, real time=783.00ms, gc time=0ns

 

 


 

Download NumbrixPuzzle.mw

Here is a classic puzzle:
You are camping, and have an 8-liter bucket which is full of fresh water. You need to share this water fairly into exactly two portions (4 + 4 liters). But you only have two empty buckets: a 5-liter and a 3-liter. Divide the 8 liters in half in as short a time as possible.

This is not an easy task and requires at least 7 transfusions to solve it. 

The procedure  Pouring  solves a similar problem for the general case. Given n buckets of known volume and the amount of water in each bucket is known. Buckets can be partially filled, be full or be empty (of course the case when all is empty or all is full is excluded). With each individual transfusion, the bucket from which it is poured must be completely free, or the bucket into which it is poured must be completely filled. It is forbidden to pour water anywhere other than the indicated buckets.

Formal parameters of the procedure: BucketVolumes  is a list of bucket volumes,  WaterVolumes  is a list of water volumes in these buckets. The procedure returns all possible states that can occur during transfusions in the form of a tree (the initial state  WaterVolumes  is its root).

restart;
Pouring:=proc(BucketVolumes::list(And(positive,{integer,float,fraction})),WaterVolumes::list(And(nonnegative,{integer,float,fraction})), Output:=graph)
local S, W, n, N, OneStep, j, v, H, G;
uses ListTools, GraphTheory;

n:=nops(BucketVolumes); 
if nops(WaterVolumes)<>n then error "The lists should be the same length" fi;
if n<2 then error "Must have at least 2 buckets" fi;
if not `or`(op(WaterVolumes>~0)) then error "There must be at least one non-empty bucket" fi;
if BucketVolumes=WaterVolumes then error "At least one bucket should not be full" fi;
if not `and`(seq(WaterVolumes[i]<=BucketVolumes[i], i=1..n)) then error "The amount of water in each bucket cannot exceed its volume" fi;
W:=[[WaterVolumes]];

OneStep:=proc(W)
local w, k, i, v, V, k1, v0;
global L;
L:=convert(Flatten(W,1), set);
k1:=0; 
for w in W do
k:=0; v:='v';
for i from 1 to n do
for j from 1 to n do
if i<>j and w[-1][i]<>0 and w[-1][j]<BucketVolumes[j] then k:=k+1; v[k]:=subsop(i=`if`(w[-1][i]<=BucketVolumes[j]-w[-1][j],0,w[-1][i]-(BucketVolumes[j]-w[-1][j])),j=`if`(w[-1][i]<=BucketVolumes[j]-w[-1][j],w[-1][j]+w[-1][i],BucketVolumes[j]),w[-1]); fi;
od; 
od; 
v:=convert(v,list);
if `and`(seq(v0 in L, v0=v)) then k1:=k1+1; V[k1]:=w else 
for v0 in v do  
if not (v0 in L) then k1:=k1+1; V[k1]:=[op(w),v0] fi;
od;
fi;
L:=L union convert(v,set);
od;
convert(V,list);
end proc:

S[0]:={};
for N from 1 do
H[N]:=(OneStep@@N)(W);
S[N]:=L;
if S[N-1]=S[N] then break fi;
od;
if Output=set then return L else
if Output=trails then interface(rtablesize=infinity);
return <H[N-1]> else
G:=Graph(seq(Trail(map(t->t[2..-2],convert~(h,string))),h=H[N-1]));
DrawGraph(G, style=tree, root=convert(WaterVolumes,string)[2..-2], stylesheet = [vertexcolor = "Yellow", vertexfont=[TIMES,20]], size=[800,500])  fi; fi;

end proc: 

Examples of use:

Here is the solution to the original puzzle above. We see that at least 7 transfusions are  
required to get equal volumes (4 + 4) in two buckets

Pouring([8,5,3], [8,0,0]);
           

 

 With an increase in the number of buckets, the number of solutions is extremely 
 increased. Here is the solution to the problem: is it possible to equalize the amount of water (7+7+7+7) in the following example? 

Pouring([14,10,9,9],[14,10,4,0]);
S:=Pouring([14,10,9,9],[14,10,4,0], set);
is([7,7,7,7] in S);
nops(S);
         

 

 

Download Pouring.mw

 

 

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

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