What course is it for? What boundaries can you describe for the scope of the project.

I can't tell yet whether you want suggestions for topics (1st yr Calculus, or 4th yr degree project?) or suggestions on how to use or learn to use Maple.

acer

Command completion (on Linux, ctl-shift-space) does seem to work inside Math Containers.

For getting the value of a Math Container,  DocumentTools:-GetProperty seems to work. But what it produces isn't necessarily digestible by MathML:-ImportModified, if the content is some particular forms of 2D Math.

It might be interesting to know exactly what subset of 2D Math in a Math Container can be handled by ImportModified.

acer

Command completion (on Linux, ctl-shift-space) does seem to work inside Math Containers.

For getting the value of a Math Container,  DocumentTools:-GetProperty seems to work. But what it produces isn't necessarily digestible by MathML:-ImportModified, if the content is some particular forms of 2D Math.

It might be interesting to know exactly what subset of 2D Math in a Math Container can be handled by ImportModified.

acer

1) I was suggesting that perhaps there was a bug in how Maple itself was collecting garbage produced during "eval" callbacks. I am not sure whether this is the case or not. But it would not be your bug.

2) In the original code, the objective would not run under evalhf due to the `int` call. But operator form can indeed run under evalhf mode. For example,

> restart:
> infolevel[GlobalOptimization]:=2:

> Data := matrix(1,1,[1]):

> G := proc(x) x^2*Data[1,1]; end proc:

> GlobalOptimization:-GlobalSolve(G,-1..1);
GlobalSolve:   calling NLP solver
SolveGeneral:   calling global optimization solver
SolveGeneral:   number of problem variables   1
SolveGeneral:   number of nonlinear inequality constraints   0
SolveGeneral:   number of nonlinear equality constraints   0
SolveGeneral:   trying evalhf mode
SolveGeneral:   total number of function evaluations   1035
SolveGeneral:   runtime in external solver   0.
SolveGeneral:   maximum constraint infeasibility   0.
                                  [0., [0.]]

Note the message about "trying evalhf mode above". Notice also this following behaviour,

> evalhf(G(2.0));
                                      4.
> H := proc(x)
>   x * int(sin(t),t=0..1);
> end proc:
> evalf( H(2.0) );
                                  0.919395388

> evalhf( H(2.0) );
Error, unable to evaluate function `int` in evalhf

3) I have to correct myself slightly, about `data`. In your original example, it's not a barrier to evalhf. But if one encapsulated the example above within a procedure, then it would be. For example,

> restart:
> infolevel[GlobalOptimization]:=2:

> f := proc()
>   local Data, G;
>   Data := matrix(1,1,[1]);
>   G := proc(x) x^2*Data[1,1]; end proc;
>   GlobalOptimization:-GlobalSolve(G,-1..1);
> end proc:

> f();
GlobalSolve:   calling NLP solver
SolveGeneral:   calling global optimization solver
SolveGeneral:   number of problem variables   1
SolveGeneral:   number of nonlinear inequality constraints   0
SolveGeneral:   number of nonlinear equality constraints   0
SolveGeneral:   trying evalhf mode
SolveGeneral:   trying evalf mode
SolveGeneral:   total number of function evaluations   1035
SolveGeneral:   runtime in external solver   0.
SolveGeneral:   maximum constraint infeasibility   0.
                                  [0., [0.]]

Notice how, in the above example, there is also a message, "trying evalf mode". The evalhf callback failed. And so would this,

> restart:

> f := proc()
>   local Data, G;
>   Data := matrix(1,1,[1]);
>   G := proc(x) x^2*Data[1,1]; end proc;
>   evalhf(G(2.0));
> end proc:

> f();
Error, (in f) lexical scoping is not supported in evalhf

You can check whether your procedure evaluates under evalhf, say by calling it on a sample point within an evalhf() call like I've done above. It can get a little tricky though. The Optimization or GlobalOptimization solvers may switch modes altogether if they encounter just one evaluation point which fails under evalhf. So success at just a single test point wouldn't assure that the solver only uses evalhf mode.

> G := proc(x)
>   if x<5 then
>     sin(x);
>   else
>     x * int(cos(t),t=Pi..x);
>   end if;
> end proc:
> evalhf(G(4.0)); # ok so far, but...
                             -0.756802495307928202

> infolevel[Optimization]:=2:
> Optimization:-Minimize(G,1..7);
NLPSolve:   calling NLP solver
NLPSolve:   using method=quadratic
NLPSolve:   number of problem variables
NLPSolve:   trying evalhf mode
NLPSolve:   trying evalf mode
attemptsolution:   number of evaluations taken   42
attemptsolution:   final value of objective   -4.79462137048670024
                 [-4.79462137048670024, [5.00000000615824103]]

> evalhf(G(7.0));
Error, unable to evaluate function `int` in evalhf

acer

Isn't InverseSurvivalFunction(Y, p) the same as Quantile(Y, 1-p), at least for continuous distribution functions?

Quantile is sometimes known as the inverse cumulative distribution function (or icdf as it's called in Matlab).

I wonder what, if anything, are the "standard" names for these things.

acer

Isn't InverseSurvivalFunction(Y, p) the same as Quantile(Y, 1-p), at least for continuous distribution functions?

Quantile is sometimes known as the inverse cumulative distribution function (or icdf as it's called in Matlab).

I wonder what, if anything, are the "standard" names for these things.

acer

I don't know for a fact that your problem is due to a memory leak. But the evidence presented so far made me consider it -- it would explain the evidence so far presented.

When one runs GlobalOptimization or Optimization routines, and set infolevel appropriately, once can see messages indicating that either evalhf or eval callbacks are being tried. A callback is a call from an external routine back to Maple, for the purpose of evaluating an objective, constraint, or integrand for example Now, evalhf is a fast, alternate Maple interpreter. It can't handle everything. But what it can handle it does quickly. It appears to first interpret most everything as a hardware double precision float, so that it can evaluate expressions very quickly while producing very little garbage.

Your code example produced a lot of intermediary garbage (objects no longer referenced or needed), with all those symbolic int calls. Using your orginal code, those symbolic integrals would all have been computed inside eval callbacks from the external GlobalSolve engine back to Maple. Normally, that garbage would get collected automatically when Maple did its periodic garbage collection. My hypothesis is that, perhaps, some of this garbage is not being properly collected. The thing that's different about this garbage in particular is that it is being formed during eval callbacks to Maple from the compiled external optimization engine. It may be this difference which is the key to the behaviour.

Using objectives/constraints/integrands which can be evaluated under evalhf will make some Maple routines faster. And it can also reduce the amount of symbolic garbage produced. Garbage is not always a bad thing. Usually is can be collected and reclaimed. But if you suspect a leak related to such symbolic garbage, then ensuring that your functions evaluate ok under evalhf may be one workaround.

A memory leak is a bug. So there is no good answer to "why isn't that memory freed up again", in the case of a leak, except to say that it is by mistake.

acer

I see, your example was deliberately contrived so as to simulate heavy memory use.

Following what John wrote, one might reasonably suspect a memory leak. I could imagine that as possible, during software float (not evalhf) callbacks from an external solver to Maple for each objective and constraint evaluation. I do not know how to get around that, except to rewrite the code similarly to what I suggested above, on a case by case basis and where possible.

As for upgrading your machine or memory, well, that will only help you so far. Whatever code is consuming memory will likely continue to do so, perhaps at the same rate, unless you can change it.

You might find this link useful, if you use a Microsoft OS.

acer

The bytes used and gc (and hence timing for the whole worksheet) is likely very high because of what goes on inside the objective function.

The objective function is this below, and it  is run for likely tens or hundreds of thousands of pairs of numeric values for parameters M1 and C1.

model_res:=proc(M1,C1)
  local res:
  res:=add((data[i,2]-
           subs({m=M1,c=C1},
                int(m*data[i,1]/(1+t+t^2)+c*t,t=0..5)))^2,
           i=1..RowDimension(data)):
  return(res):
end:

The data object is being set up a lowercase "a" array (lowercase "m" matrix) and is being passed to RowDimension which is a command for uppercase "M" Matrices. So probably a needless conversion occurs with each invocation. And it should be LinearAlgebra:-RowDimension and not rely upon LinearAlgebra being loaded at the top-level.

The most serious memory usage issue, I suspect, is that the symbolic integral is being done for each M1,C1 pair. And for each such given pair it is being done for each i. That integral can be done once, for generic C1,M1, and i, to provide a formula. Maple can evaluate such a formula efficiently.

If the objective function were also "fixed" to not rely on lexical scoping of 'data' then it might even run under fast evalhf mode. GlobalOptimization would try that, and handle it, automagically. But even if that were not so, evaluating a formula would be leaner than doing nested `int` calls.

> int(m*data[i,1]/(1+t+t^2)+c*t,t=0..5);
                                                                1/2
                      1/2                        1/2        11 3       25 c
   -1/9 m data[i, 1] 3    Pi + 2/3 m data[i, 1] 3    arctan(-------) + ----
                                                               3        2

I suggest rewriting `model_res` so that it accepts parameters m,c, and data. And have data be a datatype=float[8] Matrix. And have it use the formula for the integral rather than call `int` even once.

acer

And more generally, if one wishes to obtain P{X>x} for many other points x then integration is not necessary at each of them. One may instead compute a symbolic integral just once, to obtain a function of x, or simply use the CDF as Doug pointed out. And it may also be done for (as yet) unknown mean and variance.

> with( Statistics ):

> X := RandomVariable( Normal(m,s) ):

> PXms := int( PDF(X,t), t=x..infinity )
>   assuming m>0, m<infinity, s>0, s<infinity;
                                      1/2
                                     2    (-x + m)
                     PXms := 1/2 erf(-------------) + 1/2
                                          2 s

> Px := unapply(eval(PXms,[m=3,s=3]),x);
                                          1/2
                  Px := x -> 1/2 erf(1/6 2    (-x + 3)) + 1/2
 
> evalf(Px(0));
                                 0.8413447460

And, of course, this is the same as what one gets from the CDF,

> X := RandomVariable( Normal(m,s) ):

> 1-CDF(X,x);
                                        1/2
                                       2    (-x + m)
                         1/2 + 1/2 erf(-------------)
                                            2 s

acer

And more generally, if one wishes to obtain P{X>x} for many other points x then integration is not necessary at each of them. One may instead compute a symbolic integral just once, to obtain a function of x, or simply use the CDF as Doug pointed out. And it may also be done for (as yet) unknown mean and variance.

> with( Statistics ):

> X := RandomVariable( Normal(m,s) ):

> PXms := int( PDF(X,t), t=x..infinity )
>   assuming m>0, m<infinity, s>0, s<infinity;
                                      1/2
                                     2    (-x + m)
                     PXms := 1/2 erf(-------------) + 1/2
                                          2 s

> Px := unapply(eval(PXms,[m=3,s=3]),x);
                                          1/2
                  Px := x -> 1/2 erf(1/6 2    (-x + 3)) + 1/2
 
> evalf(Px(0));
                                 0.8413447460

And, of course, this is the same as what one gets from the CDF,

> X := RandomVariable( Normal(m,s) ):

> 1-CDF(X,x);
                                        1/2
                                       2    (-x + m)
                         1/2 + 1/2 erf(-------------)
                                            2 s

acer

That reply is much more useful and instructive for the student.

There could be more pedagogical difference between VectorCalculus:-ArcLength and Student:-VectorCalculus:-ArcLength.

acer

That reply is much more useful and instructive for the student.

There could be more pedagogical difference between VectorCalculus:-ArcLength and Student:-VectorCalculus:-ArcLength.

acer

I would probably do it like this,

> F:=proc(L)
>   subsop(1=NULL,L);
> end proc:
> L:=[x,y,z,d];
                               L := [x, y, z, d]
 
> L:= F(L);
                                L := [y, z, d]
 
> L:= F(L);
                                  L := [z, d]
 
> L:= F(L);
                                   L := [d]
 
> L:= F(L);
                                    L := []

That keeps it as simple as possible, is easiest to understand, and is less confusing.

Doing it the other way, to bring about the change in L as a so-called side-effect on the name, can cause confusion. In a sense, it's deliberately subverting any protection Maple may otherwise give you about not having your procedure overwrite the formal parameter. As a general rule, you may not wish to do that unless it's necessary. For your example, the above way to overwrite L (explicitly, by the call to the proc) is more straightforward.

Way back in the day, before Maple allowed multiple assignments, some routines would use side-effects as a mechanism to fake multiple assignment. For example,

> f := proc(a,b)
>   b := 2*a;
>   3*a;
> end proc:

> x:=11:

> x := f(x,'y'):
> x,y;
                                    33, 22
 
> x := f(x,'y'):
> x,y;
                                    99, 66

The above procedure changes both x and y. There was a time, long ago, when that was one of the easier ways to do get that effect. Nowadays, it can be done much more simply with a multiple assignment,

> restart:

> f := proc(a)
>  3*a, 2*a;
> end proc:

> x:=11:

> x,y := f(x);
                                x, y := 33, 22
 
> x,y := f(x);
                                x, y := 99, 66

When using side effects, one can get confused the by quotes, by either forgetting them or inserting them when they are not wanted.

> restart:

> f := proc(a,b)
>   b := 2*a;
>   3*a;
> end proc:

> x:=11:

> x := f(11,y): # it works once, without quotes...
> x,y;
                                    33, 22
 
> x := f(11,y): # and now, the second time...
Error, (in f) illegal use of a formal parameter

And, doing it another way,

> f := proc(x::evaln)
>   x := 2*eval(x);
> end proc:

> x:=11:

> f(x):
> x;
                                      22
 
> f('x'):
Error, illegal use of an object as a name

And also, documenting its use for others is more involved if it's implemented in these ways.

acer

I would probably do it like this,

> F:=proc(L)
>   subsop(1=NULL,L);
> end proc:
> L:=[x,y,z,d];
                               L := [x, y, z, d]
 
> L:= F(L);
                                L := [y, z, d]
 
> L:= F(L);
                                  L := [z, d]
 
> L:= F(L);
                                   L := [d]
 
> L:= F(L);
                                    L := []

That keeps it as simple as possible, is easiest to understand, and is less confusing.

Doing it the other way, to bring about the change in L as a so-called side-effect on the name, can cause confusion. In a sense, it's deliberately subverting any protection Maple may otherwise give you about not having your procedure overwrite the formal parameter. As a general rule, you may not wish to do that unless it's necessary. For your example, the above way to overwrite L (explicitly, by the call to the proc) is more straightforward.

Way back in the day, before Maple allowed multiple assignments, some routines would use side-effects as a mechanism to fake multiple assignment. For example,

> f := proc(a,b)
>   b := 2*a;
>   3*a;
> end proc:

> x:=11:

> x := f(x,'y'):
> x,y;
                                    33, 22
 
> x := f(x,'y'):
> x,y;
                                    99, 66

The above procedure changes both x and y. There was a time, long ago, when that was one of the easier ways to do get that effect. Nowadays, it can be done much more simply with a multiple assignment,

> restart:

> f := proc(a)
>  3*a, 2*a;
> end proc:

> x:=11:

> x,y := f(x);
                                x, y := 33, 22
 
> x,y := f(x);
                                x, y := 99, 66

When using side effects, one can get confused the by quotes, by either forgetting them or inserting them when they are not wanted.

> restart:

> f := proc(a,b)
>   b := 2*a;
>   3*a;
> end proc:

> x:=11:

> x := f(11,y): # it works once, without quotes...
> x,y;
                                    33, 22
 
> x := f(11,y): # and now, the second time...
Error, (in f) illegal use of a formal parameter

And, doing it another way,

> f := proc(x::evaln)
>   x := 2*eval(x);
> end proc:

> x:=11:

> f(x):
> x;
                                      22
 
> f('x'):
Error, illegal use of an object as a name

And also, documenting its use for others is more involved if it's implemented in these ways.

acer