What is the data assigned to? What is assigned to x, y, and z?

I'm guessing that you have assigned to `x` the Vector of the first independent variable's values, and to `y` the Vector of the second independent variable's values, and to `z` the Vector of the dependent variable's values.

Did you instead mean something like this?

Fit(a+b.x+c.y+d.x^2+e.x.y+f.y^2, Matrix([x,y]),z,[x,y]);

acer

Try this,

> with(student): with(RealDomain):

> F := int(sqrt(u),u=0..t^10);
                                            5   15
                             F := 2/3 csgn(t ) t

> simplify(F) assuming t::real;
                                           15
                                  2/3 | t |

> df := diff(F,t);
                                     5   15            5   14
                  df := 2/3 csgn(1, t ) t   + 10 csgn(t ) t
 
> evalc(df);
                                           5  14
                               10 signum(t)  t

> simplify(%) assuming t::real;
                                         13
                                 10 | t |   t

Maple also came up with the following, though it took much longer than it ought to have done.

> F := int(sqrt(u),u=0..t^10) assuming t::real;
                                      15
                            F := 2/3 t   signum(t)

I originally gave a flippant answer to this question. But, under the RealDomain environment, the assumption that t is real doesn't seem strained. Also, Maple didn't seem able to simplify(df) nicely and immediately under the assumption that t was real. Hence I made intermediate use of evalc. Having to use `assuming` and evalc underneath RealDomain doesn't look so good.

acer

What happens when x<0 ? What about when x>=0 and x<1?

> limit(f,s=0) assuming x>=1;
                                       0

acer

It could be that the Float(infinity) result you see is due to the working precision. In such a case, it might help to increase the working precision. For example,

> 1/(1.000000000001 - 1.0);
                                Float(infinity)
 
> Digits:=20:
> 1/(1.000000000001 - 1.0);
                                                   13
                          0.10000000000000000000 10

Or it may be that the result is a good representation of the answer, in some limiting concept. Consider,

> evalf(Int(1/(2-x),x=1..2)); # or int()
                                Float(infinity)

Note that in the above, the integrand cannot be evaluated at the end-point value of 2. You might think of that input as a shorthand notation for,

> evalf( limit(Int(1/(2-x),x=1..b), b=2, left) );
                                Float(infinity)

That's a mathematical way of saying that the area under the curve 1/(2-x) from 1 to b tends to infinity as b tends to 2 from the left. To replace maple's result of infinity by some large finite value M would be wrong, because there is always some b~ closer to 2 for which the area from 1 to b~ would be greater than M.

Nobody can write down explicitly a number with an infinite number of digits, or draw an infinite area (unscaled). But there are lots of mathematical problems for which infinity is nevertheless the best way to represent the answer (as a limiting concept).

Without knowing your example, it's hard to say more that may be of help to you.

acer

> V := Vector(10,i->3*i+11);
                                        [14]
                                        [  ]
                                        [17]
                                        [  ]
                                        [20]
                                        [  ]
                                        [23]
                                        [  ]
                                        [26]
                                   V := [  ]
                                        [29]
                                        [  ]
                                        [32]
                                        [  ]
                                        [35]
                                        [  ]
                                        [38]
                                        [  ]
                                        [41]

> V1,V2 := V[1..5],V[6..-1];
                                       [14]  [29]
                                       [  ]  [  ]
                                       [17]  [32]
                                       [  ]  [  ]
                             V1, V2 := [20], [35]
                                       [  ]  [  ]
                                       [23]  [38]
                                       [  ]  [  ]
                                       [26]  [41]

See ?interface and the item on rtablesize, for controlling the dimension above which Vectors/Matrix get printed compactly. Eg,

interface(rtablesize=9);
Vector[row](10);
interface(rtablesize=10);
Vector[row](10);

See ?LinearAlgebra,General,MVextract for more about Vector "extraction".

Sorry, I don't have any ideas about problems specific (only) to printing. You might consider contacting support@maplesoft.com to ask Technical Support.

acer

> simplify(rationalize(d));
                                       1/2
                                      2
                                    - ----
                                       2

> radnormal(d);
                                       1/2
                                      2
                                    - ----
                                       2

acer

You create a module, which you give the name you want for the package. And you give it 'option package'.

For example,

foo := module()
option package;
export f1, f2;
  f1 := proc() "stuff" end proc;
  f2 := proc() "other stuff" end proc;
end module:

with(foo);
f1();
f2();

acer

The awful results from this particular posted use of the GlobalOptimization toolbox are due to the fact that, as the problem was originally set up, the solver would only be able to do a relatively small number of global search evaluations before it timed out. The reason for this is because dsolve(...,numeric) was being called an enormous number of times. The act of setting up the dsolve/numeric solver was itself being done inside the objective function, at each evaluation of the objective.

Joe's responses show how the dsolve calls may be brought outside the evaluation routines, to be done just once as a set-up step.

Once the objective is made fast enough, and can run enough times to to the global optimizer a chance to do enough evaluations, it can find decent results.

I consider Joe's responses to be excellent.  I also feel that these techniques could be better documented in the product. A worksheet could illustrate the ideas nicely. And that might be usefully referenced from the help-pages of both dsolve/numeric, Optimization, and GlobalOptimization. One could construct a problem for which the numeric DE solving was involved within a root-finding problem, and thus the parametrized numeric DE solver needed to be called from fsolve.

acer

If you don't mind waiting for the entire worksheet to re-execute, then the "!!!" button or top menubar's "Edit"->"Execute"->"Worksheet" actions should suffice.

However, some key parts of a worksheet might take a great deal of time to execute, if heavy computations are involved. And re-execution of the whole worksheet might then be undesirable. Sometimes one would like to be able to change certain values, or plot options, etc, and re-run only some subset of an existing worksheet.

So, there's a distinction, between expensive parts that you don't want to re-execute, and interactive parts that you might like to call once again. Re-executing the entire worksheet doesn't allow for this distinction.

The commands `save` and `anames` can be used to save Maple's "state" (or parts of it). The state can be saved into a file, and then re-read quickly back into a new session with an explicit `read` command. But such state files are distinct from worksheets. The worksheet itself doesn't have the complete state bundled up within it.

The Standard GUI can prompt you about whether to bundle large amounts of plot data (that sit behind the scenes) into a worksheet. That's a bit of an exception to the rule.

Hmm. I wonder whether one could create an auto-opening .mla archive, whose load action reads a .m state file as well as opens a worksheet. Maybe it'd have to rely on -kp or a shared kernel.

acer

It's pretty straightforward to toggle on units support at the interative top-level. Note that the SI system is the default. In SI, the Newton (N) is the base unit for representing force. So output with that dimension comes out in terms on Unit(N).

In 1D Maple input, you can just type it in directly, as in the example below.

> with(Units[Standard]):
> 35*Unit(kN) + 7*Unit(MN);
                                  7035000 [N]

If you instead want kN to be the default base for output, so that the above comes back as 7035 [[kN]], then that too can be done with a customized system. To see an example of that, enter this Maple command to see a help-system page,

?Task,ChangeDefaultUnitInSystem

Alternatively, you can convert output like 7035000 [[N]] to 7035 [[kN]] using the context-menus. Select the output with the mouse, and when the context-menu pops up choose "Units" -> "Replace Units".

In 2D Math output, the Unit(N) object should get typeset as something that looks like [[N]]. To get that in 2D Math input, you can use the "Units (SI)" palette. Or use command-completion.

See the top menubar's "View" -> "Palettes" -> "Arrange Palettes" for a popup menu to configure which palettes show up in the left sidebar.

acer

It'd help if you posted all the Maple code, and your exact commands used to convert it. You can upload it here, as a worksheet.

acer

What you asked for was this,

(quote) but ideally what I would like to be able to do is

plot( Re(subs(deltaE=0.355,M=sqrt(5.4e5),Eig1(R^3))), R^3=1^(1/3)..10^(1/3));

(end quote)

That is, you asked to plot f(R^3) vs R^3, and that is just a question of labelling.

You didn't actually ask how to plot f(R) vs R^3. Is this what you would like?

Eig1 := R -> (1/9*(-27*deltaE^3*R^9-10*M^6+3*
  sqrt(-189*deltaE^4*R^12*M^4-147*deltaE^2*R^6*M^8
  -27*M^12+60*deltaE^3*R^9*M^6))^(1/3)
  -(-deltaE^2*R^6-7/9*M^4)/(-27*deltaE^3*R^9-
  10*M^6+3*sqrt(-189*deltaE^4*R^12*M^4
  -147*deltaE^2*R^6*M^8-27*M^12+60*deltaE^3*R^9*M^6))^(1/3)
  +2/3*deltaE*R^3+2/9*M^2)/R^3;

plot( Re(subs(deltaE=0.355,M=sqrt(5.4e5),Eig1(t^(1/3)))),
      t=1..10,labels=[R^3,'Eig1'(R)]);

acer

I too got the hypergeom form from Maple 12 (and 11.02, and 10.02). Which is a pity.

> In := Int(sqrt(t^5+6*t)*(5*t^4+6), t=0..1):

> value(In);
10  1/2
-- 6    hypergeom([-1/2, 11/8], [19/8], -1/6)
11
 
          1/2
     + 4 6    hypergeom([-1/2, 3/8], [11/8], -1/6)

The example (from a textbook it seems) is of course easy to do by hand. And it can be forced from Maple as well, by doing something as simple as the obvious change of variables.

> student[changevar](u=t^5+6*t,In,u);
                                    7
                                   /
                                  |    1/2
                                  |   u    du
                                  |
                                 /
                                   0

> value(%);
                                        1/2
                                    14 7
                                    -------
                                       3

Of course, as delta7 noted, the indefinite integration result could also be used. And for fun,

> sol := value(In);
       10  1/2
sol := -- 6    hypergeom([-1/2, 11/8], [19/8], -1/6)
       11
 
          1/2
     + 4 6    hypergeom([-1/2, 3/8], [11/8], -1/6)

> sol2 := value(simplify(combine(convert(convert(sol,MeijerG),Sum))));
                 54  1/2                       169    81
         sol2 := -- 6    hypergeom([-1/2, 3/8, ---], [--, 19/8], -1/6)
                 11                            88     88
 

Now, how can one get back to 13/3*sqrt(7) from that? Should we try to look them up? Better still if `int` could use a simpler method (such as easy substitutions) and get the simpler answer directly.

acer

I don't not have access to the Global Optimization Toolbox right this moment (weekend). I'll try to look, when I can. Not everyone here will have access to that add-on.

Remember that the discipline of global optimization is very difficult in general. It is an art. I have seen that toolbox work wonderfully and outperform other products.

I believe that the default method for the toolbox is branchandbound. Try changing it to multistart (or maybe even singlestart). Each method will work better than the rest for some different class of problems. For highly nonlinear functions, the multistart often does better.

No method works best for all global optimization problems. No optimizer (that I have seen) works better than all other on all problems.

Also, in global optimization the "best" solution is hard to prove, if one cannot do it analytically and has thus resorted to numerical optimizers. There might always be another (better) local minimum just over the next ridge. So it can often be a question of, "how long can you afford to wait, until you have to accept the curent best". In consequence of this you may wish to adjust the timelimit.

Try to keep Digits to at least 14 (which trunc(evalhf(Digits)) on Windows) so that Maple can evaluate objective and constraints better (in case of numerical instability of the problem).

If speed is an issue, try to ensure that objective and constraints are "evalhf'able".

Also, try setting infolevel[GlobalOptimization]:=1 or 2 to try and get an idea of how well it's doing.

You mentioned another product (scop) producing a better fit in 2 minutes. Did you have to compile the objective and constraints, when using that package?

acer

Your angles1_s from the .m file is a sequence of 10 lists. It appears that you are trying to build a new sequence, by concatenating the contents of the 10 lists in angles1_s except for their last elements.

But you have a nested set of calls. The "problem" is that your inner `seq` call is walking all the elements of the inner lists simply to extract each element. I guess that it's just not efficient to do it like that.

Does this not do the same thing? (It's fast.)

seq(op(angles1_s[j][1..-2]),j=1..10):

Now, you mentioned that you will eventually want to do something more. Maybe that means acting on these values in some way, while you form the concatenated sequence. This should do that (replace x->x with your function of choice),

seq( op(map(x->x,angles1_s[j][1..-2])), j=1..10):

acer