isCarmichael.m

function y = isCarmichael(x)
    y = false
    
    if isprime(x)
        return
    end

    for b = 2:x-1
        if gcd(b,x) ~= 1; continue; end
        
        p = 1;

        for i = 1:x-1
            p = mod(p*b,x);
        end
        
        if p ~= 1; return; end
    end
 
    y = true
end