4
$\begingroup$

I'm currently interested in looking at all the 6-digit numbers in base 10 (i.e. 000000 to 99999), and checking if they satisfy at least one of 5 equations: either 0=0+0+0+0+0, 0+0=0+0+0+0,..., 0+0+0+0+0=0. In terms of visualization, I created a 1000x1000 square, starting with 000000 to 000999 in the top left, then 001000 to 001999 a row below, etc., until the last row, which is 999000 to 999999. The coloring rule is simple: if there are any solutions to the 6-digit number, I color the square black, and if there are none, I color the square white.

I had used Excel when I was looking at 4-digit numbers, but it has quickly become unwieldy and difficult to visualize when extending to 6-digit numbers (it seems there is a limit to the size of the pixels). Is there a way I can do this in Mathematica? I am also looking at solutions modulo n=1,2,3,..., so I am not sure if it is possible to include a changing mod in the code.

$\endgroup$
3
  • 2
    $\begingroup$ ArrayPlot or MatrixPlot? $\endgroup$
    – Michael E2
    Commented Jun 20, 2022 at 16:48
  • 2
    $\begingroup$ To get you started: MatrixPlot[Partition[Range[0, 1*^6 - 1], 1000]]; what you now need to do is to write a function that tests your criterion on all those integers (IntegerDigits[] should be helpful). I get an image like this when I tried it out. $\endgroup$ Commented Jun 20, 2022 at 16:51
  • 2
    $\begingroup$ A trick that circumvents Mathematica graphics and goes straight to pixel-exact image storage (for easy external visualization): Export["something.png", Image[Raster[Reverse[A]]]], assuming that A is a matrix with values 0 and 1. $\endgroup$
    – Roman
    Commented Jun 20, 2022 at 16:52

1 Answer 1

2
$\begingroup$

This looks like a fun problem, from an optimisation point of view in particular.

It should highlight the how you can make your code faster by looking at different ways to solve your problem and also making use of Compile

The first "optimisation" is that we do not need to work on a square yet, we operate on a list and convert to a square later with ArrayReshape. ie: we define the "square" using Range.

I am also looking at solutions modulo n=1,2,3,...,

I'm, not sure where modulo n comes in so I will ignore this for now.


The naive approach:

(* init *)
$base = 10;
$squareSize = 1000;

(* define square *)
$lastNumber = $squareSize^2 - 1;
$digitPadLength = IntegerLength[$lastNumber, $base];
$numberOfEquations = IntegerLength[$lastNumber, $base] - 1;

(* create the square *)
numbers = Range[0, $lastNumber]; (* the "square" *)
digitsList = IntegerDigits[#, $base, $digitPadLength] & /@ numbers;

(* equations *)
equationsToSatisfy =
  Function[
   digits,
   Table[
    Total[digits[[;; i]]] == Total[digits[[i + 1 ;;]]]
    ,
    {i, $numberOfEquations}
    ]];

(* run equations *)
equationResult = equationsToSatisfy /@ digitsList;
isAtLeastOneEquationSatisfied = Or @@@ equationResult;
isAtLeastOneEquationSatisfiedBoolean = Boole@isAtLeastOneEquationSatisfied;

(* plot result *)
satisfactionMatrix = ArrayReshape[isAtLeastOneEquationSatisfiedBoolean, {$squareSize, $squareSize}];
ArrayPlot@satisfactionMatrix

looks similar to @J.M.'s which is nice


Looking at the timings it is clear that the function equationsToSatisfy is slow

AbsoluteTiming[
 equationResult = equationsToSatisfy /@ digitsList;
 isAtLeastOneEquationSatisfied = Or @@@ equationResult;
 ] (*{38.3839, Null}*)

Given the logical and numerical nature of the problem, we should be able to get quicker. First let's see if we can reformulate equationsToSatisfy by examining what the equations are checking.

Using 60231 as an example:

digits = {0, 6, 0, 2, 3, 1};

Total[digits] (* 12 *)

Table[
 {Total[digits[[;; i]]], Total[digits[[i + 1 ;;]]]}
 ,
 {i, $numberOfEquations}
 ] (* {{0, 12}, {6, 6}, {6, 6}, {8, 4}, {11, 1}} *)

equationsToSatisfy[digits] (* {False, True, True, False, False} *)

Some observations

  1. If the total of the digits its an odd number we can never satisfy the fact that the LHS is equal to RHS. We can then use OddQ to terminate the code early.

  2. For the LHS to be equal to RHS it must be half the total.

    • This implies we don't need to calculate to RHS and the LHS is enough, ie if the LHS is half the total then the RHS is too by symmetry
    • We can use the function Accumulate to help us
accumulation = Accumulate[digits] (* {0, 6, 6, 8, 11, 12} *)

Hence we just need to check if 6 is a member of this list then at least one of the equations has been satisfied.

(Strictly speaking we should check if 6 is a member of {0,6,6,8,11} ie Most[accumulation], however we can save one step of computation by noting that the last element of accumulation is the total of digits, hence checking for MemberQ[accumulation,totalOfDigits/2] is equivalent to checking for MemberQ[Most[accumulation],totalOfDigits/2).

Using these two observations we rewrite equationsToSatisfy and get a 3-4x speed up:

The smarter approach:

equationsToSatisfy2[digits_] :=
  Module[
   {
    accumulation = Accumulate@digits,
    totalOfDigits
    },
   totalOfDigits = Last@accumulation;
   If[OddQ[totalOfDigits], Return[False]]; (* from observation 1 *)
   MemberQ[accumulation, totalOfDigits/2] (* from observation 2 *)
   ];


AbsoluteTiming[
 equationResult2 = equationsToSatisfy2 /@ digitsList;
 ] (* {11.6489, Null} *)

isAtLeastOneEquationSatisfied == equationResult2 (* True *)

Finally noting that the problem is highly numerical and embarrassingly parallel hints that it would benefit from Compile and being made {Listable}. In fact we get a 60x speed up from the approach above and 200x from our initial naive approach

equationsToSatisfy2C =
 Compile[
  {{digits, _Real, 1}},
  Module[
   {
    accumulation = Accumulate@digits,
    totalOfDigits
    },
   totalOfDigits = Last@accumulation;
   If[OddQ[totalOfDigits], Return[False]];
   MemberQ[accumulation, totalOfDigits/2]
   ],
  CompilationTarget -> "C",
  RuntimeAttributes -> {Listable},
  RuntimeOptions -> "Speed"
  ];

AbsoluteTiming[
 equationResult2C = equationsToSatisfy2C[digitsList];
 ] (* {0.192282, Null} *)

isAtLeastOneEquationSatisfied == equationResult2C (* True *)

I'm pretty sure we can solve this problem even faster and would welcome other suggestions :)

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.