The following generates a 3dpoint plot of the maxima on the 3dplot.  If you rotate the plot you see that some of the points are just above the surface and some just below. I suppose this is due to slight numerical inacuracy finding the maxima.

spq:=proc(w) # simple procedure to ease the generation of points to plot
local q,prq;
q:=diff(subs(ga = w, TSales),delta):  
prq:=proc(x) return eval(q,delta=x); end proc:
fsolve(prq(z)=0,z);
end proc :


pts:={seq([m,spq(m),eval(TSales,{ga=m,delta=spq(m)})],m=0.0..0.5,0.1)}:
pp:=plots:-pointplot3d(pts,symbol=solidsphere ,color=black,symbolsize=10):
p1:=plot3d(TSales, delta = 0 .. 2,ga = 0 .. .5, axes = boxed):
plots:-display(pp,p1);

@hisasugi8338 

Yes, you should be able to. You need to create a set of points, and then use plots:-pointplot3d to plot them. Then use plots:-display() to display both the 3dplot and the 3dpointplot together,

@hisasugi8338 

Yes, you should be able to. You need to create a set of points, and then use plots:-pointplot3d to plot them. Then use plots:-display() to display both the 3dplot and the 3dpointplot together,

Excuse me. I forgot to explain how to actually find the maxima !

So, the problem is that the expression involves an integral in anover variable, and we can't just use fsolve(diff(abc,delta)=0,delta). There may be a better way to do it but one way is to wrap the expression in a procedure, like this:

For a(t)

qq:=diff(a(t),delta):  
pq:=proc(x) global qq;return eval(qq,delta=x); end proc:
fsolve(pq(z)=0,z);
                          1.062593448

Excuse me. I forgot to explain how to actually find the maxima !

So, the problem is that the expression involves an integral in anover variable, and we can't just use fsolve(diff(abc,delta)=0,delta). There may be a better way to do it but one way is to wrap the expression in a procedure, like this:

For a(t)

qq:=diff(a(t),delta):  
pq:=proc(x) global qq;return eval(qq,delta=x); end proc:
fsolve(pq(z)=0,z);
                          1.062593448

I don't know why the 2nd plot is missing on the above post. I'll try it again with a classic mode worksheet

Download TruncNormal1-1.mw

I don't know why the 2nd plot is missing on the above post. I'll try it again with a classic mode worksheet

Download TruncNormal1-1.mw

Yes, the truncation is effectively shifted by t. To see the attached document with 2 plots

The mathematical principle is the concept of "function of a random variable" - I've appended it, as I understand it, briefly in the document also. I'm not sure of the best way to present "nice format" math in maple as this is the first time I tried it in a document. Hope it's readable, but if you just google for "function of a random variable" you should find plenty of references.

Edit1: Note that f(x) and f'(x) are hard to distinguish above :(
Edit2: f must be monotonic for this to work, otherwise it doesn't have a unique inverse. If f is not monotonic, then the usual way to proceed is to split it pieces that are monotonic 

Download TruncNormal1.mw

Yes, the truncation is effectively shifted by t. To see the attached document with 2 plots

The mathematical principle is the concept of "function of a random variable" - I've appended it, as I understand it, briefly in the document also. I'm not sure of the best way to present "nice format" math in maple as this is the first time I tried it in a document. Hope it's readable, but if you just google for "function of a random variable" you should find plenty of references.

Edit1: Note that f(x) and f'(x) are hard to distinguish above :(
Edit2: f must be monotonic for this to work, otherwise it doesn't have a unique inverse. If f is not monotonic, then the usual way to proceed is to split it pieces that are monotonic 

Download TruncNormal1.mw

@epostma , I've followed the instructions in the link you gave and now have MS Visual C++ 2010 Express installed. I've also followed the instructions regarding the bat file, and run Compiler:-Setup() to point to it.

However, I still get this error

Error, (in Compiler:-Compile) compiler exited with nonzero status 1: 

 This occurs with the simplest example,

p := proc( x :: float ) :: float; 2.3 * x end proc:
cp:=Compiler:-Compile(p);

so I can only assume there is something wrong in my setup.

Any ideas how to track down the problem ?

Thanks

LR

Could you provide some more details ? In particular, what is T ? Perhaps you can post all the relevant code.

@epostma : It seems I don't have the compiler installed (64 bit windows). I shall look through the link you gave soon. Thanks again. 

@epostma 

Thanks again. But in version 4 I get the following error:

Error, (in Compiler:-Compile) compiler exited with nonzero status 1: 

@Christopher2222 

I don't know about graphically, I can think of two ways to test this mathematically: first, a test based on counts: for example the χ² dispersion test; and second, a test based on distances: for example Hopkins' test.  

If I have some spare time this weekend I might try to implement these for the case in hand.

@epostma thank you very much for your comments and the improved code. That's a terrific efficiency gain :-)

I didn't really understand version 3. Can you enlighten me ? And why is it so much faster that version 2 ?

I have another question: how do you determine the time spent in different parts of the procedure ?

And just a trivial point:  You seem to have left out 

uses S = Statistics;

from version 3 : I guess you had S declared at the global level in your worksheet ?

Thanks again

LR