Write program which verifies Erdős–Straus conjecture

Ruby, 119 106 characters

f=->s,c,a{m=s.to_i;c<2?m<s||(p a+[m];exit):(1+m...c*s).map{|k|f[s/(1-s/k),c-1,a+[k]]}}
f[gets.to_r/4,3,[]]

The code uses minimal bounds for each variable, e.g. n/4<x<3n/4, similarly for y. Even the last example returns instantaneous (try here).

Examples:

> 12
[4, 13, 156]

> 123
[31, 3814, 14542782]

> 1234
[309, 190654, 36348757062]

> 40881241801
[10220310451, 139272994276206121600, 22828913614743204775214996005450198400]

Mathematica 62

This plain-vanilla solution works fine--most of the time.

f@n_ := FindInstance[4/n == 1/x + 1/y + 1/z && 0 < x < y < z, {x, y, z}, Integers]

Examples and Timings (in secs)

AbsoluteTiming[f[63]]
AbsoluteTiming[f[123]]
AbsoluteTiming[f[1003]]
AbsoluteTiming[f[3003]]
AbsoluteTiming[f[999999]]
AbsoluteTiming[f[1000000]]

{0.313671, {{x -> 16, y -> 1009, z -> 1017072}}}
{0.213965, {{x -> 31, y -> 3814, z -> 14542782}}}
{0.212016, {{x -> 251, y -> 251754, z -> 63379824762}}}
{0.431834, {{x -> 751, y -> 2255254, z -> 5086168349262}}}
{1.500332, {{x -> 250000, y -> 249999750052, z -> 1201920673328124750000}}}
{1.126821, {{x -> 375000, y -> 1125000, z -> 2250000}}}

---

But it does not constitute a complete solution. There are a some numbers that it cannot solve for. For example,

AbsoluteTiming[f[30037]]
AbsoluteTiming[f[130037]]

{2.066699, FindInstance[4/30037 == 1/x + 1/y + 1/z && 0 < x < y < z, {x, y, z}, Integers]}
{1.981802, FindInstance[4/130037 == 1/x + 1/y + 1/z && 0 < x < y < z, {x, y, z}, Integers]}

C#

Disclamer: this is not a serious answer

This just bruteforces all the possibilities from 1 to 1<<30. It's huge, it's slow, I don't even know if it works correctly, but it follows the specifications quite literally, as it checks the condition every single time, so that's nice. I haven't tested this because ideone has a 5 second time limit for programs and therefore this won't finish executing.

(In case anyone was wondering: this is a whopping 308 bytes long)

static double[]f(double n)
{
    for(double x=1;x<1<<30;x++)
    {
        for(double y=1;y<1<<30;y++)
        {
            for(double z=1;z<1<<30;z++)
            {
                if(4/n==1/x+1/y+1/z)
                    return new[]{x,y,z};
            }
        }
    }
    return null;
}

Update: fixed it so it actually works