如何利用matlab求出10000之内所有的水仙花数（如153=1^3+5^3+3^3，这种数称为

水仙花数只是自幂数的一种，严格来说三位数的3次幂数才成为水仙花数。
如153=1^3+5^3+3^3，这种数称为水仙花数。

若是10000以内，则必有四位数，四位自幂数不叫水仙花数，而叫做四叶玫瑰数。
现给一种简单的求水仙花数的代码，自创的。
x=0:9;
y=x;z=x;
[x,y,z]=meshgrid(x,y,z);
u=x.^3+y.^3+z.^3-100*x-10*y-z;
k=find(u==0&x~=0);
w=[x(k)';y(k)';z(k)']
sxhs=[100*w(1,:)+10*w(2,:)+w(3,:)]

w =

3 3 1 4
7 7 5 0
0 1 3 7

sxhs =

370 371 153 407

当然你也可以用循环来做
for m=100:999
m1=fix(m/100);
m2=rem(fix(m/10),10);
m3=rem(m,10);
if m==m1*m1*m1+m2*m2*m2+m3*m3*m3
disp(m)
end
end
输出结果为：
153 370 371 407

如果非要求四位数的四叶玫瑰数，有
for m=1000:9999
m1=fix(m/1000);
m2=rem(fix(m/100),10);
m3=rem(fix(m/10),10);
m4=rem(m,10);
if m==m1^4+m2^4+m3^4+m4^4
disp(m)
end
end

输出结果为： 1634 8208 9474

你也可以将后面的两步统一，即
for m=100:9999
m1=fix(m/1000);
m2=rem(fix(m/100),10);
m3=rem(fix(m/10),10);
m4=rem(m,10);
if (m==m1^4+m2^4+m3^4+m4^4)|m==m2^3+m3^3+m4^3
disp(m)
end
end

输出结果为：
153

370

371

407

1634

8208

9474

水仙花数（Narcissistic number），也被称为超完全数字不变数（pluperfect digital invariant, PPDI）、自恋数、自幂数、阿姆斯壮数或阿姆斯特朗数（Armstrong number） ，用来描述一个N位非负整数，其各位数字的N次方和等于该数本身。

若将条件放宽，一个N位数，其各个数之M次方和等于该数，M和N不一定相等，这样的数称为完全数字不变量（perfect digital invariant, PDI），例如数字4150等于各位数字的5次方。

水仙花数一定是完全数字不变量，但完全数字不变量不一定是水仙花数。在十进制意义下，所有的个位数字0~9都是水仙花数，但没有二位的水仙花数。

用MATLAB求水仙花数可以用下面的简单程序：

N=1:10000;
n=floor(log10(max(N(:))))+1;
isnarcissistic=@(N)sum(mod(fix(N./10.^(n:-1:0)),10).^(floor(log10(max(N(:))))+1))==N;
inx=arrayfun(isnarcissistic,N);
N(inx)

程序的核心是构造一个判断是否为水仙花数的匿名函数isnarcissistic，使用的时候只需要修改最前面的N，其它都不需要改动。

输出结果：

ans =
  Columns 1 through 5
           1           2           3           4           5
  Columns 6 through 10
           6           7           8           9         153
  Columns 11 through 15
         370         371         407        1634        8208
  Column 16
        9474

clear
clc
y=[];
N=9999;
a=[0 0 0 0];
for i=1:N
a(1)=floor(i/1000); %取第一个数字
a(2)=mod(floor(i/100),10);%第二
a(3)=floor(mod(i,100)/10);%第三
a(4)=mod(i,10);%第四

if i<=9
t=1;%取几次方

elseif 10<=i&i<=99
t=2;
elseif 100<=i&i<=999

t=3;
elseif 1000<=i&i<=9999
t=4;
end
if i==a(1)^t+a(2)^t+a(3)^t+a(4)^t %如果相同 就放入y中
y=[y,i];
end
end
y%输出y

%%%楼主采纳 ！！！谢谢 不懂再问！！！！！！！