#Haskell: 87 bytes

Performs modular exponentiation x^y mod n in a number of operations which is logarithmic in the exponent y. The code

(x!y)m|odd y=z q|y>0=q|0<1=1where z 0=0;z n=until(<m)(+(-m))$x+z(n-1);q=(z x!div y 2)m

defines an operator (!) to do modular exponentiation, so you can use it like

(2!3) 5 == 3

The program works by using the two relations x^{2 n} = (x²)ⁿ and x^{2 n + 1} = x (x²)ⁿ to accumulate the result, halving the exponent at each step. Intermediate results are taken mod m by repeated subtraction. Here's an ungolfed version:

(x ^% y) m
    | y == 0     =  1
    | odd y      =  x `times` (squaredX ^% halfY) m
    | otherwise  =  (squaredX ^% halfY) m
    where squaredX = x `times` x
          halfY = y `div` 2
          a `times` 0 = 0
          a `times` b = reduce $ a + a `times` (b - 1)
          reduce n = if n < m then n else reduce (n - m)

Note: div is used to carry out a right shift operation only.

Haskell: 87 bytes

Performs modular exponentiation x^y mod n in a number of operations which is logarithmic in the exponent y. The code

(x!y)m|odd y=z q|y>0=q|0<1=1where z 0=0;z n=until(<m)(+(-m))$x+z(n-1);q=(z x!div y 2)m

defines an operator (!) to do modular exponentiation, so you can use it like

(2!3) 5 == 3

The program works by using the two relations x^{2 n} = (x²)ⁿ and x^{2 n + 1} = x (x²)ⁿ to accumulate the result, halving the exponent at each step. Intermediate results are taken mod m by repeated subtraction. Here's an ungolfed version:

(x ^% y) m
    | y == 0     =  1
    | odd y      =  x `times` (squaredX ^% halfY) m
    | otherwise  =  (squaredX ^% halfY) m
    where squaredX = x `times` x
          halfY = y `div` 2
          a `times` 0 = 0
          a `times` b = reduce $ a + a `times` (b - 1)
          reduce n = if n < m then n else reduce (n - m)

Note: div is used to carry out a right shift operation only.

#Haskell: 87 bytes

Performs modular exponentiation x^y mod n in a number of operations which is logarithmic in the exponent y. The code

(x!y)m|odd y=z q|y>0=q|0<1=1where z 0=0;z n=until(<m)(+(-m))$x+z(n-1);q=(z x!div y 2)m

defines an operator (!) to do modular exponentiation, so you can use it like

(2!3) 5 == 3

The program works by using the two relations x^{2 n} = (x²)ⁿ and x^{2 n + 1} = x (x²)ⁿ to accumulate the result, halving the exponent at each step. Intermediate results are taken mod m by repeated subtraction. Here's an ungolfed version:

(x ^% y) m
    | y == 0     =  1
    | odd y      =  x `times` (squaredX ^% halfY) m
    | otherwise  =  (squaredX ^% halfY) m
    where squaredX = x `times` x
          halfY = y `div` 2
          a `times` 0 = 0
          a `times` b = reduce $ a + a `times` (b - 1)
          reduce n = if n < m then n else reduce (n - m)

#Haskell: 87 bytes

Performs modular exponentiation x^y mod n in a number of operations which is logarithmic in the exponent y. The code

(x!y)m|odd y=z q|y>0=q|0<1=1where z 0=0;z n=until(<m)(+(-m))$x+z(n-1);q=(z x!div y 2)m

defines an operator (!) to do modular exponentiation, so you can use it like

(2!3) 5 == 3

The program works by using the two relations x^{2 n} = (x²)ⁿ and x^{2 n + 1} = x (x²)ⁿ to accumulate the result, halving the exponent at each step. Intermediate results are taken mod m by repeated subtraction. Here's an ungolfed version:

(x ^% y) m
    | y == 0     =  1
    | odd y      =  x `times` (squaredX ^% halfY) m
    | otherwise  =  (squaredX ^% halfY) m
    where squaredX = x `times` x
          halfY = y `div` 2
          a `times` 0 = 0
          a `times` b = reduce $ a + a `times` (b - 1)
          reduce n = if n < m then n else reduce (n - m)

Note: div is used to carry out a right shift operation only.