additive_primes: procedure options (main);
  %replace
    search_limit by 500,
    true by '1'b,
    false by '0'b;
  dcl (i, count) fixed bin;
  put skip edit('Searching up to ', search_limit,
     ' for additive primes') (a,f(3),a);
  put skip;
  count = 0;
  do i = 2 to search_limit;
     if isprime(i) then
        do;
          if isprime(sumdigits(i)) then
             do;
                put edit(i) (f(5));
                count = count + 1;
                if mod(count,8) = 0 then put skip;
             end;
        end;
  end;
  put skip edit(count, ' were found') (f(3), a);

/* return true if n is prime */
isprime: proc(n) returns (bit(1));
   dcl
     (n, i, limit) fixed bin;
   if n < 2 then return (false);
   if mod(n, 2) = 0 then return (n = 2);
   limit = floor(sqrt(n));
   i = 3;
   do while ((i <= limit) & (mod(n, i) ^= 0));
     i = i + 2;
   end;
   return (i > limit);
end isprime;

/* return the sum of the digits of n */
sumdigits: proc(n) returns (fixed bin);
   dcl
     (n, nn, sum) fixed bin;
   /* use copy, since n is passed by reference */
   nn = n;
   sum = 0;
   do while (nn > 0);
     sum = sum + mod(nn, 10);
     nn = nn / 10;
   end;
   return (sum);
end sumdigits;

end additive_primes;