####################################################################################### with(ListTools): Help := proc() if args=NULL then print(" Part 1: minAllST, minAllRecursiveST "); print(" Part 2: Ord, AllTrip, EachTrip, Trip1, Trip2, Trip3, Trip4"); print(" Part 3: Zeil, RunColor, FindLength, Case1, Case2 "); print(" Case3, Case4, FindMin, FindT1T2, FindCase3 "); print(" For more info: type 'Help(Part1)' or 'Help(Part2)' or 'Help(Part3)' "); elif nargs=1 and args[1]=Part1 then print(" Part 1 is directly search for the pattern of the minimum Schur Triples"); print(" The main function is minAllST(n,r)"); print(" input : length of interval, n and color, r"); print(" output : list of r-coloring of length n with minimum Schur Triples"); print(" Try: minAllST(10,3)"); elif nargs=1 and args[1]=Part2 then print(" Part 2 is the numerical counting number of Schur Triples "); print(" The main function is Ord(C,L,n) and AllTrip(C,L)"); print(" Ord(C,L,n)"); print(" input : list of color, list of lengths correspond to color, symbol n"); print(" output : the number of Schur Triples of the list above in term of n"); print(" Try: Ord([1,2,1],[4,6,1],n)"); print(" AllTrip(C,L)"); print(" input : list of color, list of lengths correspond to color"); print(" output : the number of Schur Triples of the list above in term of k"); print(" Try: Ord([1,2,1],[4,6,1])"); print(" Note that [4,6,1] means the lengths are 4k , 6k and k repectively"); elif nargs=1 and args[1]=Part3 then print(" Part 3 is calculating the minimum number of monochromatic Schur Triples"); print(" Using greedy algorithm and calculus"); print(" The main function is Zeil(r)"); print(" input : color, r"); print(" output : list and order of r-coloring with minimum Schur Triples"); print(" Try: Zeil(2)"); else print("Sorry no help for ", args[1]); fi: end: ############################################################################### ############################################################################### # Part 1) # Main Function - Min All Schur Triple # Input: length of interval, number of color # Output: minimum value, set of minimum interval # minAllST := proc(n, r) option remember; global minValue, minInterval; minValue := n*(n-1)/2; minInterval := []; minAllRecursiveST(n, r, [1], 0); return [minValue, minInterval]; end: # # Min All Recursive Schur Triple # minAllRecursiveST := proc(n, r, interval, count) global minValue, minInterval; local i, size, count2; size := nops(interval); count2 := count; for i from 1 to size/2 do if interval[i] = interval[size] and interval[i] = interval[size-i] then count2 := count2+1; fi; od; if count2 > minValue then return; elif size = n then if count2 < minValue then minValue := count2; minInterval := [interval]; elif count2 = minValue then minInterval := [op(minInterval), interval]; fi; else for i from 1 to r do minAllRecursiveST(n, r, [op(interval), i], count2); od; fi; end: ############################################################################# ############################################################################# ## Part 2) ## Count Triple and give order. # Ord # input : list of color and length of each color # output : the order of minimum number of Schur Triples. Ord := proc(color,length,n) local nT, L; nT := AllTrip(color, length); L := add(length[i] , i = 1..nops(length)); return simplify(subs( k = n/L, nT)); end: # input: color := [1,2,1] , length := [4,6,1] AllTrip := proc(color, length) local i,j,k,numcol,Atrip, temp; numcol := max(op(color)); Atrip := 0; for i from 1 to numcol do temp := []; for j from 1 to nops(color) do if color[j] = i then temp := [op(temp), [add(length[k], k = 1..j-1),add(length[k], k = 1..j)]]; fi: od: Atrip := Atrip + EachTrip(temp); od: return(Atrip); end: # EachTrip :: # input list of interval # output number of triple EachTrip := proc(AA) local i,j,m,Etrip; Etrip := 0; for i from 1 to nops(AA) do Etrip := Etrip + Trip1(AA[i]); for j from i+1 to nops(AA) do Etrip := Etrip + Trip2(AA[i],AA[j])+ Trip3(AA[i],AA[j]); for m from j+1 to nops(AA) do Etrip := Etrip + Trip4(AA[i],AA[j],AA[m]); od: od: od: return simplify(Etrip); end: # Trip1: three points in T1 Trip1 := proc(T1) local a1,a2; a1 := T1[1]; a2 := T1[2]; return (max(a2-2*a1,0)^2/4*k^2); end: # Trip2 : one point in T1 , two point in T2. Trip2 := proc(T1,T2) local a1,a2,b1,b2,L; a1 := T1[1]; a2 := T1[2]; b1 := T2[1]; b2 := T2[2]; if not(a1>=0) or not(a1 <= a2) or not(a2 <= b1) or not(b1 <= b2) then ERROR(BadInput); fi: L := a2-a1; if b2 > b1 + a2 then return ((b2-a2-b1)*L + L^2/2)*k^2; else return max(b2-a1-b1,0)^2/2*k^2; fi: end: # Trip3 : two points in T1 , one point in T2. Trip3 := proc(T1,T2) local a1,a2,b1,b2,x,y,z,w; a1 := T1[1]; a2 := T1[2]; b1 := T2[1]; b2 := T2[2]; if not(a1>=0) or not(a1 <= a2) or not(a2 <= b1) or not(b1 <= b2) then ERROR(BadInput); fi: x := b1-a2; y := b1/2; z := b2-a2; w := b2/2; if y <= x or 2*a1 >= b2 then return(0); elif z < y then if a1 < z then return (((z-max(a1,x))*(max(a1-x,0)+z-x)/2 + (y-z)*(z-x)+ (z-x)^2/4)*k^2); elif a1 < y then return( ((y-a1)*(z-x)+ (z-x)^2/4)*k^2); elif a1 < w then return( (1/2*(w-a1)*2*(w-a1))*k^2); fi: elif z < w then if a1 < y then return ( ((y-max(a1,x))*(max(a1-x,0)+y-x)/2 +1/2*(z-y)*(3*y-2*x-z) + (2*y-x-z)^2/4)*k^2); elif a1 < z then return ( (1/2*(z-a1)*(2*y-x-z+2*y-x-a1)+(2*y-x-z)^2/4)*k^2); elif a1 < w then return ( ((w-a1)*2*(w-a1)/2)*k^2); fi: elif z >= w then if a1 < y then return (((y-max(a1,x))*(max(a1-x,0)+y-x)/2+(y-x)^2/2)*k^2); elif a1 < w then return ((2*y-x-a1)^2/2*k^2); fi: else ERROR("No case in Trip3."); fi: end: # Trip4 : one point in T1, one point in T2, one point in T3. Trip4 := proc(T1,T2,T3) local a1,a2,b1,b2,c1,c2,x,y,z,w; a1 := T1[1]; a2 := T1[2]; b1 := T2[1]; b2 := T2[2]; c1 := T3[1]; c2 := T3[2]; if not(a1>=0) or not(a1 <= a2) or not(a2 <= b1) or not(b1 <= b2) or not(b2 <= c1) or not(c1 <= c2) then ERROR(BadInput); fi: x := c1-a2; y := min(c1-a1,c2-a2); z := max(c1-a1,c2-a2); w := c2-a1; if w <= b1 or x >= b2 then return(0); elif b1 < y then if b2 < y then return(((b2-max(b1,x))*(max(b1-x,0)+b2-x)/2)*k^2); elif b2 < z then return(((y-max(b1,x))*(max(b1-x,0)+y-x)/2 + (y-x)*(b2-y))*k^2); else return(((y-max(b1,x))*(max(b1-x,0)+y-x)/2 + (y-x)*(z-y) + 1/2*(min(b2,w)-z)*(y-x+max(w-b2,0)))*k^2); fi: elif b1 < z then if b2 < z then return( (b2-b1)*(y-x)*k^2); else return( ((z-b1)*(y-x) + 1/2*(min(b2,w)-z)*(y-x+max(w-b2,0)))*k^2); fi: elif b1 < w then return( (1/2*(min(b2,w)-b1)*(w-b1+max(w-b2,0)))*k^2); else ERROR("No case in Trip4"); fi: end: ############################################################################ ############################################################################ # Part 3) # Use calculus and greedy algorithm to find the coloring that give the # minimum number of Schur Triples. # Zeil # input : color # output : the minimum coloring, length of each coloring, the order. Zeil := proc(r) local color , length , temp; color := [1]; length := [1]; temp := RunColor(color,length,r); while temp <> [color,length] do color := temp[1]; length := temp[2]; temp := RunColor(color,length,r); od: return(color, length, Ord(color, length,n)); end: # RunColor # input : color, length, r # output : the new minimum coloring, new length of each coloring. # of the smallest possible color (according to greedy algorithm). RunColor := proc(color, length, r) local i, temp; for i from 1 to r do temp := FindLength(color, length, i); if temp <> 0 then return([[op(color),i], [op(length),temp]]); fi: od: return([color,length]); end: # FindLength # input : color, length, index # output : the new optimal length (attaching at the end) # of the index coloring. FindLength := proc(color, length, index) local n,ret, more, newL; if nops(color) = 1 and index = 1 then return 0; fi: more := Case1(color, length, index, var) + Case2(color, length, index,var) + Case3(color, length, index,var) + Case4(color, length, index,var); newL := FindMin(color, length, index,more , var); ############################################## #check output, I check the output in FindMin; ############################################## return(newL); end: ############################################################################## ############################################################################## # Case 1 Case1 := proc(color, length, index, X) local long; long := add(length[i], i = 1..nops(length)); return(piecewise(X > long, (X-long)^2/4*k^2)); end: #Case 2 Case2 := proc(color, length, index, X) local a1,a2,b1,L; if not(member(index, {op(color)})) then return 0; fi: a1 := add(length[i], i = 1..Search(index, color)-1); a2 := add(length[i], i = 1..Search(index, color)); b1 := add(length[i], i = 1..nops(length)); L := a2-a1; return piecewise( X > a2,((X-a2)*L + L^2/2)*k^2, X >= 0, max(X-a1,0)^2/2*k^2); end: #Case 3 Case3 := proc(color, length, index, var) local a1,a2,b1,b2 ,x,y,z,w, endPos, listTT; listTT := FindCase3(color, length, index); if listTT = [] then return(0); fi: endPos := add(length[i], i = 1..nops(length)); a1 := listTT[1]; a2 := listTT[2]; b1 := endPos; b2 := endPos+var; x := b1-a2; y := b1/2; z := b2-a2; w := b2/2; if y <= x then return(0); fi: return piecewise( a1 <= x and z <= y, (1/2*(z-x)^2+(y-z)*(z-x)+ (z-x)^2/4)*k^2, a1 <= z and z <= y, (1/2*(z-a1)*(z-x+a1-x)+(y-z)*(z-x)+ (z-x)^2/4)*k^2, a1 <= y and z <= y, ((y-a1)*(z-x)+(z-x)^2/4)*k^2, a1 <= w and z <= y, ((w-a1)*2*(w-a1)/2)*k^2, a1 <= x and z <= w, ((y-x)^2/2 + 1/2*(z-y)*(y-x+2*y-x-z)+ 1/4*(2*y-x-z)^2)*k^2, a1 <= y and z <= w, (1/2*(y-a1)*(y-x+a1-x) +1/2*(z-y)*(y-x+2*y-x-z)+1/4*(2*y-x-z)^2)*k^2, a1 <= z and z <= w, (1/2*(z-a1)*(y-x-(a1-y)+2*y-x-z)+1/4*(2*y-x-z)^2)*k^2, a1 <= w and z <= w, ((w-a1)*2*(w-a1)/2)*k^2, a1 <= x and w < z, ((y-x)^2)*k^2, a1 <= y and w < z, (1/2*(y-a1)*(y-x+a1-x)+1/2*(y-x)^2)*k^2, a1 <= w and w < z, ((y+y-x-a1)^2/2)*k^2, 0); end: #Case4 Case4 := proc(color, length, index, var) local a1,a2,b1,b2,c1,c2,x,y,z,w ,i, listTT, f,g; listTT := FindT1T2( color, length, index); f := 0 ; for i from 1 to nops(listTT) do a1 := listTT[i][1][1]; a2 := listTT[i][1][2]; b1 := listTT[i][2][1]; b2 := listTT[i][2][2]; c1 := add(length[j], j = 1..nops(length)); c2 := c1+var; x := c1-a2; y := c2-a2; z := c1-a1; w := c2-a1; g := piecewise( b2 < x , 0, b2 < y and y < z, ((b2-x)^2/2)*k^2, b2 < z and y < z, ((y-x)^2/2 + (b2-y)*(y-x))*k^2 , b2 < w and y < z, ((y-x)^2/2 + (z-y)*(y-x) + (b2-z)*((y-x)+w-b2)/2)*k^2, y < z, ((y-x)^2/2 + (z-y)*(y-x) + (w-z)^2/2)*k^2, b2 < z , ((b2-x)^2/2)*k^2, b2 < y , ((z-x)^2/2 + (b2-z)*(z-x))*k^2 , b2 < w , ((z-x)^2/2 + (y-z)*(z-x) + (b2-y)*((z-x)+w-b2)/2)*k^2 , ((z-x)^2/2 + (y-z)*(z-x) + (w-y)^2/2)*k^2 ) - piecewise( b1 < x , 0, b1 < y and y < z, ((b1-x)^2/2)*k^2, b1 < z and y < z, ((y-x)^2/2 + (b1-y)*(y-x))*k^2 , b1 < w and y < z, ((y-x)^2/2 + (z-y)*(y-x) + (b1-z)*((y-x)+w-b1)/2)*k^2, y < z, ((y-x)^2/2 + (z-y)*(y-x) + (w-z)^2/2 )*k^2, b1 < z , ((b1-x)^2/2)*k^2, b1 < y , ((z-x)^2/2 + (b1-z)*(z-x))*k^2 , b1 < w , ((z-x)^2/2 + (y-z)*(z-x) + (b1-y)*((z-x)+w-b1)/2)*k^2 , ((z-x)^2/2 + (y-z)*(z-x) + (w-y)^2/2)*k^2); f := f + g; od: return f; end: # FindT1T2 for case 4. FindT1T2 := proc(color,length,index) local i, j, listCol, endPos, ret; endPos := add(length[i], i = 1..nops(length)); listCol := []; ret := []; #make list color for i from 1 to nops(color) do if color[i] = index then listCol := [op(listCol), [add(length[m], m = 1..i-1),add(length[m], m = 1..i)]]; fi: od: for i from 1 to nops(listCol) do for j from i +1 to nops(listCol) do if ( listCol[i][1] + listCol[j][1] <= endPos and listCol[i][2] + listCol[j][2] >= endPos) then ret := [op(ret),[listCol[i],listCol[j]]]; fi: od: od: return ret; end: FindCase3 := proc(color, length, index) local i, endPos ; endPos := add(length[i], i = 1..nops(length)); i := 1; while i < nops(color) and add(length[j],j=1..i) < endPos/2 do i := i +1; od: if color[i] = index then return( [add(length[j],j=1..i-1) ,add(length[j],j=1..i)] ); else return( [] ); fi: end: #FindMin FindMin := proc(color, length, index,more, var) local n, long, nT ,EQ, temp ; long := add(length[i], i = 1..nops(length)); nT := AllTrip(color,length) + more; assume(n >= 0); EQ := simplify(subs(k = n/(long + var),nT)) ; temp := minimize( EQ , var >= 0, location); #check output if (evalb( temp[1] = Ord( [op(color),index],[op(length),rhs(temp[2][1][1][1])],n))) = false then print("Counting triples meet something not expected"); ERROR(); fi: return rhs(temp[2][1][1][1]); end: