@chomchom The code follows:
 

restart;
schro := {diff(psi(x), x, x)-(alpha*x^4+x^2-energy)*psi(x) = 0};
ic := {psi(3) = 0, (D(psi))(3) = 1};
Ic := [{psi(3) = 0, (D(psi))(3) = -1}, %$2];
E := [1.06538, 3.306, 5.74796];
schro1 := [seq(subs(energy = e, alpha = .1, schro), e = E)];
##
soln1 := [seq(dsolve(schro1[i] union Ic[i], {psi(x)}, numeric,output=listprocedure), i = 1 .. nops(E))]; #output=listprocedure is used instead of default.
##
for i from 1 to nops(E) do N[i]:=evalf(Int(subs(soln1[i],psi(x))(x)^2,x=-3..3,epsilon=1e-8))^(1/2) end do;
##
with(plots):
display(seq(odeplot(soln1[i], [x, psi(x)/N[i]], -3 .. 3, color = [red, blue, green][i]), i = 1 .. nops(E)));
## And extending to -7..7 (without the normalization) as I commented on later:
##
display(seq(odeplot(soln1[i], [x, psi(x)], -7 .. 7, color = [red, blue, green][i]), i = 1 .. nops(E)));

My guess is that in your quantum mechanical model the Hilbert space is L²(-3..3) and not L²(-infinity..infinity).
## As it turns out dsolve can find the solutions without resorting to numerics.
Here the first one as an example.
sol:=dsolve(schro1[1] union Ic[1]);
nm:=evalf(Int(rhs(sol)^2,x=-3..3))^(1/2); #Using numerical integration for the norm.
plot(rhs(sol)/nm,x=-3..3);
## I'm not necessarily implying that you should avoid the numerical approach since the expression in rhs(sol) is rather complicated.

This extrapolation problem to which I have given an answer below, reminds me that I never got any reaction from you to my answer to your question here:

http://www.mapleprimes.com/questions/219344-How-To-Solve-Odeproblem

I would like to hear your reaction to the answer I gave there.

@9009134 In my recollection pdf-files should be no problem.

Here is a test file. It has only one line.
test.pdf

Notice that when you try to upload, you are told that valid extensions are:
doc, xls, txt, zip, pdf, jpg, jpeg, gif, png, mw, mws, mla, ppt, html, htm, msim, docx, xlsx, pptx, qu

@Mac Dude Yes, that arctan(0,0) = 0 follows from the definition of argument(0) as 0, which is written in the help page for argument (here z:=x+I*y):

"For the case where y = 0 , argument(z) = 0 if  0 <= x."

@Carl Love Yes, I just noticed. My first deletion after the new interface appeared (you guys must have been busy deleting).
When you don't see that the deleted item is gone and if there are several other items basically looking the same this could be annoying.

@9009134 I try to keep all my dealings with MaplePrimes matters entirely in this forum. Furthermore, if the reference material you have take up that much space I'm afraid I wouldn't take the time to read it or even just browse it.

@loramina123 You give us h(0)=0 D(h)(0)=tehta (theta, I suppose), D(D(h))(infinity)=1/R and say that infinity might be taken as 2 or 3 micron. By a micron I suppose that (in this context with no units) you mean 1e-6.
But then you also say that you have h(2e-6)=200e-9.  But if infinity can be taken as e.g. 2e-6 then you have a conflict of enormous proportions since R = 1e-7, so 1/R=1e7.
Finally you mention D(D(H)(-infinity)=0 and that -infinity can be equal to 2e-6. I suppose that you meant -2e-6?

In these matters one has to be more careful or else people give up.

@chomchom Are you sure that you want to integrate from -infinity to infinity?

If you are, then try e.g.

display(seq(odeplot(soln1[i], [x, psi(x)], -7 .. 7, color = [red, blue, green][i]), i = 1 .. nops(E)))

If you think that is fine then replace -7..7 by -19 .. 19.

@loramina123 So what are the boundary conditions? You have presented more than one set of conditions. The ode appears to be the same though.
The original 3-point boundary condition problem may be doable, but you do have to use 4 conditions for a 4th order ode. So two conditions could be given at 0, e.g. values of h(0) and D(h)(0). Besides one condition at each of the other two points. By using smoothness at the middle point this can (sometimes) be done.
I did as a rather trivial example the following, where an exact solution can also be found. Thus we can compare the numerical result with the exact.

ode:=diff(y(x),x$4)+4*y(x)=0;
sol:=dsolve({ode,y(0)=0,D(y)(0)=1,y(Pi/2)=0,y(Pi)=0}); #This gives you a simple formula for the solution

By doing what I sketched above I got very good agreement between numerical and exact results.

Before attempting anything similar to your problem I need the actual boundary conditions if you are still interested in the 3-point problem.

@loramina123 As Tom Leslie has examined it is difficult if not impossible to solve the revised problem by dsolve/numeric/bvp.

The original problem involved conditions at 3 points, say a,b,c. Was that the real problem you wanted to solve?
Certainly it cannot be done in a straightforward way by dsolve/numeric/bvp, but maybe you could split the problem into two boundary value problems for the same ode but on the two intervals a..b and on b..c.
Extra conditions at b would be the requirement that the solution be smooth (also) at b.


Do your revised conditions also reflect a change in the actual (physical) problem?

Your code is extremely difficult to read since it is probably some copy of 2D-math input.

But you cannot use psi(infinity)= whatever as a boundary condition in dsolve/numeric/bvp. You must replace infinity by a real number.

@YasH I tried the following variant of your code:

for k from 1 to 4 do [Transpose(Eigenvalues(M)),Transpose(op(1, [Eigenvectors(M)] ))] end do;

and found that the eigenvalues in the first version of the eigenvalues vector (the one from Eigenvalues) always came in the same order for the given matrix M whereas this was certainly not the case for the second version.
I have no idea why.

@9009134 Do you have a reference for this problem? I should like to see how SYS was derived.
I notice that the coefficients in SYS are either 1, -1, 1/2, or 3/2, i.e. there aren't any values that come from physical data as opposed to quite a lot of bvp problems presented on MaplePrimes.

It is safe to say that none of us like images if we have to examine some problem. We want code in text or an uploaded worksheet: Use the big fat green arrow in the MaplePrimes editor to upload a worksheet.

@9009134 My point in trying revised boundary conditions was actually to make you examine the physical situation of which (presumably) this ode bvp is a mathematical model. Could the boundary conditions be different; why impose f''(1)=0? Is there a good physical reason?
Obviously, I wouldn't have tried revising the bcs if we hadn't run into problems like hints at nonuniqueness with the original bcs.