Jump to content

Multiplicatively perfect numbers

From Rosetta Code
(Redirected from Special numbers)
Multiplicatively perfect numbers is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

If the product of the divisors of an integer n (including n itself) is equal to n^2, then n is a multiplicatively perfect number. Alternatively: the product of the proper divisors of n (i.e. excluding n) is equal to n.

Note that the integer '1' qualifies under the first definition (1 = 1 x 1), but not under the second since '1' has no proper divisors. Given this ambiguity, it is optional whether '1' is included or not in solutions to this task.

Task

Find and show on this page the multiplicatively perfect numbers below 500.

Stretch

Find and show the number of multiplicatively perfect numbers under 500, 5,000, 50,000 and 500,000 and for each of these limits deduce (avoid counting separately) and show the number of semi-primes (numbers which are the product of exactly two primes) under that limit.

Related (and near duplicate) task
See also



ALGOL 68

Library: ALGOL 68-primes

Uses sieves, counts 1 as a multiplicatively perfect number.
As with the Phix and probably other samples, uses the fact that the multiplicatively perfect numbers (other than 1) must have three proper divisors - see OIES A007422. 

BEGIN # find multiplicatively perfect numbers - numbers whose proper         #
      # divisor product is the number itself ( includes 1 )                  #
    PR read "primes.incl.a68" PR          # include prime utilties           #
    INT max number = 500 000;             # largest number we wlil consider  #
    []BOOL prime = PRIMESIEVE max number; # sieve the primes to max number   #
    [ 1 : max number ]INT pdc;            # table of proper divisor counts   #
    [ 1 : max number ]INT pfc;            # table of prime factor counts     #
    # count the proper divisors and prime factors                            #
    FOR n TO UPB pdc DO pdc[ n ] := 1; pfc[ n ] := 0 OD;
    FOR i FROM 2 TO UPB pdc DO
        INT m := 1;                       # j will be the mth multiple of i  #
        FOR j FROM i + i BY i TO UPB pdc DO
            pdc[ j ] +:= 1;
            IF prime[ m +:= 1 ] THEN      # j is a prime multiple of i       #
               pfc[ j ] +:= 1;
               IF i = m THEN              # j is i squared                   #
                  pfc[ j ] +:= 1
               FI
            ELSE                          # j is not a prime multiple of i   #
               pfc[ j ] +:= pfc[ m ]      # add the prime factor count of m  #
            FI
        OD
    OD;
    # find the mutiplicatively perfect numbers                               #
    INT mp count     :=   0;
    INT sp count     :=   0;
    INT next to show := 500;
    FOR n TO UPB pdc DO
        IF n = 1 OR pdc[ n ] = 3 THEN
            # n is 1 or has 3 proper divisors so is multiplicatively perfect #
            # number - see OEIS A007422                                      #
            mp count +:= 1;
            IF n < 500 THEN
                print( ( " ", whole( n, -3 ) ) );
                IF mp count MOD 10 = 0 THEN print( ( newline ) ) FI
            FI
        FI;
        IF pfc[ n ] = 2 THEN
            # 2 prime factors, n is semi-prime                               #
            sp count +:= 1
        FI;
        IF n = next to show THEN
            print( ( "Up to ", whole( next to show, -8 )
                   , " there are ", whole( mp count, -6 )
                   , " multiplicatively perfect numbers and ", whole( sp count, -6 )
                   , " semi-primes", newline
                   )
                 );
            next to show *:= 10
        FI
    OD
END
Output:
   1   6   8  10  14  15  21  22  26  27
  33  34  35  38  39  46  51  55  57  58
  62  65  69  74  77  82  85  86  87  91
  93  94  95 106 111 115 118 119 122 123
 125 129 133 134 141 142 143 145 146 155
 158 159 161 166 177 178 183 185 187 194
 201 202 203 205 206 209 213 214 215 217
 218 219 221 226 235 237 247 249 253 254
 259 262 265 267 274 278 287 291 295 298
 299 301 302 303 305 309 314 319 321 323
 326 327 329 334 335 339 341 343 346 355
 358 362 365 371 377 381 382 386 391 393
 394 395 398 403 407 411 413 415 417 422
 427 437 445 446 447 451 453 454 458 466
 469 471 473 478 481 482 485 489 493 497
Up to      500 there are    150 multiplicatively perfect numbers and    153 semi-primes
Up to     5000 there are   1354 multiplicatively perfect numbers and   1365 semi-primes
Up to    50000 there are  12074 multiplicatively perfect numbers and  12110 semi-primes
Up to   500000 there are 108223 multiplicatively perfect numbers and 108326 semi-primes

BASIC

BASIC256

Translation of: FreeBASIC
arraybase 1
limit = 500
dim Divisors(1)
c = 0

print "Special numbers under "; limit; ":"

for n = 1 to limit
	pro = 1
	for m = 2 to ceil(n / 2)
		if n mod m = 0 then
			pro *= m
			c += 1
			redim Divisors(c)
			Divisors[c] = m
		end if
	next m
	ub = Divisors[?]
	if n = pro and ub > 2 then
		print rjust(string(n), 3); " = "; rjust(string(Divisors[ub-1]), 2); " x "; rjust(string(Divisors[ub]), 3)
	end if
next n
Output:
Same as FreeBASIC entry.

True BASIC

Translation of: FreeBASIC
LET limit = 50
DIM divisors(0)
LET c = 1

PRINT "Special numbers under"; limit; ":"

FOR n = 1 TO limit
    LET pro = 1
    FOR m = 2 TO ceil(n/2)
        IF remainder(n, m) = 0 THEN
           LET pro = pro*m
           MAT REDIM divisors(c)
           LET divisors(c) = m
           LET c = c+1
        END IF
    NEXT m
    LET ub = UBOUND(divisors)
    IF n = pro AND ub > 1 THEN PRINT  USING "### = ## x ###": n, divisors(ub-1), divisors(ub)
NEXT n
END
Output:
Same as FreeBASIC entry.

PureBasic

Translation of: FreeBASIC
Macro Ceil (_x_)
  Round(_x_, #PB_Round_Up)
EndMacro

OpenConsole()
Define.i limit = 500
Dim Divisors.i(1)
Define.i n, pro, m, c = 0, ub

PrintN("Special numbers under" + Str(limit) + ":")

For n = 1 To limit
  pro = 1
  For m = 2 To Ceil(n / 2)
    If Mod(n, m) = 0:
      pro * m
      ReDim Divisors(c) : Divisors(c) = m
      c + 1
    EndIf
  Next m
  ub = ArraySize(Divisors())
  If n = pro And ub > 1:
    PrintN(RSet(Str(n),3) + " = " + RSet(Str(Divisors(ub-1)),2) + " x " + RSet(Str(Divisors(ub)),3))
  EndIf
Next n
Input()
CloseConsole()
Output:
Same as FreeBASIC entry.

Yabasic

Translation of: FreeBASIC
limit = 500
dim Divisors(1)
c = 1

print "Special numbers under", limit, ":"

for n = 1 to limit
    pro = 1
    for m = 2 to ceil(n / 2)
        if mod(n, m) = 0 then
            pro = pro * m
            dim Divisors(c) : Divisors(c) = m
            c = c + 1
        fi
    next m
    ub = arraysize(Divisors(), 1)
    if n = pro and ub > 1 then
        print n using "###", " = ", Divisors(ub-1) using "##", " x ", Divisors(ub) using "###"
    fi
next n
Output:
Same as FreeBASIC entry.

C

This includes '1' as an MPN. Run time around 2.3 seconds. 

#include <stdio.h>
#include <stdbool.h>
#include <locale.h>

bool isPrime(int n) {
    if (n < 2) return false;
    if (n%2 == 0) return n == 2;
    if (n%3 == 0) return n == 3;
    int d = 5;
    while (d*d <= n) {
        if (n%d == 0) return false;
        d += 2;
        if (n%d == 0) return false;
        d += 4;
    }
    return true;
}

void divisors(int n, int *divs, int *length) {
    int i, j, k = 1, c = 0;
    if (n%2) k = 2;
    for (i = 1; i*i <= n; i += k) {
        if (i == 1) continue; // exclude 1 and n
        if (!(n%i)) {
            divs[c++] = i;
            if (c > 2) break; // not eligible if has > 2 divisors 
            j = n / i;
            if (j != i) divs[c++] = j;
        }
    }
    *length = c;
}

int main() {
    int i, d, j, k, t, length, prod;
    int divs[4], count = 0, limit = 500, s = 3, c = 3, squares = 1, cubes = 1;
    printf("Multiplicatively perfect numbers under %d:\n", limit);
    setlocale(LC_NUMERIC, "");
    for (i = 1; ; ++i) {
        if (i != 1) {
            divisors(i, divs, &length);
        } else {
            divs[1] = divs[0] = 1;
            length = 2;
        }
        if (length == 2 && divs[0] * divs[1] == i) {
            ++count;
            if (i < 500) {
                printf("%3d  ", i);
                if (!(count%10)) printf("\n");
            }
        }
        if (i == 499) printf("\n");
        if (i >= limit - 1) {
            for (j = s; j * j < limit; j += 2) if (isPrime(j)) ++squares;
            for (k = c; k * k * k < limit; k +=2 ) if (isPrime(k)) ++cubes;
            t = count + squares - cubes - 1;
            printf("Counts under %'9d: MPNs = %'7d  Semi-primes = %'7d\n", limit, count, t);
            if (limit == 5000000) break;
            s = j;
            c = k;
            limit *= 10;
        }
    }
    return 0;
}
Output:
Multiplicatively perfect numbers under 500:
  1    6    8   10   14   15   21   22   26   27  
 33   34   35   38   39   46   51   55   57   58  
 62   65   69   74   77   82   85   86   87   91  
 93   94   95  106  111  115  118  119  122  123  
125  129  133  134  141  142  143  145  146  155  
158  159  161  166  177  178  183  185  187  194  
201  202  203  205  206  209  213  214  215  217  
218  219  221  226  235  237  247  249  253  254  
259  262  265  267  274  278  287  291  295  298  
299  301  302  303  305  309  314  319  321  323  
326  327  329  334  335  339  341  343  346  355  
358  362  365  371  377  381  382  386  391  393  
394  395  398  403  407  411  413  415  417  422  
427  437  445  446  447  451  453  454  458  466  
469  471  473  478  481  482  485  489  493  497  

Counts under       500: MPNs =     150  Semi-primes =     153
Counts under     5,000: MPNs =   1,354  Semi-primes =   1,365
Counts under    50,000: MPNs =  12,074  Semi-primes =  12,110
Counts under   500,000: MPNs = 108,223  Semi-primes = 108,326
Counts under 5,000,000: MPNs = 978,983  Semi-primes = 979,274

C++ 

#include <iostream>
#include <vector>
#include <numeric>

std::vector<int> divisors( int n ) {
   std::vector<int> divisors ;
   for ( int i = 1 ; i < n + 1 ; i++ ) {
      if ( n % i == 0 ) 
         divisors.push_back( i ) ;
   }
   return divisors ;
}

int main( ) {
   std::vector<int> multi_perfect ;
   for ( int i = 1 ; i < 501 ; i++ ) {
      std::vector<int> divis { divisors( i ) } ;
      if ( std::accumulate( divis.begin( ) , divis.end( ) , 1 , 
               std::multiplies<int>() ) == (i * i ) ) 
         multi_perfect.push_back( i ) ;
   }
   std::cout << '(' ;
   int count = 1 ;
   for ( int i : multi_perfect ) {
      std::cout << i << ' ' ;
      if ( count % 15 == 0 )
         std::cout << '\n' ;
      count++ ;
   }
   std::cout << ")\n" ;
   return 0 ;
}
Output:
(1 6 8 10 14 15 21 22 26 27 33 34 35 38 39 
46 51 55 57 58 62 65 69 74 77 82 85 86 87 91 
93 94 95 106 111 115 118 119 122 123 125 129 133 134 141 
142 143 145 146 155 158 159 161 166 177 178 183 185 187 194 
201 202 203 205 206 209 213 214 215 217 218 219 221 226 235 
237 247 249 253 254 259 262 265 267 274 278 287 291 295 298 
299 301 302 303 305 309 314 319 321 323 326 327 329 334 335 
339 341 343 346 355 358 362 365 371 377 381 382 386 391 393 
394 395 398 403 407 411 413 415 417 422 427 437 445 446 447 
451 453 454 458 466 469 471 473 478 481 482 485 489 493 497 
)

EasyLang 

for n = 1 to 500
   pro = 1
   for m = 2 to n div 2
      if n mod m = 0
         pro *= m
      .
   .
   if n = pro
      write n & " "
   .
.
Output:
1 6 8 10 14 15 21 22 26 27 33 34 35 38 39 46 51 55 57 58 62 65 69 74 77 82 85 86 87 91 93 94 95 106 111 115 118 119 122 123 125 129 133 134 141 142 143 145 146 155 158 159 161 166 177 178 183 185 187 194 201 202 203 205 206 209 213 214 215 217 218 219 221 226 235 237 247 249 253 254 259 262 265 267 274 278 287 291 295 298 299 301 302 303 305 309 314 319 321 323 326 327 329 334 335 339 341 343 346 355 358 362 365 371 377 381 382 386 391 393 394 395 398 403 407 411 413 415 417 422 427 437 445 446 447 451 453 454 458 466 469 471 473 478 481 482 485 489 493 497

FreeBASIC 

#define ceil(x) (-((-x*2.0-0.5) Shr 1))
Dim As Integer limit = 500
Dim As Integer n, pro, Divisors(), m, c = 0, ub
Print "Special numbers under"; limit; ":"

For n = 1 To limit
    pro = 1
    For m = 2 To ceil(n / 2)
        If n Mod m = 0 Then            
            pro *= m
            Redim Preserve Divisors(c) : Divisors(c) = m
            c += 1
        End If
    Next m
    ub = Ubound(Divisors)
    If n = pro And ub > 1 Then
        Print Using "### = ## x ###"; n; Divisors(ub-1); Divisors(ub)
    End If
Next n

Sleep
Output:
Similar to Ring entry.

Go

Translation of: C
Library: Go-rcu

Run time around 9.2 seconds. 

package main

import (
    "fmt"
    "rcu"
)

// library method customized for this task
func Divisors(n int) []int {
    var divisors []int
    i := 1
    k := 1
    if n%2 == 1 {
        k = 2
    }
    for ; i*i <= n; i += k {
        if i > 1 && n%i == 0 { // exclude 1 and n
            divisors = append(divisors, i)
            if len(divisors) > 2 { // not eligible if has > 2 divisors
                break
            }
            j := n / i
            if j != i {
                divisors = append(divisors, j)
            }
        }
    }
    return divisors
}

func main() {
    count := 0
    limit := 500
    s := 3
    c := 3
    squares := 1
    cubes := 1
    fmt.Printf("Multiplicatively perfect numbers under %d:\n", limit)
    var divs []int
    for i := 0; ; i++ {
        if i != 1 {
            divs = Divisors(i)
        } else {
            divs = []int{1, 1}
        }
        if len(divs) == 2 && divs[0]*divs[1] == i {
            count++
            if i < 500 {
                fmt.Printf("%3d  ", i)
                if count%10 == 0 {
                    fmt.Println()
                }
            }
        }
        if i == 499 {
            fmt.Println()
        }
        if i >= limit-1 {
            var j, k int
            for j = s; j*j < limit; j += 2 {
                if rcu.IsPrime(j) {
                    squares++
                }
            }
            for k = c; k*k*k < limit; k += 2 {
                if rcu.IsPrime(k) {
                    cubes++
                }
            }
            t := count + squares - cubes - 1
            slimit := rcu.Commatize(limit)
            scount := rcu.Commatize(count)
            st := rcu.Commatize(t)
            fmt.Printf("Counts under %9s: MPNs = %7s  Semi-primes = %7s\n", slimit, scount, st)
            if limit == 5000000 {
                break
            }
            s, c = j, k
            limit *= 10
        }
    }
}
Output:
Same as C example.

Haskell 

divisors :: Int -> [Int]
divisors n = [d | d <- [1..n] , mod n d == 0 ] 

isMultiplicativelyPerfect :: Int -> Bool
isMultiplicativelyPerfect n = (product $ divisors n) == n ^ 2

solution :: [Int]
solution = filter isMultiplicativelyPerfect [1..500]
Output:
[1,6,8,10,14,15,21,22,26,27,33,34,35,38,39,46,51,55,57,58,62,65,69,74,77,82,85,86,87,91,93,94,95,106,111,115,118,119,122,123,125,129,133,134,141,142,143,145,146,155,158,159,161,166,177,178,183,185,187,194,201,202,203,205,206,209,213,214,215,217,218,219,221,226,235,237,247,249,253,254,259,262,265,267,274,278,287,291,295,298,299,301,302,303,305,309,314,319,321,323,326,327,329,334,335,339,341,343,346,355,358,362,365,371,377,381,382,386,391,393,394,395,398,403,407,411,413,415,417,422,427,437,445,446,447,451,453,454,458,466,469,471,473,478,481,482,485,489,493,497]

J

Implementation: 

factors=: [: /:~@, */&>@{@((^ i.@>:)&.>/)@q:~&__
isMPerfect=: *: = */@(factors ::_:)"0

Task example: 

   #(#~ isMPerfect) i.500
150
   10 15$(#~ isMPerfect) i.500
  1   6   8  10  14  15  21  22  26  27  33  34  35  38  39
 46  51  55  57  58  62  65  69  74  77  82  85  86  87  91
 93  94  95 106 111 115 118 119 122 123 125 129 133 134 141
142 143 145 146 155 158 159 161 166 177 178 183 185 187 194
201 202 203 205 206 209 213 214 215 217 218 219 221 226 235
237 247 249 253 254 259 262 265 267 274 278 287 291 295 298
299 301 302 303 305 309 314 319 321 323 326 327 329 334 335
339 341 343 346 355 358 362 365 371 377 381 382 386 391 393
394 395 398 403 407 411 413 415 417 422 427 437 445 446 447
451 453 454 458 466 469 471 473 478 481 482 485 489 493 497

For the stretch goal, we need to determine the number of semi-primes, given the number of multiplicatively perfect numbers less than N: 

adjSemiPrime=: + _1 + %: -&(p:inv) 3&%:

Thus (first number in following results is count of multiplicatively perfect numbers, second is count of semiprimes): 

   {{ (, adjSemiPrime&y) +/isMPerfect i.y}} 500
150 153
   {{ (, adjSemiPrime&y) +/isMPerfect i.y}} 5000
1354 1365
   {{ (, adjSemiPrime&y) +/isMPerfect i.y}} 50000
12074 12110
   {{ (, adjSemiPrime&y) +/isMPerfect i.y}} 500000
108223 108326

Java 

public final class MultiplicativelyPerfectNumbers {

	public static void main(String[] args) {
		int count = 0;
		for ( int i = 1; i < 500; i++ ) {
			if ( isMPN(i) ) {
				count += 1;
				System.out.print(String.format("%3d%s", i, ( count % 10 == 0 ? "\n" : " " )));
			}
		}
		System.out.println();		
		
		int mpnCount = 0; 
		int limit = 500;
		int ns = 3;
		int nc = 3;
		int squares = 1;
		int cubes = 1;
		int n = 0;
		while ( true ) {
			n += 1;
			if ( n == limit ) {
				while ( ns * ns < limit ) {
					if ( isPrime(ns) ) {
						squares += 1;
					}
					ns += 2;
				}
				while ( nc * nc * nc < limit ) {
					if ( isPrime(nc) ) { 
						cubes += 1;
					}
					nc += 2;
				}
			System.out.println("Under " + limit + " there are " + mpnCount + " MPNs and "
				+ ( mpnCount - cubes + squares ) + " semi-primes.");
		    if ( limit == 500_000 ) {
		    	break;
		    }
		    limit *= 10;
		}
		if ( isMPN(n) ) {
			mpnCount += 1;
		}
	}

}
	
	private static boolean isMPN(int number) {
		if ( number == 1 ) { // Counting 1 as an MPN in accordance with www.oeis.net/A007422
			return true; 
		} 
		
		final int delta = 1 + ( number & 1 );		
		int d = delta + 1;
		int firstDivisor = 0;
		int secondDivisor = 0;
		while ( d * d <= number ) {
			if ( number % d == 0 ) {
				if ( secondDivisor != 0 ) { // From www.oeis.net/A007422, an MPN
					return false;           // cannot have more than two proper divisors 
				} 
				firstDivisor = d;                           
				final int q = number / d;
				if ( q != d ) { 
					secondDivisor = q;
				}
			}
			d += delta;
		}
		return firstDivisor * secondDivisor == number;
	}

	private static boolean isPrime(int n) {
		if ( n < 2 ) { 
			return false;
		}
		if ( n % 2 == 0 ) { 
			return n == 2;
		}
		if ( n % 3 == 0 ) {
			return n == 3;
		}
		
		int k = 5;
		int delta = 2;
		while ( k * k <= n ) {
			if ( n % k == 0 ) {
				return false;
			}
			k += delta;
			delta = 6 - delta;
		}
		return true;
	}

}
Output:
  1   6   8  10  14  15  21  22  26  27
 33  34  35  38  39  46  51  55  57  58
 62  65  69  74  77  82  85  86  87  91
 93  94  95 106 111 115 118 119 122 123
125 129 133 134 141 142 143 145 146 155
158 159 161 166 177 178 183 185 187 194
201 202 203 205 206 209 213 214 215 217
218 219 221 226 235 237 247 249 253 254
259 262 265 267 274 278 287 291 295 298
299 301 302 303 305 309 314 319 321 323
326 327 329 334 335 339 341 343 346 355
358 362 365 371 377 381 382 386 391 393
394 395 398 403 407 411 413 415 417 422
427 437 445 446 447 451 453 454 458 466
469 471 473 478 481 482 485 489 493 497

Under 500 there are 150 MPNs and 154 semi-primes.
Under 5000 there are 1354 MPNs and 1366 semi-primes.
Under 50000 there are 12074 MPNs and 12111 semi-primes.
Under 500000 there are 108223 MPNs and 108327 semi-primes.

jq

Adapted from Wren

Works with jq, the C implementation of jq

Works with jaq, the Rust implementation of jq

Works with gojq, the Go implementation of jq (*)

(*) gojq's space requirements for computing the sieve are extravagent. 

### Generic utilities

def lpad($len): tostring | ($len - length) as $l | (" " * $l) + .;

def sum(s): reduce s as $x (0; .+$x);

# like while/2 but emit the final term rather than the first one
def whilst(cond; update):
     def _whilst:
         if cond then update | (., _whilst) else empty end;
     _whilst;

### Prime number functions

# Input:  a positive integer
# Output: an array, $a, of length .+1 such that
#         $a[$i] is $i if $i is prime, and false otherwise.
def primeSieve:
  # erase(i) sets .[i*j] to false for integral j > 1
  def erase($i):
    if .[$i] then
      reduce (range(2*$i; length; $i)) as $j (.; .[$j] = false) 
    else .
    end;
  (. + 1) as $n
  | (($n|sqrt) / 2) as $s
  | [null, null, range(2; $n)]
  | reduce (2, 1 + (2 * range(1; $s))) as $i (.; erase($i)) ;

# Number of primes up to and including $n,
# assuming . is a sufficiently large
# array as produced by primeSieve 
def primeCount($n):
  sum(.[] | select(. and . <= $n) | 1);
  

### Specialized functions

# collect at most 3 divisors
def eligible_divisors:
  . as $n
  | ($n|sqrt) as $sqrt
  | {divisors: [],
     i: 1,
     k: (if $n%2 == 0 then 1 else 2 end)
    }
  | until (.i > $sqrt;
        if .i > 1 and ($n%.i == 0)  # exclude 1 and n 
        then .divisors += [.i]
        | if (.divisors|length > 2) then .i = $sqrt + 1 # i.e. break
          else (($n/.i)|floor) as $j
          | if $j != .i then .divisors += [$j] end
          end
        end
        | .i += .k )
  | .divisors ;

# First give details for the multiplicatively-perfect numbers under $n
# then give counts for 500 * 10^k <= $upto
def task($n; $upto):
  def p: lpad(4);
  "Multiplicatively perfect numbers under \($n):",
  ( ($upto | primeSieve) as $sieved
    | { limit: 500,
       count: 0,
       line: "",
       i: 1 }
    | whilst( .limit <= $upto;
        .emit = null
        | (.i | if (. != 1) then eligible_divisors else [1, 1] end ) as $pd
        | if ($pd|length) == 2 and ($pd[0] * $pd[1]) == .i
          then .count += 1
          | if .i < $n
            then .line +=  "\(.i|lpad(3))  =\($pd[0]|p) x\($pd[1]|p)   "
            | if .count % 5 == 0
              then .emit = .line
              | .line = ""
              end
            end
          end 
          | if .i == $n - 1
            then .emit = .line + "\n"
            | .line = ""
            end
          | (.limit - 1) as $l
          | if .i == $l
            then ($sieved | primeCount($l|sqrt|floor)) as $squares
            |    ($sieved | primeCount($l|cbrt|floor)) as $cubes
            |     (.count + $squares - $cubes - 1) as $count2
            | .emit += "Counts under \(.limit|lpad(8)): MPNs = \(.count|lpad(6))  Semi-primes = \($count2|lpad(6))"
            | .limit *= 10
            end
         | .i += 1 )
    | select(.emit).emit
  );

task(500; 500000)
Output:
Multiplicatively perfect numbers under 500:
  1  =   1 x   1     6  =   2 x   3     8  =   2 x   4    10  =   2 x   5    14  =   2 x   7   
 15  =   3 x   5    21  =   3 x   7    22  =   2 x  11    26  =   2 x  13    27  =   3 x   9   
 33  =   3 x  11    34  =   2 x  17    35  =   5 x   7    38  =   2 x  19    39  =   3 x  13   
 46  =   2 x  23    51  =   3 x  17    55  =   5 x  11    57  =   3 x  19    58  =   2 x  29   
 62  =   2 x  31    65  =   5 x  13    69  =   3 x  23    74  =   2 x  37    77  =   7 x  11   
 82  =   2 x  41    85  =   5 x  17    86  =   2 x  43    87  =   3 x  29    91  =   7 x  13   
 93  =   3 x  31    94  =   2 x  47    95  =   5 x  19   106  =   2 x  53   111  =   3 x  37   
115  =   5 x  23   118  =   2 x  59   119  =   7 x  17   122  =   2 x  61   123  =   3 x  41   
125  =   5 x  25   129  =   3 x  43   133  =   7 x  19   134  =   2 x  67   141  =   3 x  47   
142  =   2 x  71   143  =  11 x  13   145  =   5 x  29   146  =   2 x  73   155  =   5 x  31   
158  =   2 x  79   159  =   3 x  53   161  =   7 x  23   166  =   2 x  83   177  =   3 x  59   
178  =   2 x  89   183  =   3 x  61   185  =   5 x  37   187  =  11 x  17   194  =   2 x  97   
201  =   3 x  67   202  =   2 x 101   203  =   7 x  29   205  =   5 x  41   206  =   2 x 103   
209  =  11 x  19   213  =   3 x  71   214  =   2 x 107   215  =   5 x  43   217  =   7 x  31   
218  =   2 x 109   219  =   3 x  73   221  =  13 x  17   226  =   2 x 113   235  =   5 x  47   
237  =   3 x  79   247  =  13 x  19   249  =   3 x  83   253  =  11 x  23   254  =   2 x 127   
259  =   7 x  37   262  =   2 x 131   265  =   5 x  53   267  =   3 x  89   274  =   2 x 137   
278  =   2 x 139   287  =   7 x  41   291  =   3 x  97   295  =   5 x  59   298  =   2 x 149   
299  =  13 x  23   301  =   7 x  43   302  =   2 x 151   303  =   3 x 101   305  =   5 x  61   
309  =   3 x 103   314  =   2 x 157   319  =  11 x  29   321  =   3 x 107   323  =  17 x  19   
326  =   2 x 163   327  =   3 x 109   329  =   7 x  47   334  =   2 x 167   335  =   5 x  67   
339  =   3 x 113   341  =  11 x  31   343  =   7 x  49   346  =   2 x 173   355  =   5 x  71   
358  =   2 x 179   362  =   2 x 181   365  =   5 x  73   371  =   7 x  53   377  =  13 x  29   
381  =   3 x 127   382  =   2 x 191   386  =   2 x 193   391  =  17 x  23   393  =   3 x 131   
394  =   2 x 197   395  =   5 x  79   398  =   2 x 199   403  =  13 x  31   407  =  11 x  37   
411  =   3 x 137   413  =   7 x  59   415  =   5 x  83   417  =   3 x 139   422  =   2 x 211   
427  =   7 x  61   437  =  19 x  23   445  =   5 x  89   446  =   2 x 223   447  =   3 x 149   
451  =  11 x  41   453  =   3 x 151   454  =   2 x 227   458  =   2 x 229   466  =   2 x 233   
469  =   7 x  67   471  =   3 x 157   473  =  11 x  43   478  =   2 x 239   481  =  13 x  37   
482  =   2 x 241   485  =   5 x  97   489  =   3 x 163   493  =  17 x  29   497  =   7 x  71   

Counts under      500: MPNs =    150  Semi-primes =    153
Counts under     5000: MPNs =   1354  Semi-primes =   1365
Counts under    50000: MPNs =  12074  Semi-primes =  12110
Counts under   500000: MPNs = 108223  Semi-primes = 108326

Julia 

using Printf
using Primes


""" Find and count multiplicatively perfect integers up to thresholds """
function testmultiplicativelyperfects(thresholds = [500, 5000, 50_000, 500_000])
    mpcount, scount = 0, 0
    pmask = primesmask(thresholds[end])
    println("Multiplicatively perfect numbers under $(thresholds[begin]):")
    for n in 1:thresholds[end]
        f = factor(n).pe
        flen = length(f)
        if flen == 2 && f[1][2] == 1 == f[2][2] || flen == 1 && f[1][2] == 3
            mpcount += 1
            if n < thresholds[begin]
                @printf("%3d * %3d = %3d   ", f[1][1], n ÷ f[1][1], n)
                mpcount % 5 == 0 && println()
            end
        end
        if n in thresholds
            cbsum, sqsum = sum(pmask[1:Int(floor(n^(1/3)))]), sum(pmask[1:isqrt(n)])
            scount = mpcount - cbsum + sqsum
            @printf("\nCounts under %d: MPNs = %7d  Semi-primes = %7d\n", n, mpcount, scount)
        end
    end
end

testmultiplicativelyperfects()
Output:

Multiplicatively perfect numbers under 500:
  2 *   3 =   6     2 *   4 =   8     2 *   5 =  10     2 *   7 =  14     3 *   5 =  15
  3 *   7 =  21     2 *  11 =  22     2 *  13 =  26     3 *   9 =  27     3 *  11 =  33
  2 *  17 =  34     5 *   7 =  35     2 *  19 =  38     3 *  13 =  39     2 *  23 =  46
  3 *  17 =  51     5 *  11 =  55     3 *  19 =  57     2 *  29 =  58     2 *  31 =  62
  5 *  13 =  65     3 *  23 =  69     2 *  37 =  74     7 *  11 =  77     2 *  41 =  82
  5 *  17 =  85     2 *  43 =  86     3 *  29 =  87     7 *  13 =  91     3 *  31 =  93
  2 *  47 =  94     5 *  19 =  95     2 *  53 = 106     3 *  37 = 111     5 *  23 = 115
  2 *  59 = 118     7 *  17 = 119     2 *  61 = 122     3 *  41 = 123     5 *  25 = 125
  3 *  43 = 129     7 *  19 = 133     2 *  67 = 134     3 *  47 = 141     2 *  71 = 142
 11 *  13 = 143     5 *  29 = 145     2 *  73 = 146     5 *  31 = 155     2 *  79 = 158
  3 *  53 = 159     7 *  23 = 161     2 *  83 = 166     3 *  59 = 177     2 *  89 = 178
  3 *  61 = 183     5 *  37 = 185    11 *  17 = 187     2 *  97 = 194     3 *  67 = 201
  2 * 101 = 202     7 *  29 = 203     5 *  41 = 205     2 * 103 = 206    11 *  19 = 209
  3 *  71 = 213     2 * 107 = 214     5 *  43 = 215     7 *  31 = 217     2 * 109 = 218
  3 *  73 = 219    13 *  17 = 221     2 * 113 = 226     5 *  47 = 235     3 *  79 = 237
 13 *  19 = 247     3 *  83 = 249    11 *  23 = 253     2 * 127 = 254     7 *  37 = 259
  2 * 131 = 262     5 *  53 = 265     3 *  89 = 267     2 * 137 = 274     2 * 139 = 278
  7 *  41 = 287     3 *  97 = 291     5 *  59 = 295     2 * 149 = 298    13 *  23 = 299
  7 *  43 = 301     2 * 151 = 302     3 * 101 = 303     5 *  61 = 305     3 * 103 = 309
  2 * 157 = 314    11 *  29 = 319     3 * 107 = 321    17 *  19 = 323     2 * 163 = 326
  3 * 109 = 327     7 *  47 = 329     2 * 167 = 334     5 *  67 = 335     3 * 113 = 339
 11 *  31 = 341     7 *  49 = 343     2 * 173 = 346     5 *  71 = 355     2 * 179 = 358
  2 * 181 = 362     5 *  73 = 365     7 *  53 = 371    13 *  29 = 377     3 * 127 = 381
  2 * 191 = 382     2 * 193 = 386    17 *  23 = 391     3 * 131 = 393     2 * 197 = 394
  5 *  79 = 395     2 * 199 = 398    13 *  31 = 403    11 *  37 = 407     3 * 137 = 411   
  7 *  59 = 413     5 *  83 = 415     3 * 139 = 417     2 * 211 = 422     7 *  61 = 427
 19 *  23 = 437     5 *  89 = 445     2 * 223 = 446     3 * 149 = 447    11 *  41 = 451
  3 * 151 = 453     2 * 227 = 454     2 * 229 = 458     2 * 233 = 466     7 *  67 = 469
  3 * 157 = 471    11 *  43 = 473     2 * 239 = 478    13 *  37 = 481     2 * 241 = 482
  5 *  97 = 485     3 * 163 = 489    17 *  29 = 493     7 *  71 = 497
Counts under 500: MPNs =     149  Semi-primes =     153

Counts under 5000: MPNs =    1353  Semi-primes =    1365

Counts under 50000: MPNs =   12073  Semi-primes =   12110

Counts under 500000: MPNs =  108222  Semi-primes =  108326

Nim

In our solution, 1 is not considered to be a multiplicatively perfect number.

Using 32 bits integers rather than 64 bits integers improves considerably the performance. With a limit set to 5_000_000 and 64 bits integers, the program runs in around 10 seconds. This time drops to 3 seconds if 32 bits integers are used. 

import std/strformat

func isMPN(n: int32): bool =
  ## Return true if "n" is a multiplicatively perfect number.
  ## We consider than 1 is not an MPN.
  var first, second = 0i32   # First and second proper divisors.
  let delta = 1 + (n and 1)
  var d = delta + 1
  while d * d <= n:
    if n mod d == 0:
      if second != 0: return false  # More than two proper divisors.
      first = d
      let q = n div d
      if q != d: second = q
    inc d, delta
  result = first * second == n

### Task ###
var count = 0
for n in 1i32..499i32:
  if n.isMPN:
    inc count
    stdout.write &"{n:3}"
    stdout.write if count mod 10 == 0: '\n' else: ' '
echo '\n'

### Stretch task ###

func isPrime(n: int32): bool =
  ## Return true if "n" is prime.
  if n < 2: return false
  if (n and 1) == 0: return n == 2
  if n mod 3 == 0: return n == 3
  var k = 5
  var delta = 2
  while k * k <= n:
    if n mod k == 0: return false
    inc k, delta
    delta = 6 - delta
  result = true

var mpnCount = 0
var limit = 500i32
var ns, nc = 3i32
var squares, cubes = 1i32
var n = 1i32
while true:
  inc n
  if n == limit:
    while ns * ns < limit:
      if ns.isPrime: inc squares
      inc ns, 2
    while nc * nc * nc < limit:
      if nc.isPrime: inc cubes
      inc nc, 2
    echo &"Under {limit} there are {mpnCount} MPNs and {mpnCount - cubes + squares} semi-primes."
    if limit == 500_000: break
    limit *= 10
  if n.isMPN: inc mpnCount

{{out]] 

  6   8  10  14  15  21  22  26  27  33
 34  35  38  39  46  51  55  57  58  62
 65  69  74  77  82  85  86  87  91  93
 94  95 106 111 115 118 119 122 123 125
129 133 134 141 142 143 145 146 155 158
159 161 166 177 178 183 185 187 194 201
202 203 205 206 209 213 214 215 217 218
219 221 226 235 237 247 249 253 254 259
262 265 267 274 278 287 291 295 298 299
301 302 303 305 309 314 319 321 323 326
327 329 334 335 339 341 343 346 355 358
362 365 371 377 381 382 386 391 393 394
395 398 403 407 411 413 415 417 422 427
437 445 446 447 451 453 454 458 466 469
471 473 478 481 482 485 489 493 497 

Under 500 there are 149 MPNs and 153 semi-primes.
Under 5000 there are 1353 MPNs and 1365 semi-primes.
Under 50000 there are 12073 MPNs and 12110 semi-primes.
Under 500000 there are 108222 MPNs and 108326 semi-primes.

Perl

Translation of: Raku
Library: ntheory
use v5.36;
use enum    <false true>;
use ntheory <is_prime nth_prime is_semiprime gcd>;

sub comma      { reverse ((reverse shift) =~ s/.{3}\K/,/gr) =~ s/^,//r }
sub table (@V) { my $t = 10 * (my $w = 5); ( sprintf( ('%'.$w.'d')x@V, @V) ) =~ s/.{1,$t}\K/\n/gr }

sub find_factor ($n, $constant = 1) {   # NB: required that n > 1
    my($x, $rho, $factor) = (2, 1, 1);
    while ($factor == 1) {
        $rho *= 2;
        my $fixed = $x;
        for (0..$rho) {
            $x = ( $x * $x + $constant ) % $n;
            $factor = gcd(($x-$fixed), $n);
            last if 1 < $factor;
        }
    }
    $factor = find_factor($n, $constant+1) if $n == $factor;
    $factor
}

# Call with range 1..$limit
sub is_mpn($n) {
    state $cube = 1; $cube = 1 if $n == 1; # set and reset
    $n == 1 ? return true : is_prime($n) ? return false : ();
    ++$cube, return true if $n == nth_prime($cube)**3;
    my $factor = find_factor($n);
    my $div    = int $n/$factor;
    return true if is_prime $factor and is_prime $div and $div != $factor;
    false
}

say "Multiplicatively perfect numbers less than 500:\n" . table grep is_mpn($_), 1..499;

say 'There are:';
for my $limit (5e2, 5e3, 5e4, 5e5, 5e6) {
    my($m,$s) = (0,0);
    is_mpn       $_ and $m++ for 1..$limit-1;
    is_semiprime $_ and $s++ for 1..$limit-1;
    printf "%8s MPNs less than %8s, %8s semiprimes\n", comma($m), $limit, comma $s
}
Output:
Multiplicatively perfect numbers less than 500:
    1    6    8   10   14   15   21   22   26   27
   33   34   35   38   39   46   51   55   57   58
   62   65   69   74   77   82   85   86   87   91
   93   94   95  106  111  115  118  119  122  123
  125  129  133  134  141  142  143  145  146  155
  158  159  161  166  177  178  183  185  187  194
  201  202  203  205  206  209  213  214  215  217
  218  219  221  226  235  237  247  249  253  254
  259  262  265  267  274  278  287  291  295  298
  299  301  302  303  305  309  314  319  321  323
  326  327  329  334  335  339  341  343  346  355
  358  362  365  371  377  381  382  386  391  393
  394  395  398  403  407  411  413  415  417  422
  427  437  445  446  447  451  453  454  458  466
  469  471  473  478  481  482  485  489  493  497

There are:
     150 MPNs less than      500,      153 semiprimes
   1,354 MPNs less than     5000,    1,365 semiprimes
  12,074 MPNs less than    50000,   12,110 semiprimes
 108,223 MPNs less than   500000,  108,326 semiprimes
 978,983 MPNs less than  5000000,  979,274 semi primes

Phix

Includes 1, to exclude just start with multiplicatively_perfect_numbers = 0 and r = {}. 

with javascript_semantics
integer multiplicatively_perfect_numbers = 1,
        semiprime_numbers = 0,
        five_e_n = 5e2
sequence r = {1}
for n=1 to 5e6 do
    sequence pn = vslice(prime_powers(n),2)
    if product(sq_add(pn,1))=4 then
        multiplicatively_perfect_numbers += 1
        if n<=500 then r &= n end if
    end if
    if sum(pn)=2 then
        semiprime_numbers += 1
    end if
    if n=five_e_n then
        if n=5e2 then
            printf(1,"%d multiplicatively perfect numbers under 500: %s\n",
                     {length(r),join(shorten(r,"",5,"%d"),",")})
        end if
        printf(1,"Counts under %,9d: MPNs = %,7d  Semi-primes = %,7d\n",
          {five_e_n,multiplicatively_perfect_numbers,semiprime_numbers})
        five_e_n *= 10
    end if
end for
Output:
150 multiplicatively perfect numbers under 500: 1,6,8,10,14,...,482,485,489,493,497
Counts under       500: MPNs =     150  Semi-primes =     153
Counts under     5,000: MPNs =   1,354  Semi-primes =   1,365
Counts under    50,000: MPNs =  12,074  Semi-primes =  12,110
Counts under   500,000: MPNs = 108,223  Semi-primes = 108,326
Counts under 5,000,000: MPNs = 978,983  Semi-primes = 979,274

PL/M

Works with: 8080 PL/M Compiler

... under CP/M (or an emulator) 

100H: /* FIND MUKTIPLICATIVELY PERFECT NUMBERS - NUMBERS WHOSE PROPER        */
      /* DIVISOR PRODUCT IS THE NUMBER ITSELF                                */
      /* NOTE THIS IS EQUIVALENT TO FINDING NUMBERS WITH 3 PROPER DIVISORS   */
      /* SEE OEIS A007422                                                    */

   /* CP/M BDOS SYSTEM CALL                                                  */
   BDOS: PROCEDURE( FN, ARG ); DECLARE FN BYTE, ARG ADDRESS; GOTO 5; END;
   /* I/O ROUTINES                                                           */
   PR$CHAR:   PROCEDURE( C ); DECLARE C BYTE;    CALL BDOS( 2, C );  END;
   PR$STRING: PROCEDURE( S ); DECLARE S ADDRESS; CALL BDOS( 9, S );  END;
   PR$NL:     PROCEDURE;   CALL PR$CHAR( 0DH ); CALL PR$CHAR( 0AH ); END;
   PR$NUMBER: PROCEDURE( N ); /* PRINTS A NUMBER IN THE MINIMUN FIELD WIDTH  */
      DECLARE N ADDRESS;
      DECLARE V ADDRESS, N$STR ( 6 )BYTE, W BYTE;
      V = N;
      W = LAST( N$STR );
      N$STR( W ) = '$';
      N$STR( W := W - 1 ) = '0' + ( V MOD 10 );
      DO WHILE( ( V := V / 10 ) > 0 );
         N$STR( W := W - 1 ) = '0' + ( V MOD 10 );
      END;
      CALL PR$STRING( .N$STR( W ) );
   END PR$NUMBER;

   /* FIND THE NUMBERS UP TO 500                                             */

   DECLARE PDC      ( 501 )ADDRESS; /* TABLE OF PROPER DIVISOR COUNTS        */
   DECLARE ( I, J, COUNT ) ADDRESS;

   /* COUNT THE PROPER DIVISORS OF NUMBERS TO 500                            */
   DO I = 0 TO LAST( PDC ); PDC( I ) = 1; END;
   DO I = 2 TO LAST( PDC );
      DO J = I + I TO LAST( PDC ) BY I;
         PDC( J ) = PDC( J ) + 1;
      END;
   END;
   PDC( 1 ) = 3;       /* PRETEND 1 HAS 3 PROPER DIVISORS SO IT IS INCLUDED  */

   /* SHOW YHE MULTIPLICATIVELY PERFECT NUMBERS                              */
   COUNT = 0;
   DO I = 1 TO LAST( PDC );
      IF PDC( I ) = 3 THEN DO;
         CALL PR$CHAR( ' ' );
         IF I < 100 THEN DO;
            IF I < 10 THEN CALL PR$CHAR( ' ' );
            CALL PR$CHAR( ' ' );
         END;
         CALL PR$NUMBER( I );
         IF ( COUNT := COUNT + 1 ) MOD 10 = 0 THEN CALL PR$NL;
      END;
   END;

EOF
Output:
   1   6   8  10  14  15  21  22  26  27
  33  34  35  38  39  46  51  55  57  58
  62  65  69  74  77  82  85  86  87  91
  93  94  95 106 111 115 118 119 122 123
 125 129 133 134 141 142 143 145 146 155
 158 159 161 166 177 178 183 185 187 194
 201 202 203 205 206 209 213 214 215 217
 218 219 221 226 235 237 247 249 253 254
 259 262 265 267 274 278 287 291 295 298
 299 301 302 303 305 309 314 319 321 323
 326 327 329 334 335 339 341 343 346 355
 358 362 365 371 377 381 382 386 391 393
 394 395 398 403 407 411 413 415 417 422
 427 437 445 446 447 451 453 454 458 466
 469 471 473 478 481 482 485 489 493 497

Raku 

use List::Divvy;
use Lingua::EN::Numbers;

constant @primes = (^∞).grep: &is-prime;
constant @cubes  = @primes.map: *³;

state $cube = 0;
sub is-mpn(Int $n ) {
    return False if $n.is-prime;
    if $n == @cubes[$cube] {
        ++$cube;
        return True
    }
    my $factor = find-factor($n);
    return True if ($factor.is-prime && ( my $div = $n div $factor ).is-prime && ($div != $factor));
    False;
}

sub find-factor ( Int $n, $constant = 1 ) {
    my $x      = 2;
    my $rho    = 1;
    my $factor = 1;
    while $factor == 1 {
        $rho *= 2;
        my $fixed = $x;
        for ^$rho {
            $x = ( $x * $x + $constant ) % $n;
            $factor = ( $x - $fixed ) gcd $n;
            last if 1 < $factor;
        }
    }
    $factor = find-factor( $n, $constant + 1 ) if $n == $factor;
    $factor;
}

constant @mpn = lazy 1, |(2..*).grep: &is-mpn;

say 'Multiplicatively perfect numbers less than 500:';
put @mpn.&upto(500).batch(10)».fmt("%3d").join: "\n";

put "\nThere are:";
for 5e2, 5e3, 5e4, 5e5, 5e6 {
   printf  "%8s MPNs less than %9s, %7s semiprimes.\n",
     comma(my $count = +@mpn.&upto($_)), .Int.&comma,
     comma $count + @primes.map(*²).&upto($_) - @cubes.&upto($_) - 1;
}
Output:
Multiplicatively perfect numbers less than 500:
  1   6   8  10  14  15  21  22  26  27
 33  34  35  38  39  46  51  55  57  58
 62  65  69  74  77  82  85  86  87  91
 93  94  95 106 111 115 118 119 122 123
125 129 133 134 141 142 143 145 146 155
158 159 161 166 177 178 183 185 187 194
201 202 203 205 206 209 213 214 215 217
218 219 221 226 235 237 247 249 253 254
259 262 265 267 274 278 287 291 295 298
299 301 302 303 305 309 314 319 321 323
326 327 329 334 335 339 341 343 346 355
358 362 365 371 377 381 382 386 391 393
394 395 398 403 407 411 413 415 417 422
427 437 445 446 447 451 453 454 458 466
469 471 473 478 481 482 485 489 493 497

There are:
     150 MPNs less than       500,     153 semiprimes.
   1,354 MPNs less than     5,000,   1,365 semiprimes.
  12,074 MPNs less than    50,000,  12,110 semiprimes.
 108,223 MPNs less than   500,000, 108,326 semiprimes.
 978,983 MPNs less than 5,000,000, 979,274 semiprimes.

Ring 

see "working..." + nl
see "Special numbers under 500:" + nl
limit = 500
Divisors = []
for n = 1 to limit
    pro = 1
    Divisors = []
    for m = 2 to ceil(n/2)
        if n % m = 0
           pro = pro * m
           add(Divisors,m)
        ok
    next
    str = ""
    if n = pro and len(Divisors) > 1
       for m = 1 to len(Divisors)
           str = str + Divisors[m] + " * "
           if m = len(Divisors)
              str = left(str,len(str)-2) 
           ok
       next
       see "" + n + " = " + str + nl
    ok
next
see "done..." + nl
Output:
working...
Special numbers under 500:
  6  =   2 x   3
  8  =   2 x   4
 10  =   2 x   5
 14  =   2 x   7
 15  =   3 x   5
 21  =   3 x   7
 22  =   2 x  11
 26  =   2 x  13
 27  =   3 x   9
 33  =   3 x  11
 34  =   2 x  17
 35  =   5 x   7
 38  =   2 x  19
 39  =   3 x  13
 46  =   2 x  23
 51  =   3 x  17
 55  =   5 x  11
 57  =   3 x  19
 58  =   2 x  29
 62  =   2 x  31
 65  =   5 x  13
 69  =   3 x  23
 74  =   2 x  37
 77  =   7 x  11
 82  =   2 x  41
 85  =   5 x  17
 86  =   2 x  43
 87  =   3 x  29
 91  =   7 x  13
 93  =   3 x  31
 94  =   2 x  47
 95  =   5 x  19
106  =   2 x  53
111  =   3 x  37
115  =   5 x  23
118  =   2 x  59
119  =   7 x  17
122  =   2 x  61
123  =   3 x  41
125  =   5 x  25
129  =   3 x  43
133  =   7 x  19
134  =   2 x  67
141  =   3 x  47
142  =   2 x  71
143  =  11 x  13
145  =   5 x  29
146  =   2 x  73
155  =   5 x  31
158  =   2 x  79
159  =   3 x  53
161  =   7 x  23
166  =   2 x  83
177  =   3 x  59
178  =   2 x  89
183  =   3 x  61
185  =   5 x  37
187  =  11 x  17
194  =   2 x  97
201  =   3 x  67
202  =   2 x 101
203  =   7 x  29
205  =   5 x  41
206  =   2 x 103
209  =  11 x  19
213  =   3 x  71
214  =   2 x 107
215  =   5 x  43
217  =   7 x  31
218  =   2 x 109
219  =   3 x  73
221  =  13 x  17
226  =   2 x 113
235  =   5 x  47
237  =   3 x  79
247  =  13 x  19
249  =   3 x  83
253  =  11 x  23
254  =   2 x 127
259  =   7 x  37
262  =   2 x 131
265  =   5 x  53
267  =   3 x  89
274  =   2 x 137
278  =   2 x 139
287  =   7 x  41
291  =   3 x  97
295  =   5 x  59
298  =   2 x 149
299  =  13 x  23
301  =   7 x  43
302  =   2 x 151
303  =   3 x 101
305  =   5 x  61
309  =   3 x 103
314  =   2 x 157
319  =  11 x  29
321  =   3 x 107
323  =  17 x  19
326  =   2 x 163
327  =   3 x 109
329  =   7 x  47
334  =   2 x 167
335  =   5 x  67
339  =   3 x 113
341  =  11 x  31
343  =   7 x  49
346  =   2 x 173
355  =   5 x  71
358  =   2 x 179
362  =   2 x 181
365  =   5 x  73
371  =   7 x  53
377  =  13 x  29
381  =   3 x 127
382  =   2 x 191
386  =   2 x 193
391  =  17 x  23
393  =   3 x 131
394  =   2 x 197
395  =   5 x  79
398  =   2 x 199
403  =  13 x  31
407  =  11 x  37
411  =   3 x 137
413  =   7 x  59
415  =   5 x  83
417  =   3 x 139
422  =   2 x 211
427  =   7 x  61
437  =  19 x  23
445  =   5 x  89
446  =   2 x 223
447  =   3 x 149
451  =  11 x  41
453  =   3 x 151
454  =   2 x 227
458  =   2 x 229
466  =   2 x 233
469  =   7 x  67
471  =   3 x 157
473  =  11 x  43
478  =   2 x 239
481  =  13 x  37
482  =   2 x 241
485  =   5 x  97
489  =   3 x 163
493  =  17 x  29
497  =   7 x  71
done...

RPL

Works with: HP version 49
≪ { 1 }
   2 ROT FOR j
      IF j SQ j DIVIS ΠLIST == THEN j + END
   NEXT
≫ 'TASK' STO

≪ (1,0) 
   2 ROT FOR j
      j DIVIS
      IF CASE 
         DUP SIZE DUP 4 > OVER " < THEN NOT END
         3 == THEN 1 END
         DUP 2 GET OVER 3 GET GCD 1 == END
      THEN SWAP (0,1) + SWAP END
      IF ΠLIST j SQ == THEN 1 + END
   NEXT
≫ 'STRETCH' STO

500 TASK
500 STRETCH
Output:
2: { 1 6 8 10 14 15 21 22 26 27 33 34 35 38 39 46 51 55 57 58 62 65 69 74 77 82 85 86 87 91 93 94 95 106 111 115 118 119 122 123 125 129 133 134 141 142 143 145 146 155 158 159 161 166 177 178 183 185 187 194 201 202 203 205 206 209 213 214 215 217 218 219 221 226 235 237 247 249 253 254 259 262 265 267 274 278 287 291 295 298 299 301 302 303 305 309 314 319 321 323 326 327 329 334 335 339 341 343 346 355 358 362 365 371 377 381 382 386 391 393 394 395 398 403 407 411 413 415 417 422 427 437 445 446 447 451 453 454 458 466 469 471 473 478 481 482 485 489 493 497 }
1: (150,153)

Rust 

fn divisors( num : u128 ) -> Vec<u128> {
   (1..= num).filter( | &d | num % d == 0 ).collect( )
}

fn main() {
   println!("{:?}" , (1..= 500).filter( | &d | {
            let divis : Vec<u128> = divisors( d ) ;
            let prod : u128 = divis.iter( ).product( ) ;
            prod == d.checked_pow( 2 ).unwrap( )
            }).collect::<Vec<u128>>( ) ) ;
}
Output:
[1, 6, 8, 10, 14, 15, 21, 22, 26, 27, 33, 34, 35, 38, 39, 46, 51, 55, 57, 58, 62, 65, 69, 74, 77, 82, 85, 86, 87, 91, 93, 94, 95, 106, 111, 115, 118, 119, 122, 123, 125, 129, 133, 134, 141, 142, 143, 145, 146, 155, 158, 159, 161, 166, 177, 178, 183, 185, 187, 194, 201, 202, 203, 205, 206, 209, 213, 214, 215, 217, 218, 219, 221, 226, 235, 237, 247, 249, 253, 254, 259, 262, 265, 267, 274, 278, 287, 291, 295, 298, 299, 301, 302, 303, 305, 309, 314, 319, 321, 323, 326, 327, 329, 334, 335, 339, 341, 343, 346, 355, 358, 362, 365, 371, 377, 381, 382, 386, 391, 393, 394, 395, 398, 403, 407, 411, 413, 415, 417, 422, 427, 437, 445, 446, 447, 451, 453, 454, 458, 466, 469, 471, 473, 478, 481, 482, 485, 489, 493, 497]


Sidef 

func multiplicatively_perfect_numbers(N) {
    [1] + semiprimes(N).grep{|n| isqrt(n**tau(n)) == n.sqr } + N.icbrt.primes.map { .cube } -> sort
}

say ("Terms <= 500: ", multiplicatively_perfect_numbers(500).join(' '), "\n")

for n in (500, 5_000, 50_000, 500_000) {
    var M = multiplicatively_perfect_numbers(n)
    say "There are #{M.len} MPNs and #{semiprime_count(n)} semiprimes <= #{n.commify}."
}
Output:
Terms <= 500: 1 6 8 10 14 15 21 22 26 27 33 34 35 38 39 46 51 55 57 58 62 65 69 74 77 82 85 86 87 91 93 94 95 106 111 115 118 119 122 123 125 129 133 134 141 142 143 145 146 155 158 159 161 166 177 178 183 185 187 194 201 202 203 205 206 209 213 214 215 217 218 219 221 226 235 237 247 249 253 254 259 262 265 267 274 278 287 291 295 298 299 301 302 303 305 309 314 319 321 323 326 327 329 334 335 339 341 343 346 355 358 362 365 371 377 381 382 386 391 393 394 395 398 403 407 411 413 415 417 422 427 437 445 446 447 451 453 454 458 466 469 471 473 478 481 482 485 489 493 497

There are 150 MPNs and 153 semiprimes <= 500.
There are 1354 MPNs and 1365 semiprimes <= 5,000.
There are 12074 MPNs and 12110 semiprimes <= 50,000.
There are 108223 MPNs and 108326 semiprimes <= 500,000.

Wren

Library: Wren-math
Library: Wren-fmt

This includes '1' as an MPN. Not very quick at around 112 seconds but reasonable for Wren. 

import "./math" for Int, Nums
import "./fmt" for Fmt

// library method customized for this task
var divisors = Fn.new { |n|
    var divisors = []
    var i = 1
    var k = (n%2 == 0) ? 1 : 2
    while (i <= n.sqrt) {
        if (i > 1 && n%i == 0) {  // exclude 1 and n 
            divisors.add(i)
            if (divisors.count > 2) break // not eligible if has > 2 divisors
            var j = (n/i).floor
            if (j != i) divisors.add(j)
        }
        i = i + k
    }
    return divisors
}

var limit = 500
var count = 0
var i = 1
System.print("Multiplicatively perfect numbers under %(limit):")
while (true) {
    var pd = (i != 1) ? divisors.call(i) : [1, 1]
    if (pd.count == 2 && pd[0] * pd[1] == i) {
        count = count + 1
        if (i < 500) {
            var pds = Fmt.swrite("$3d x $3d", pd[0], pd[1])
            Fmt.write("$3d  = $s   ", i, pds)
            if (count % 5 == 0) System.print()
        }
    }
    if (i == 499) System.print()
    if (i >= limit - 1) {
        var squares = Int.primeCount((limit - 1).sqrt.floor)
        var cubes   = Int.primeCount((limit - 1).cbrt.floor)
        var count2 = count + squares - cubes - 1
        Fmt.print("Counts under $,9d: MPNs = $,7d  Semi-primes = $,7d", limit, count, count2)
        if (limit == 5000000) return
        limit = limit * 10
    }
    i = i + 1
}
Output:
Multiplicatively perfect numbers under 500:
  1  =   1 x   1     6  =   2 x   3     8  =   2 x   4    10  =   2 x   5    14  =   2 x   7   
 15  =   3 x   5    21  =   3 x   7    22  =   2 x  11    26  =   2 x  13    27  =   3 x   9   
 33  =   3 x  11    34  =   2 x  17    35  =   5 x   7    38  =   2 x  19    39  =   3 x  13   
 46  =   2 x  23    51  =   3 x  17    55  =   5 x  11    57  =   3 x  19    58  =   2 x  29   
 62  =   2 x  31    65  =   5 x  13    69  =   3 x  23    74  =   2 x  37    77  =   7 x  11   
 82  =   2 x  41    85  =   5 x  17    86  =   2 x  43    87  =   3 x  29    91  =   7 x  13   
 93  =   3 x  31    94  =   2 x  47    95  =   5 x  19   106  =   2 x  53   111  =   3 x  37   
115  =   5 x  23   118  =   2 x  59   119  =   7 x  17   122  =   2 x  61   123  =   3 x  41   
125  =   5 x  25   129  =   3 x  43   133  =   7 x  19   134  =   2 x  67   141  =   3 x  47   
142  =   2 x  71   143  =  11 x  13   145  =   5 x  29   146  =   2 x  73   155  =   5 x  31   
158  =   2 x  79   159  =   3 x  53   161  =   7 x  23   166  =   2 x  83   177  =   3 x  59   
178  =   2 x  89   183  =   3 x  61   185  =   5 x  37   187  =  11 x  17   194  =   2 x  97   
201  =   3 x  67   202  =   2 x 101   203  =   7 x  29   205  =   5 x  41   206  =   2 x 103   
209  =  11 x  19   213  =   3 x  71   214  =   2 x 107   215  =   5 x  43   217  =   7 x  31   
218  =   2 x 109   219  =   3 x  73   221  =  13 x  17   226  =   2 x 113   235  =   5 x  47   
237  =   3 x  79   247  =  13 x  19   249  =   3 x  83   253  =  11 x  23   254  =   2 x 127   
259  =   7 x  37   262  =   2 x 131   265  =   5 x  53   267  =   3 x  89   274  =   2 x 137   
278  =   2 x 139   287  =   7 x  41   291  =   3 x  97   295  =   5 x  59   298  =   2 x 149   
299  =  13 x  23   301  =   7 x  43   302  =   2 x 151   303  =   3 x 101   305  =   5 x  61   
309  =   3 x 103   314  =   2 x 157   319  =  11 x  29   321  =   3 x 107   323  =  17 x  19   
326  =   2 x 163   327  =   3 x 109   329  =   7 x  47   334  =   2 x 167   335  =   5 x  67   
339  =   3 x 113   341  =  11 x  31   343  =   7 x  49   346  =   2 x 173   355  =   5 x  71   
358  =   2 x 179   362  =   2 x 181   365  =   5 x  73   371  =   7 x  53   377  =  13 x  29   
381  =   3 x 127   382  =   2 x 191   386  =   2 x 193   391  =  17 x  23   393  =   3 x 131   
394  =   2 x 197   395  =   5 x  79   398  =   2 x 199   403  =  13 x  31   407  =  11 x  37   
411  =   3 x 137   413  =   7 x  59   415  =   5 x  83   417  =   3 x 139   422  =   2 x 211   
427  =   7 x  61   437  =  19 x  23   445  =   5 x  89   446  =   2 x 223   447  =   3 x 149   
451  =  11 x  41   453  =   3 x 151   454  =   2 x 227   458  =   2 x 229   466  =   2 x 233   
469  =   7 x  67   471  =   3 x 157   473  =  11 x  43   478  =   2 x 239   481  =  13 x  37   
482  =   2 x 241   485  =   5 x  97   489  =   3 x 163   493  =  17 x  29   497  =   7 x  71   

Counts under       500: MPNs =     150  Semi-primes =     153
Counts under     5,000: MPNs =   1,354  Semi-primes =   1,365
Counts under    50,000: MPNs =  12,074  Semi-primes =  12,110
Counts under   500,000: MPNs = 108,223  Semi-primes = 108,326
Counts under 5,000,000: MPNs = 978,983  Semi-primes = 979,274

XPL0 

func Special(N);
int  N, D, P;
[D:= 2;  P:= 1;
while D < N do
    [if rem(N/D) = 0 then P:= P*D;
    D:= D+1;
    ];
return P = N;
];

int N, C;
[C:= 0;
Format(4, 0);
for N:= 2 to 500-1 do
    if Special(N) then
        [RlOut(0, float(N));
        C:= C+1;
        if rem(C/20) = 0 then CrLf(0);
        ];
]
Output:
   6   8  10  14  15  21  22  26  27  33  34  35  38  39  46  51  55  57  58  62
  65  69  74  77  82  85  86  87  91  93  94  95 106 111 115 118 119 122 123 125
 129 133 134 141 142 143 145 146 155 158 159 161 166 177 178 183 185 187 194 201
 202 203 205 206 209 213 214 215 217 218 219 221 226 235 237 247 249 253 254 259
 262 265 267 274 278 287 291 295 298 299 301 302 303 305 309 314 319 321 323 326
 327 329 334 335 339 341 343 346 355 358 362 365 371 377 381 382 386 391 393 394
 395 398 403 407 411 413 415 417 422 427 437 445 446 447 451 453 454 458 466 469
 471 473 478 481 482 485 489 493 497
Retrieved from "https://rosettacode.org/wiki/Multiplicatively_perfect_numbers?oldid=372209"
Cookies help us deliver our services. By using our services, you agree to our use of cookies.
More information