3
$\begingroup$

I want to do this,

Solve[ { n1 + n2 + n3 + .... = 10, n1 + 2 n2 + 3 n3 + .... = 100}, 
{n1,n2,n3,....}, Assumptions -> {n1,n2,n3,.... are Integers >=0} ]

Note that here .... means "ad infinitum" . Can someone write a genuine Mathematica code that solves this system?

$\endgroup$
  • $\begingroup$ Are you asking if Mathematica can calculate values for $n_1$, $n_2$, $n_3$, ... , $n_{\infty - 3}$, $n_{\infty - 2}$, $n_{\infty - 1}$, $n_{\infty}$? How exactly would you display the solutions or store them in memory? $\endgroup$ – MassDefect Feb 22 at 16:37
  • 1
    $\begingroup$ Have you seen FrobeniusSolve[]? $\endgroup$ – J. M. is away Feb 24 at 2:26
  • $\begingroup$ FrobeniusSolve[] is indeed a nice shortcut. Thanks. Now that I think about it the header is somewhat misleading. The reason why I felt I may have to declare an infinite array of integers is because I wanted to make the numbers 10 and 100 variables. Even then one can get by without having to declare an infinite array I guess. It is hard to guess that the problem I have described arises naturally while trying to compute the entropy of a physical system. $\endgroup$ – Quasar Supernova Feb 24 at 5:46
  • $\begingroup$ But I am unable to see how to Solve two Frobenius Equations Simultaneously.. $\endgroup$ – Quasar Supernova Feb 24 at 13:33
3
$\begingroup$

To solve the full problem, I don't know how to coax Mathematica into running 100 nested Do's in a reasonable amount of time. Instead, I would use Mathematica to generate C-code that is then executed externally. In this way, the the full problem can be solved in about 40 seconds.

Here's the Mathematica code that generates the C-code. Assume that we want $\sum_{i=1}^{i_{\text{max}}}n_i=\text{sum}_0=10$ and $\sum_{i=1}^{i_{\text{max}}}in_i=\text{sum}_1=100$, and go to the full problem with $i_{\text{max}}=100$:

sum0 = 10;
sum1 = 100;
imax = 100;
f = OpenWrite["~/Desktop/loop.c"]; (* or whereever you want to save it *)
WriteString[f, "#include <stdio.h>\n"];
WriteString[f, "int main() {\n"];
WriteString[f, "int n[" <> ToString[imax] <> "];\n"];
WriteString[f, "int s0[" <> ToString[imax + 1] <> "];\n"];
WriteString[f, "int s1[" <> ToString[imax + 1] <> "];\n"];
WriteString[f, "s0[0]=s1[0]=0;\n"];
Do[WriteString[f, "for (n[" <> ToString[i - 1] <> "]=0; (s0[" <> ToString[i - 1] <>
  "]+n[" <> ToString[i - 1] <> "]<=" <> ToString[sum0] <> ") && (s1[" <> 
  ToString[i - 1] <> "]+" <> ToString[i] <> "*n[" <> ToString[i - 1] <> "]<=" <> 
  ToString[sum1] <> "); n[" <> ToString[i - 1] <> "]++) {\ns0[" <> ToString[i] <> 
  "]=s0[" <> ToString[i - 1] <> "]+n[" <> ToString[i - 1] <> "];\ns1[" <> 
  ToString[i] <> "]=s1[" <> ToString[i - 1] <> "]+" <> ToString[i] <> "*n[" <> 
  ToString[i - 1] <> "];\n"], {i, imax}];
WriteString[f, "if ((s0[" <> ToString[imax] <> "]==" <> ToString[sum0] <>
  ") && (s1[" <> ToString[imax] <> "]==" <> ToString[sum1] <> "))\n"];
WriteString[f, "printf(\""];
Do[WriteString[f, "%d "], {imax}];
WriteString[f, "\\n\", n[0]"];
Do[WriteString[f, ",n[" <> ToString[i - 1] <> "]"], {i, 2, imax}];
WriteString[f, ");\n"];
Do[WriteString[f, "}"], {imax + 1}];
WriteString[f, "\n"];
Close[f];

Compile this code in a terminal with

gcc -O3 loop.c -o loop

and run it with

./loop > loop.dat

On my laptop the code runs in about 40 seconds and generates 2'977'866 solutions. These solutions can be read into Mathematica with

data = ReadList["~/Desktop/loop.dat", Number, RecordLists -> True];
Dimensions[data]

{2977866, 100}

To determine how many solutions you get as a function of $i_{\text{max}}$, first count how many nonzero $n_i$ are involved in each solution:

lengths = Replace[data, {x___, Except[0], 0 ...} :> Length[{x}] + 1, {1}];

Then count how many solutions have a length $\le i_{\text{max}}$ as a function of $i_{\text{max}}$:

Transpose[{Range[100], Accumulate[Lookup[Counts[lengths], Range[100], 0]]}]

{{1, 0}, {2, 0}, {3, 0}, {4, 0}, {5, 0}, {6, 0}, {7, 0}, {8, 0}, {9, 0}, {10, 1}, {11, 42}, {12, 463}, {13, 2507}, {14, 8861}, {15, 23601}, {16, 51376}, {17, 96314}, {18, 161073}, {19, 246448}, {20, 351344}, {21, 473259}, {22, 608704}, {23, 753813}, {24, 904675}, {25, 1057740}, {26, 1209868}, {27, 1358546}, {28, 1501764}, {29, 1638097}, {30, 1766535}, {31, 1886522}, {32, 1997755}, {33, 2100247}, {34, 2194143}, {35, 2279771}, {36, 2357509}, {37, 2427846}, {38, 2491252}, {39, 2548263}, {40, 2599369}, {41, 2645085}, {42, 2685871}, {43, 2722199}, {44, 2754470}, {45, 2783097}, {46, 2808428}, {47, 2830808}, {48, 2850528}, {49, 2867882}, {50, 2883106}, {51, 2896444}, {52, 2908092}, {53, 2918248}, {54, 2927072}, {55, 2934729}, {56, 2941344}, {57, 2947052}, {58, 2951956}, {59, 2956162}, {60, 2959751}, {61, 2962811}, {62, 2965403}, {63, 2967597}, {64, 2969442}, {65, 2970991}, {66, 2972282}, {67, 2973358}, {68, 2974245}, {69, 2974977}, {70, 2975575}, {71, 2976063}, {72, 2976456}, {73, 2976774}, {74, 2977026}, {75, 2977227}, {76, 2977384}, {77, 2977507}, {78, 2977601}, {79, 2977674}, {80, 2977728}, {81, 2977769}, {82, 2977799}, {83, 2977821}, {84, 2977836}, {85, 2977847}, {86, 2977854}, {87, 2977859}, {88, 2977862}, {89, 2977864}, {90, 2977865}, {91, 2977866}, {92, 2977866}, {93, 2977866}, {94, 2977866}, {95, 2977866}, {96, 2977866}, {97, 2977866}, {98, 2977866}, {99, 2977866}, {100, 2977866}}


For reference, here's the generated C-code:

#include <stdio.h>
int main() {
int n[100];
int s0[101];
int s1[101];
s0[0]=s1[0]=0;
for (n[0]=0; (s0[0]+n[0]<=10) && (s1[0]+1*n[0]<=100); n[0]++) {
s0[1]=s0[0]+n[0];
s1[1]=s1[0]+1*n[0];
for (n[1]=0; (s0[1]+n[1]<=10) && (s1[1]+2*n[1]<=100); n[1]++) {
s0[2]=s0[1]+n[1];
s1[2]=s1[1]+2*n[1];
for (n[2]=0; (s0[2]+n[2]<=10) && (s1[2]+3*n[2]<=100); n[2]++) {
s0[3]=s0[2]+n[2];
s1[3]=s1[2]+3*n[2];
for (n[3]=0; (s0[3]+n[3]<=10) && (s1[3]+4*n[3]<=100); n[3]++) {
s0[4]=s0[3]+n[3];
s1[4]=s1[3]+4*n[3];
for (n[4]=0; (s0[4]+n[4]<=10) && (s1[4]+5*n[4]<=100); n[4]++) {
s0[5]=s0[4]+n[4];
s1[5]=s1[4]+5*n[4];
for (n[5]=0; (s0[5]+n[5]<=10) && (s1[5]+6*n[5]<=100); n[5]++) {
s0[6]=s0[5]+n[5];
s1[6]=s1[5]+6*n[5];
for (n[6]=0; (s0[6]+n[6]<=10) && (s1[6]+7*n[6]<=100); n[6]++) {
s0[7]=s0[6]+n[6];
s1[7]=s1[6]+7*n[6];
for (n[7]=0; (s0[7]+n[7]<=10) && (s1[7]+8*n[7]<=100); n[7]++) {
s0[8]=s0[7]+n[7];
s1[8]=s1[7]+8*n[7];
for (n[8]=0; (s0[8]+n[8]<=10) && (s1[8]+9*n[8]<=100); n[8]++) {
s0[9]=s0[8]+n[8];
s1[9]=s1[8]+9*n[8];
for (n[9]=0; (s0[9]+n[9]<=10) && (s1[9]+10*n[9]<=100); n[9]++) {
s0[10]=s0[9]+n[9];
s1[10]=s1[9]+10*n[9];
for (n[10]=0; (s0[10]+n[10]<=10) && (s1[10]+11*n[10]<=100); n[10]++) {
s0[11]=s0[10]+n[10];
s1[11]=s1[10]+11*n[10];
for (n[11]=0; (s0[11]+n[11]<=10) && (s1[11]+12*n[11]<=100); n[11]++) {
s0[12]=s0[11]+n[11];
s1[12]=s1[11]+12*n[11];
for (n[12]=0; (s0[12]+n[12]<=10) && (s1[12]+13*n[12]<=100); n[12]++) {
s0[13]=s0[12]+n[12];
s1[13]=s1[12]+13*n[12];
for (n[13]=0; (s0[13]+n[13]<=10) && (s1[13]+14*n[13]<=100); n[13]++) {
s0[14]=s0[13]+n[13];
s1[14]=s1[13]+14*n[13];
for (n[14]=0; (s0[14]+n[14]<=10) && (s1[14]+15*n[14]<=100); n[14]++) {
s0[15]=s0[14]+n[14];
s1[15]=s1[14]+15*n[14];
for (n[15]=0; (s0[15]+n[15]<=10) && (s1[15]+16*n[15]<=100); n[15]++) {
s0[16]=s0[15]+n[15];
s1[16]=s1[15]+16*n[15];
for (n[16]=0; (s0[16]+n[16]<=10) && (s1[16]+17*n[16]<=100); n[16]++) {
s0[17]=s0[16]+n[16];
s1[17]=s1[16]+17*n[16];
for (n[17]=0; (s0[17]+n[17]<=10) && (s1[17]+18*n[17]<=100); n[17]++) {
s0[18]=s0[17]+n[17];
s1[18]=s1[17]+18*n[17];
for (n[18]=0; (s0[18]+n[18]<=10) && (s1[18]+19*n[18]<=100); n[18]++) {
s0[19]=s0[18]+n[18];
s1[19]=s1[18]+19*n[18];
for (n[19]=0; (s0[19]+n[19]<=10) && (s1[19]+20*n[19]<=100); n[19]++) {
s0[20]=s0[19]+n[19];
s1[20]=s1[19]+20*n[19];
for (n[20]=0; (s0[20]+n[20]<=10) && (s1[20]+21*n[20]<=100); n[20]++) {
s0[21]=s0[20]+n[20];
s1[21]=s1[20]+21*n[20];
for (n[21]=0; (s0[21]+n[21]<=10) && (s1[21]+22*n[21]<=100); n[21]++) {
s0[22]=s0[21]+n[21];
s1[22]=s1[21]+22*n[21];
for (n[22]=0; (s0[22]+n[22]<=10) && (s1[22]+23*n[22]<=100); n[22]++) {
s0[23]=s0[22]+n[22];
s1[23]=s1[22]+23*n[22];
for (n[23]=0; (s0[23]+n[23]<=10) && (s1[23]+24*n[23]<=100); n[23]++) {
s0[24]=s0[23]+n[23];
s1[24]=s1[23]+24*n[23];
for (n[24]=0; (s0[24]+n[24]<=10) && (s1[24]+25*n[24]<=100); n[24]++) {
s0[25]=s0[24]+n[24];
s1[25]=s1[24]+25*n[24];
for (n[25]=0; (s0[25]+n[25]<=10) && (s1[25]+26*n[25]<=100); n[25]++) {
s0[26]=s0[25]+n[25];
s1[26]=s1[25]+26*n[25];
for (n[26]=0; (s0[26]+n[26]<=10) && (s1[26]+27*n[26]<=100); n[26]++) {
s0[27]=s0[26]+n[26];
s1[27]=s1[26]+27*n[26];
for (n[27]=0; (s0[27]+n[27]<=10) && (s1[27]+28*n[27]<=100); n[27]++) {
s0[28]=s0[27]+n[27];
s1[28]=s1[27]+28*n[27];
for (n[28]=0; (s0[28]+n[28]<=10) && (s1[28]+29*n[28]<=100); n[28]++) {
s0[29]=s0[28]+n[28];
s1[29]=s1[28]+29*n[28];
for (n[29]=0; (s0[29]+n[29]<=10) && (s1[29]+30*n[29]<=100); n[29]++) {
s0[30]=s0[29]+n[29];
s1[30]=s1[29]+30*n[29];
for (n[30]=0; (s0[30]+n[30]<=10) && (s1[30]+31*n[30]<=100); n[30]++) {
s0[31]=s0[30]+n[30];
s1[31]=s1[30]+31*n[30];
for (n[31]=0; (s0[31]+n[31]<=10) && (s1[31]+32*n[31]<=100); n[31]++) {
s0[32]=s0[31]+n[31];
s1[32]=s1[31]+32*n[31];
for (n[32]=0; (s0[32]+n[32]<=10) && (s1[32]+33*n[32]<=100); n[32]++) {
s0[33]=s0[32]+n[32];
s1[33]=s1[32]+33*n[32];
for (n[33]=0; (s0[33]+n[33]<=10) && (s1[33]+34*n[33]<=100); n[33]++) {
s0[34]=s0[33]+n[33];
s1[34]=s1[33]+34*n[33];
for (n[34]=0; (s0[34]+n[34]<=10) && (s1[34]+35*n[34]<=100); n[34]++) {
s0[35]=s0[34]+n[34];
s1[35]=s1[34]+35*n[34];
for (n[35]=0; (s0[35]+n[35]<=10) && (s1[35]+36*n[35]<=100); n[35]++) {
s0[36]=s0[35]+n[35];
s1[36]=s1[35]+36*n[35];
for (n[36]=0; (s0[36]+n[36]<=10) && (s1[36]+37*n[36]<=100); n[36]++) {
s0[37]=s0[36]+n[36];
s1[37]=s1[36]+37*n[36];
for (n[37]=0; (s0[37]+n[37]<=10) && (s1[37]+38*n[37]<=100); n[37]++) {
s0[38]=s0[37]+n[37];
s1[38]=s1[37]+38*n[37];
for (n[38]=0; (s0[38]+n[38]<=10) && (s1[38]+39*n[38]<=100); n[38]++) {
s0[39]=s0[38]+n[38];
s1[39]=s1[38]+39*n[38];
for (n[39]=0; (s0[39]+n[39]<=10) && (s1[39]+40*n[39]<=100); n[39]++) {
s0[40]=s0[39]+n[39];
s1[40]=s1[39]+40*n[39];
for (n[40]=0; (s0[40]+n[40]<=10) && (s1[40]+41*n[40]<=100); n[40]++) {
s0[41]=s0[40]+n[40];
s1[41]=s1[40]+41*n[40];
for (n[41]=0; (s0[41]+n[41]<=10) && (s1[41]+42*n[41]<=100); n[41]++) {
s0[42]=s0[41]+n[41];
s1[42]=s1[41]+42*n[41];
for (n[42]=0; (s0[42]+n[42]<=10) && (s1[42]+43*n[42]<=100); n[42]++) {
s0[43]=s0[42]+n[42];
s1[43]=s1[42]+43*n[42];
for (n[43]=0; (s0[43]+n[43]<=10) && (s1[43]+44*n[43]<=100); n[43]++) {
s0[44]=s0[43]+n[43];
s1[44]=s1[43]+44*n[43];
for (n[44]=0; (s0[44]+n[44]<=10) && (s1[44]+45*n[44]<=100); n[44]++) {
s0[45]=s0[44]+n[44];
s1[45]=s1[44]+45*n[44];
for (n[45]=0; (s0[45]+n[45]<=10) && (s1[45]+46*n[45]<=100); n[45]++) {
s0[46]=s0[45]+n[45];
s1[46]=s1[45]+46*n[45];
for (n[46]=0; (s0[46]+n[46]<=10) && (s1[46]+47*n[46]<=100); n[46]++) {
s0[47]=s0[46]+n[46];
s1[47]=s1[46]+47*n[46];
for (n[47]=0; (s0[47]+n[47]<=10) && (s1[47]+48*n[47]<=100); n[47]++) {
s0[48]=s0[47]+n[47];
s1[48]=s1[47]+48*n[47];
for (n[48]=0; (s0[48]+n[48]<=10) && (s1[48]+49*n[48]<=100); n[48]++) {
s0[49]=s0[48]+n[48];
s1[49]=s1[48]+49*n[48];
for (n[49]=0; (s0[49]+n[49]<=10) && (s1[49]+50*n[49]<=100); n[49]++) {
s0[50]=s0[49]+n[49];
s1[50]=s1[49]+50*n[49];
for (n[50]=0; (s0[50]+n[50]<=10) && (s1[50]+51*n[50]<=100); n[50]++) {
s0[51]=s0[50]+n[50];
s1[51]=s1[50]+51*n[50];
for (n[51]=0; (s0[51]+n[51]<=10) && (s1[51]+52*n[51]<=100); n[51]++) {
s0[52]=s0[51]+n[51];
s1[52]=s1[51]+52*n[51];
for (n[52]=0; (s0[52]+n[52]<=10) && (s1[52]+53*n[52]<=100); n[52]++) {
s0[53]=s0[52]+n[52];
s1[53]=s1[52]+53*n[52];
for (n[53]=0; (s0[53]+n[53]<=10) && (s1[53]+54*n[53]<=100); n[53]++) {
s0[54]=s0[53]+n[53];
s1[54]=s1[53]+54*n[53];
for (n[54]=0; (s0[54]+n[54]<=10) && (s1[54]+55*n[54]<=100); n[54]++) {
s0[55]=s0[54]+n[54];
s1[55]=s1[54]+55*n[54];
for (n[55]=0; (s0[55]+n[55]<=10) && (s1[55]+56*n[55]<=100); n[55]++) {
s0[56]=s0[55]+n[55];
s1[56]=s1[55]+56*n[55];
for (n[56]=0; (s0[56]+n[56]<=10) && (s1[56]+57*n[56]<=100); n[56]++) {
s0[57]=s0[56]+n[56];
s1[57]=s1[56]+57*n[56];
for (n[57]=0; (s0[57]+n[57]<=10) && (s1[57]+58*n[57]<=100); n[57]++) {
s0[58]=s0[57]+n[57];
s1[58]=s1[57]+58*n[57];
for (n[58]=0; (s0[58]+n[58]<=10) && (s1[58]+59*n[58]<=100); n[58]++) {
s0[59]=s0[58]+n[58];
s1[59]=s1[58]+59*n[58];
for (n[59]=0; (s0[59]+n[59]<=10) && (s1[59]+60*n[59]<=100); n[59]++) {
s0[60]=s0[59]+n[59];
s1[60]=s1[59]+60*n[59];
for (n[60]=0; (s0[60]+n[60]<=10) && (s1[60]+61*n[60]<=100); n[60]++) {
s0[61]=s0[60]+n[60];
s1[61]=s1[60]+61*n[60];
for (n[61]=0; (s0[61]+n[61]<=10) && (s1[61]+62*n[61]<=100); n[61]++) {
s0[62]=s0[61]+n[61];
s1[62]=s1[61]+62*n[61];
for (n[62]=0; (s0[62]+n[62]<=10) && (s1[62]+63*n[62]<=100); n[62]++) {
s0[63]=s0[62]+n[62];
s1[63]=s1[62]+63*n[62];
for (n[63]=0; (s0[63]+n[63]<=10) && (s1[63]+64*n[63]<=100); n[63]++) {
s0[64]=s0[63]+n[63];
s1[64]=s1[63]+64*n[63];
for (n[64]=0; (s0[64]+n[64]<=10) && (s1[64]+65*n[64]<=100); n[64]++) {
s0[65]=s0[64]+n[64];
s1[65]=s1[64]+65*n[64];
for (n[65]=0; (s0[65]+n[65]<=10) && (s1[65]+66*n[65]<=100); n[65]++) {
s0[66]=s0[65]+n[65];
s1[66]=s1[65]+66*n[65];
for (n[66]=0; (s0[66]+n[66]<=10) && (s1[66]+67*n[66]<=100); n[66]++) {
s0[67]=s0[66]+n[66];
s1[67]=s1[66]+67*n[66];
for (n[67]=0; (s0[67]+n[67]<=10) && (s1[67]+68*n[67]<=100); n[67]++) {
s0[68]=s0[67]+n[67];
s1[68]=s1[67]+68*n[67];
for (n[68]=0; (s0[68]+n[68]<=10) && (s1[68]+69*n[68]<=100); n[68]++) {
s0[69]=s0[68]+n[68];
s1[69]=s1[68]+69*n[68];
for (n[69]=0; (s0[69]+n[69]<=10) && (s1[69]+70*n[69]<=100); n[69]++) {
s0[70]=s0[69]+n[69];
s1[70]=s1[69]+70*n[69];
for (n[70]=0; (s0[70]+n[70]<=10) && (s1[70]+71*n[70]<=100); n[70]++) {
s0[71]=s0[70]+n[70];
s1[71]=s1[70]+71*n[70];
for (n[71]=0; (s0[71]+n[71]<=10) && (s1[71]+72*n[71]<=100); n[71]++) {
s0[72]=s0[71]+n[71];
s1[72]=s1[71]+72*n[71];
for (n[72]=0; (s0[72]+n[72]<=10) && (s1[72]+73*n[72]<=100); n[72]++) {
s0[73]=s0[72]+n[72];
s1[73]=s1[72]+73*n[72];
for (n[73]=0; (s0[73]+n[73]<=10) && (s1[73]+74*n[73]<=100); n[73]++) {
s0[74]=s0[73]+n[73];
s1[74]=s1[73]+74*n[73];
for (n[74]=0; (s0[74]+n[74]<=10) && (s1[74]+75*n[74]<=100); n[74]++) {
s0[75]=s0[74]+n[74];
s1[75]=s1[74]+75*n[74];
for (n[75]=0; (s0[75]+n[75]<=10) && (s1[75]+76*n[75]<=100); n[75]++) {
s0[76]=s0[75]+n[75];
s1[76]=s1[75]+76*n[75];
for (n[76]=0; (s0[76]+n[76]<=10) && (s1[76]+77*n[76]<=100); n[76]++) {
s0[77]=s0[76]+n[76];
s1[77]=s1[76]+77*n[76];
for (n[77]=0; (s0[77]+n[77]<=10) && (s1[77]+78*n[77]<=100); n[77]++) {
s0[78]=s0[77]+n[77];
s1[78]=s1[77]+78*n[77];
for (n[78]=0; (s0[78]+n[78]<=10) && (s1[78]+79*n[78]<=100); n[78]++) {
s0[79]=s0[78]+n[78];
s1[79]=s1[78]+79*n[78];
for (n[79]=0; (s0[79]+n[79]<=10) && (s1[79]+80*n[79]<=100); n[79]++) {
s0[80]=s0[79]+n[79];
s1[80]=s1[79]+80*n[79];
for (n[80]=0; (s0[80]+n[80]<=10) && (s1[80]+81*n[80]<=100); n[80]++) {
s0[81]=s0[80]+n[80];
s1[81]=s1[80]+81*n[80];
for (n[81]=0; (s0[81]+n[81]<=10) && (s1[81]+82*n[81]<=100); n[81]++) {
s0[82]=s0[81]+n[81];
s1[82]=s1[81]+82*n[81];
for (n[82]=0; (s0[82]+n[82]<=10) && (s1[82]+83*n[82]<=100); n[82]++) {
s0[83]=s0[82]+n[82];
s1[83]=s1[82]+83*n[82];
for (n[83]=0; (s0[83]+n[83]<=10) && (s1[83]+84*n[83]<=100); n[83]++) {
s0[84]=s0[83]+n[83];
s1[84]=s1[83]+84*n[83];
for (n[84]=0; (s0[84]+n[84]<=10) && (s1[84]+85*n[84]<=100); n[84]++) {
s0[85]=s0[84]+n[84];
s1[85]=s1[84]+85*n[84];
for (n[85]=0; (s0[85]+n[85]<=10) && (s1[85]+86*n[85]<=100); n[85]++) {
s0[86]=s0[85]+n[85];
s1[86]=s1[85]+86*n[85];
for (n[86]=0; (s0[86]+n[86]<=10) && (s1[86]+87*n[86]<=100); n[86]++) {
s0[87]=s0[86]+n[86];
s1[87]=s1[86]+87*n[86];
for (n[87]=0; (s0[87]+n[87]<=10) && (s1[87]+88*n[87]<=100); n[87]++) {
s0[88]=s0[87]+n[87];
s1[88]=s1[87]+88*n[87];
for (n[88]=0; (s0[88]+n[88]<=10) && (s1[88]+89*n[88]<=100); n[88]++) {
s0[89]=s0[88]+n[88];
s1[89]=s1[88]+89*n[88];
for (n[89]=0; (s0[89]+n[89]<=10) && (s1[89]+90*n[89]<=100); n[89]++) {
s0[90]=s0[89]+n[89];
s1[90]=s1[89]+90*n[89];
for (n[90]=0; (s0[90]+n[90]<=10) && (s1[90]+91*n[90]<=100); n[90]++) {
s0[91]=s0[90]+n[90];
s1[91]=s1[90]+91*n[90];
for (n[91]=0; (s0[91]+n[91]<=10) && (s1[91]+92*n[91]<=100); n[91]++) {
s0[92]=s0[91]+n[91];
s1[92]=s1[91]+92*n[91];
for (n[92]=0; (s0[92]+n[92]<=10) && (s1[92]+93*n[92]<=100); n[92]++) {
s0[93]=s0[92]+n[92];
s1[93]=s1[92]+93*n[92];
for (n[93]=0; (s0[93]+n[93]<=10) && (s1[93]+94*n[93]<=100); n[93]++) {
s0[94]=s0[93]+n[93];
s1[94]=s1[93]+94*n[93];
for (n[94]=0; (s0[94]+n[94]<=10) && (s1[94]+95*n[94]<=100); n[94]++) {
s0[95]=s0[94]+n[94];
s1[95]=s1[94]+95*n[94];
for (n[95]=0; (s0[95]+n[95]<=10) && (s1[95]+96*n[95]<=100); n[95]++) {
s0[96]=s0[95]+n[95];
s1[96]=s1[95]+96*n[95];
for (n[96]=0; (s0[96]+n[96]<=10) && (s1[96]+97*n[96]<=100); n[96]++) {
s0[97]=s0[96]+n[96];
s1[97]=s1[96]+97*n[96];
for (n[97]=0; (s0[97]+n[97]<=10) && (s1[97]+98*n[97]<=100); n[97]++) {
s0[98]=s0[97]+n[97];
s1[98]=s1[97]+98*n[97];
for (n[98]=0; (s0[98]+n[98]<=10) && (s1[98]+99*n[98]<=100); n[98]++) {
s0[99]=s0[98]+n[98];
s1[99]=s1[98]+99*n[98];
for (n[99]=0; (s0[99]+n[99]<=10) && (s1[99]+100*n[99]<=100); n[99]++) {
s0[100]=s0[99]+n[99];
s1[100]=s1[99]+100*n[99];
if ((s0[100]==10) && (s1[100]==100))
printf("%d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d \n", n[0],n[1],n[2],n[3],n[4],n[5],n[6],n[7],n[8],n[9],n[10],n[11],n[12],n[13],n[14],n[15],n[16],n[17],n[18],n[19],n[20],n[21],n[22],n[23],n[24],n[25],n[26],n[27],n[28],n[29],n[30],n[31],n[32],n[33],n[34],n[35],n[36],n[37],n[38],n[39],n[40],n[41],n[42],n[43],n[44],n[45],n[46],n[47],n[48],n[49],n[50],n[51],n[52],n[53],n[54],n[55],n[56],n[57],n[58],n[59],n[60],n[61],n[62],n[63],n[64],n[65],n[66],n[67],n[68],n[69],n[70],n[71],n[72],n[73],n[74],n[75],n[76],n[77],n[78],n[79],n[80],n[81],n[82],n[83],n[84],n[85],n[86],n[87],n[88],n[89],n[90],n[91],n[92],n[93],n[94],n[95],n[96],n[97],n[98],n[99]);
}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}}

Surely someone more skilled than me can find a C++ template programming way of generating the same code without needing Mathematica in the first place.

$\endgroup$
3
$\begingroup$

All the $n_i$ with $i>100$ must be zero, so it's enough to look at $n_1\ldots n_{100}$. The below code will crash for this large system though.

However, if you only look at $i\le i_{\text{max}}$ (and you may make $i_{\text{max}}$ as large as your computer allows, though probably not 100 as truly needed), then you could first generate the list of all tuples of length imax that sum to 10:

imax = 15;
S = Join @@ Permutations /@ IntegerPartitions[10, {imax}, Range[0, 10]];

Careful with memory, there are $\binom{i_{\text{max}} + 9}{10}$ of them. For $i_{\text{max}}=100$ this will generate 42'634'215'112'710 tuples of length 100.

Then find the positions of those that sum to 100 when dotted with the vector $\{1,2,3,\ldots,i_{\text{max}}\}$:

P = Position[S.Range[imax], 100] // Flatten

{10, 39, 276, 359, 398, 573, ...}

Show the results:

S[[P]]

{{0, 0, 0, 0, 0, 0, 0, 0, 0, 10, 0, 0, 0, 0, 0}, {1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 9, 0, 0, 0, 0}, {0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 8, 0, 0, 0}, {0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 8, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 0, 0, 0, 8, 0, 0, 0, 0, 2, 0}, {1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 8, 0, 0, 0}, ...}

To calculate the number of solutions as a function of $i_{\text{max}}$:

num[imax_Integer /; 1 <= imax <= 100] := Module[{S},
  S = Join @@ Permutations /@ IntegerPartitions[10, {imax}, Range[0, 10]];
  P = Position[S.Range[imax], 100] // Flatten;
  Length[P]]

Table[num[imax], {imax, 1, 20}]

{0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 42, 463, 2507, 8861, 23601, 51376, 96314, 161073, 246448, 351344}

$\endgroup$
  • $\begingroup$ Thanks. Can you add a code fragment that computes the number of such tuples as a function of imax $\endgroup$ – Quasar Supernova Feb 22 at 16:53
2
$\begingroup$

One approach is to name your variables as elements of an array (instead of having separate names). For example, allA defines 100 variables called a[1], a[2], etc. You can use these in Solve directly, and force them all to be positive integers:

 n = 100;
 allA = Array[a, n]
 Solve[{Sum[a[i], {i, 1, n}] == 100 && Sum[i a[i], {i, 1, n}] == 1000 &&
        Thread[allA >= ConstantArray[0, n]] //.List -> And}, allA, Integers]

We can test that it works by using a simple case where n=3

n = 3;
allA = Array[a, n]
Solve[{Sum[a[i], {i, 1, n}] == 6 && Sum[i a[i], {i, 1, n}] == 14 && 
    Thread[allA >= ConstantArray[0, n]] //.List -> And}, allA, Integers]

which returns three answers:

{{a[1] -> 0, a[2] -> 4, a[3] -> 2}, {a[1] -> 1, a[2] -> 2, a[3] -> 3}
 {a[1] -> 2, a[2] -> 0, a[3] -> 4}}

For larger n, if you don't want to wait for all the answers, you can replace Solve with FindInstance and just get a few answers.

$\endgroup$
1
$\begingroup$

I used FrobeniusSolve but in the most wasteful way imaginable. I solved these two equations separately and found the intersection of the solutions and counted how many there are as a function of U. The Log of this number is the entropy of the system. The physical model is a set of N identical marbles occupying an infinite staircase where ascending each step implies gaining a unit of energy. The total energy of the system is U.

n1 + n2 + ... + nK = N
n1 + 2 n2 + .... + K nK = U

It is easy to see that

K = U - N + 1

The Code I wrote is,

funit[x_] := 1
Arr[max_] := Array[funit, max]

FSolveTotN[ N_, U_] 
:= FrobeniusSolve[ Arr[U - N + 1], N]

FSolveEgy[N_, U_] 
:= FrobeniusSolve[ Range[U - N + 1], U]

MarbleList[N_, U_] 
:= Intersection[FSolveTotN[N, U],FSolveEgy[N, U]]

Entropy of N marbles on an infinite staircase as a function of the total energy U

S[N_, U_] := Log[ Length[MarbleList[N, U]] ]
TASK[N_] := Table[{U, S[N, U]},{U, N+1, N+15}]

ListPlot[TASK[3], AxesLabel -> {"U", "Entropy"}]

TASK[3] is

{{4, 0}, {5, Log[2]}, {6, Log[3]}, {7, Log[4]}, 
{8, Log[5]}, {9,Log[7]}, {10, Log[8]}, 
{11, Log[10]}, {12, Log[12]}, {13,Log[14]},
{14, Log[16]}, {15, Log[19]}, {16, Log[21]},
{17, Log[24]}, {18, Log[27]}}    
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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