h=1E-34;
 m=9.1E-31;
 V=0;
 En=0.7E-17;
 dx=0.000000000001;
 x=100;
 y=100;
 A=2*m*(En-V)/h^2*dx^2;
 %g = @(x,shi)[shi(1); -A*shi(2)];
 %[x, shix] = ode45(g, [0:100*dx], [1 0]);
 shix = zeros(100);
 shiy = zeros(100);
 shi = zeros(100,100);
 shix(1)=0;
 shix(2)=0.1;
 shiy(1)=0;
 shiy(2)=0.1;
 %pertamax
 for i=3:x
     shix(i)=(2-A)*shix(i-1)-shix(i-2);
     shiy(i)=(2-A)*shiy(i-1)-shiy(i-2);
 end
 for i=1:x
     for j=1:y
         shi(i,j) = shix(i)*shiy(j);
     end
 end
 %keduax
 akir=shi(x,y);
 Eni=En;
 En=En*0.9;
 A=2*m*(En-V)/h^2*dx^2;
 for i=3:x
     shix(i)=(2-A)*shix(i-1)-shix(i-2);
     shiy(i)=(2-A)*shiy(i-1)-shiy(i-2);
 end
 for i=1:x
     for j=1:y
         shi(i,j) = shix(i)*shiy(j);
     end
 end 
 %cari energi
 while shi(x,y)~=akir
     endi=En;
     En=En-shi(x,y)*(Eni-En)/(akir-shi(x,y));
     Eni=endi;
     akir=shi(x,y);
     A=2*m*(En-V)/h^2*dx^2;
    for i=3:x
        shix(i)=(2-A)*shix(i-1)-shix(i-2);
        shiy(i)=(2-A)*shiy(i-1)-shiy(i-2);
    end
    for i=1:x
        for j=1:y
            shi(i,j) = shix(i)*shiy(j);
        end
    end
 end
 surf(shi)
 En;