witam pisze program w Maple, który metodą strzałów ma liczyć równanie schrodingera, cząstka umieszczona w kwadratowej studni potencjałów, niestety program w kilku miejscach i nie wysypuje, czy ktoś by pomógł?

 > restart;with(plots):  hbar:=6.6260755E-34;  Me:=9.1093987E-31;  ec:=1.60217733E-19;  ck:=evalf((hbar*1E9/(2*Pi))^2 /(2*Me * ec));
                                    -34
                        6.6260755 10   
                                    -31
                        9.1093987 10   
                                     -19
                        1.60217733 10   
                         0.03809980339
> 
> depth := 4; width := 2;
                               4
                               2
> V := -depth*(1-Heaviside(abs(x)-width/2));  g1:=plot(V,x=-width/2-1..width/2+1,color=blue,axes=framed,title=`V vs x`,labels=[nm,eV]): #`see if V looks ok`  g1;
                   -4 + 4 Heaviside(|x| - 1)

> eq1 := -ck*(diff(y(x), x, x))+V*y(x) = eV*y(x);
                       / d  / d      \\
        -0.03809980339 |--- |--- y(x)||
                       \ dx \ dx     //

           + (-4 + 4 Heaviside(|x| - 1)) y(x) = eV y(x)
> s_odd:= k -> dsolve({subs(eV=k,eq1),D(y)(0)=1,y(0)=0},y(x),numeric);  s_even:= k -> dsolve({subs(eV=k,eq1),D(y)(0)=0,y(0)=1},y(x),numeric);
 k -> dsolve({subs(eV = k, eq1), y(0) = 0, D(y)(0) = 1}, y(x), 

   numeric)
 k -> dsolve({subs(eV = k, eq1), y(0) = 1, D(y)(0) = 0}, y(x), 

   numeric)
> test_point:=width/2+4;  r_odd:= p -> (rhs(s_odd(p)(test_point)[2]));  r_even:= p -> (rhs(s_even(p)(test_point)[2]));
                               5
p -> rhs(s_odd(p)(test_point)[2])
p -> rhs(s_even(p)(test_point)[2])
> _oddl:= p -> 3*sign(rhs(s_odd(p)(test_point)[2]));  r_evenl:= p -> 3*sign(rhs(s_even(p)(test_point)[2]));
p -> 3 sign(rhs(s_odd(p)(test_point)[2]))
p -> 3 sign(rhs(s_even(p)(test_point)[2]))
> plot(r_oddl,-depth..0,y=-2..2,color=blue,labels=[eV,` `],ytickmarks=0,title=`odd energy levels`);  plot(r_evenl,-depth..0,y=-2..2,color=blue,labels=[eV,` `],ytickmarks=0,title=`even energy levels`);
Warning, unable to evaluate the function to numeric values in the region; see the plotting command's help page to ensure the calling sequence is correct


> levels := [even, -3.9], [odd, -3.7], [even, -3.3], [odd, -2.8], [even, -2.1], [odd, -1.3], [even, -.4];
     [even, -3.9], [odd, -3.7], [even, -3.3], [odd, -2.8], 

       [even, -2.1], [odd, -1.3], [even, -0.4]
> numlevels := nops([levels]);
                               7
> dr := proc (x) options operator, arrow; (r(x+0.1e-5)-r(x))/0.1e-5 end proc;
     r(x + 0.000001) - r(x)
x -> ----------------------
            0.000001       
> for j to numlevels do  end do;
> level:=levels[j][1];  mid:=levels[j][2];
Error, invalid subscript selector
Error, invalid subscript selector
> if level='odd'  then r:=r_odd; s:=s_odd;  else r:=r_even; s:=s_even;  fi;
                             r_even
                             s_even
> mid0:=0: counter:=0:  while (abs(mid-mid0)>0.000000002) and (counter < 10) do  mid0:=mid:  mid:=mid0-r(mid0)/dr(mid0):  print(mid);  counter:=counter + 1:  od:
Error, cannot determine if this expression is true or false: 0.2e-8 < abs(mid)
> energy:=convert(mid,string);  g2:=odeplot(s(mid),[x,y(x)],-width/2-1..width/2+1,axes=framed,title=`Ev = `.energy,labels=[nm,eV]):  print(display([g2,g1]));
                             "mid"
Warning, The use of global variables in numerical ODE problems is deprecated, and will be removed in a future release. Use the 'parameters' argument instead (see ?dsolve,numeric,parameters)
Warning, cannot evaluate the solution past the initial point, problem may be complex, initially singular or improperly set up
Warning, cannot evaluate the solution past the initial point, problem may be complex, initially singular or improperly set up

> print(s(mid.j));  mid.j:=mid;  s.j:=s(mid.j);  print(Ev=mid.j);  od:
Error, unable to parse
printApplyFunction(sApplyFunction(midperiodj))semi  midperiodj

  Assignmidsemi  speriodjAssignsApplyFunction(midperiodj)semi  

  printApplyFunction(Evequalsmidperiodj)semi  

  Typesetting:-mambiguous(od, Typesetting:-merror(unable to parse)

  )colon
> for j from 1 to numlevels do  s_abs.j:=(x) -> rhs(s.j(x)[2])^2;  energy:=convert(mid.j,string);  plot(s_abs.j,-width/2-1/2..width/2+1/2,title=`original |psi|^2 Ev=`.energy);  area.j:=evalf(Int(s_abs.j,x=-width/2-1/2..width/2+1/2));  od;
                                                  2
x -> rhs(Typesetting:-delayDotProduct(s, j(x)[2])) 
                             "mid"
Warning, unable to evaluate the function to numeric values in the region; see the plotting command's help page to ensure the calling sequence is correct

Error, (in evalf/int) invalid arguments
> for j from 1 to numlevels do  s_norm.j:= (x) -> rhs(s.j(x)[2])/sqrt(area.j);  energy:=convert(mid.j,string):  g2:=plot(s_norm.j,evalf(-width/2-1..width/2+1),labels=[nm,eV],title=`normalized psi Ev=`.energy,axe  s=framed):  print(display([g2,g1]));  od:
Error, (in plot) unexpected option: axe*s_even = framed