# Code Producing Power::Infty

I tried computing the following

Unprotect[C]
\[Psi] = Sqrt[2] - 1;
\[Epsilon] = .001;
\[Omega] = 50000;
k = 13;
j = 3;
A1 = Flatten[
Table[(s1^j)/(s2^k), {s2, 1, Floor[\[Omega]^(1/k)]}, {s1,Ceiling[((s2^k) (\[Psi] - \[Epsilon]))^(1/j)],
Floor[((s2^k) (\[Psi] + \[Epsilon]))^(1/j)] }]];
A2 = Flatten[
Table[(s1^j)/(s2^k), {s1, 1, Floor[\[Omega]^(1/j)]}, {s2,
Ceiling[((s1^j)/(\[Psi] + \[Epsilon]))^(1/k)],
Floor[((s1^j)/(\[Psi] - \[Epsilon]))^(1/k)]}]];
A = DeleteDuplicates[Flatten[Intersection[A1, A2]]];
C1 = DeleteDuplicates[
Flatten[Table[(s1^j)/(2 s2 + 1)^k, {s2, 1,
Floor[(10 \[Omega])^(1/k)]}, {s1,
Ceiling[(((2 s2 + 1)^k) (\[Psi] - \[Epsilon]))^(1/j)],
Floor[(((2 s2 + 1)^k) (\[Psi] + \[Epsilon]))^(1/j)]}]]];
C2 = DeleteDuplicates[
Flatten[Table[(s1^j)/(2 s2 + 1)^k, {s1, 1,
Floor[(10 \[Omega])^(1/j)]}, {s2,
Ceiling[(((s1^j)/(\[Psi] + \[Epsilon]))^(1/k) - 1)/2],
Floor[(((s1^j)/(\[Psi] - \[Epsilon]))^(1/k) - 1)/2]}]]];
C = DeleteDuplicates[Flatten[Intersection[C1, C2]]];
g1 = N[Length[Intersection[A, C]]/Length[A]];


When j=1 the code works fine but when j is an integer greater than 1, it returns an error sign:

Power::infy: Infinite expression 1/0 encountered.

Infinity::indet: Indeterminate expression 0 ComplexInfinity encountered.


I tried fixing it using DeleteCases to delete the zeros (see t1 and t2)

Unprotect[C]
\[Psi] = Sqrt[2] - 1;
\[Epsilon] = .001;
\[Omega] = 50000;
k = 13;
j = 3;
t1 = DeleteCases[
Table[s1, {s1, Ceiling[((s2^k) (\[Psi] - \[Epsilon]))^(1/j)],
Floor[((s2^k) (\[Psi] + \[Epsilon]))^(1/j)]}], 0]
t2 = DeleteCases[
Table[s2, {s2, Ceiling[((s1^j)/(\[Psi] + \[Epsilon]))^(1/k)],
Floor[((s1^j)/(\[Psi] - \[Epsilon]))^(1/k)]}], 0]
A1 = Flatten[
Table[(s1^j)/(s2^k), {s2, 1, Floor[\[Omega]^(1/k)]}, {s1,t1 }]];
A2 = Flatten[
Table[(s1^j)/(s2^k), {s1, 1, Floor[\[Omega]^(1/j)]}, {s2,
t2}]];
A = DeleteDuplicates[Flatten[Intersection[A1, A2]]];
C1 = DeleteDuplicates[
Flatten[Table[(s1^j)/(2 s2 + 1)^k, {s2, 1,
Floor[(10 \[Omega])^(1/k)]}, {s1,
Ceiling[(((2 s2 + 1)^k) (\[Psi] - \[Epsilon]))^(1/j)],
Floor[(((2 s2 + 1)^k) (\[Psi] + \[Epsilon]))^(1/j)]}]]];
C2 = DeleteDuplicates[
Flatten[Table[(s1^j)/(2 s2 + 1)^k, {s1, 1,
Floor[(10 \[Omega])^(1/j)]}, {s2,
Ceiling[(((s1^j)/(\[Psi] + \[Epsilon]))^(1/k) - 1)/2],
Floor[(((s1^j)/(\[Psi] - \[Epsilon]))^(1/k) - 1)/2]}]]];
C = DeleteDuplicates[Flatten[Intersection[C1, C2]]];
g1 = N[Length[Intersection[A, C]]/Length[A]];


However, I continue to get Power::Infty.

How do we fix my code?

• The way to solve this is to make the code smaller and smaller until you find the exact location of the problem. In other words: please create a minimal example. While this does take some effort, it is the way to solve the problem, and anyone answering your question would have to do it. Therefore, please try to do it on your own first. See here: mathematica.meta.stackexchange.com/q/2126/12 Jun 25, 2021 at 12:54
• Unprotect[C] <- This is a very bad idea and it does break stuff. Even if it is not responsible for the specific issue you are asking about, it will prevent some people from evaluating your code and try to help you. Jun 25, 2021 at 13:06
• @Szabolcs I get that A1 and A2 are empty when j is an integer greater than one. I believe this is since the lower bound of A1 and A2 is greater than the upper bound of A1 and A2. Jun 25, 2021 at 13:11
• Crossposted here. Jun 25, 2021 at 17:45
• What do you want to happen when Length[A]==0? You have N[Length[Intersection[A, C]]/Length[A]] so of course you can expect that error when A is empty. Jun 25, 2021 at 18:19

Your code is working as intended. As others have mentioned in the comments, your code could use some cleaning up. But the basic issue is when $$j>1$$ both A1 and A2 are empty (the bounds on your table are trivial). That means Intersection[A1, A2] is empty so A is empty and Length[A]==0. This gives g1 = N[Length[Intersection[A, C]]/Length[A]] $$=$$ Something$$/0$$ which gives your error. So the question is, what do you want g1 to be when Length[A]==0? Whatever you want to do, you can accomplish it simply by doing
g1 = If[Length[A]>0, N[Length[Intersection[A,C]]/Length[A]],