acer

32587 Reputation

29 Badges

20 years, 37 days
Ontario, Canada

Social Networks and Content at Maplesoft.com

MaplePrimes Activity


These are answers submitted by acer

Quite often a good way to do this is to set it up with dsolve/numeric's events functionality. The solver itself is in charge of tracking its own tolerances internally, and adjusting stepping dynamically, and letting it do that with the goal of noticing a particular event (such as hitting a zero) makes a lot of sense. Slightly less good (sometimes, b/c the accuracy of the procs emitted by dsolve/numeric are only as good as the piecewise polynomial interpolants that get set up) is to do rootfinding (fsolve, DirectSearch, etc) on option listprocedure output. And interpolating an array of data points (output from dsolve/numeric at steps chosen by the programmer) can be worse, since that data can lack any finer grained stepping that dsolve/numerics solver might have had to do near a root.

By the way, you should be able to find several older posts that relate; just use the search field at the top of this site's pages. (eg. searching for events gets here, here, etc.)

acer

If you're willing to write out the `legend` entries by hand (or do a little coding to manipulate given values of the parameter), you might try it like,

  legend=[typeset(Beta[0]=10^(-``*9))]

That should give you 9 in an upright font. Using 10^(-`9`) would show the 9 in italics, a mismatch with the upright 10. And 10^``(-9) would display with surrounding brackets.

acer

If I understand the issue rightly, then you might be able to use fsolve's `avoid` option to find (potentially) up to all three roots. And then you could just pick whichever of those had the least `t` value.

I'm not sure that I understand why Digits=180 is necessary.

The method I've outlined will have to be allowed to fail out, by which I mean that at every iteration if will have to be allowed to search until it hits its own internal number-of-attempts-limit. This means that it will run slower as a whole.

When using a type-check to test the return value from fsolve computed at the top-level, you should use eval(soln,1), If you don't then in the case that fsolve failed the returned unevaluated `fsolve` call will get re-evaluated by virtue of its being an argument in the call to `type`.

billiard_simulatio_m.mw

Let me know if I've misunderstood.

acer

Here's one way, without changing your defns for `dist` and `v`. (This is 1D input, but you should be able to do it in 2D input as well.)

restart:

v:=3.9*Unit(km/h):
dist:=t->v*t;

plot( combine(dist(t*Unit(s)),units), t = 0*Unit(s)..3*Unit(s));

acer

You've forgotten the colon in := when trying to assign to `Rounding`.

restart;

ans:=fsolve(sin(x)=0, x=2..4);

                          3.141592654

evalf[5](ans);

                             3.1416

Rounding:=-infinity:

evalf[5](ans);

                             3.1415

acer

Your subscripted thetac is a table reference. You sheet is doing the 2D Math input equivalent of this 1D input,

`θc`:=1:

`θc`;

                               1

`θc`[deg]:=4:

`θc`;
                            θc

By assigning to the table entry you clobber the previously assigned value of the base name. You can get around this by making the 2D Math input of the subscripted thetac[deg] instead be a unique name (so that it does not collide with base name thetac). To do this, highlight all of your subscripted thetac[deg] input, on the left hand side of that assignment statement, and use the right-click context-menu action 2-D Math -> Convert To -> Atomic Identifier.

If you do this then you'll have to make sure that you also use the same atomic identifier whenever you use it or refer to it later. You could repeat the context-menu action each time, or cut & paste it from earlier.

acer

This is a cross-platform regression in the Standard GUI's 2D plot renderer, introduced it seems in 16.02.

You could try reinstalling 16.01, if you are able. Or you could try Preben's CURVE splitting fix-up procedure.

It's been reported several times, in various plotting commands, since December. (eg, 1, 2, 3, 4, 5)

acer

It's not clear, right now, that you know for sure that there is in fact a root with a1:=1.514 and a2:=5.049. Perhaps you are willing to allow a small violation of the equations, in which case those could be constraints of an optimization problem for which you'd likely need a global solver.

It might be interesting to see how the other variables change, w.r.t. a1 and a2, near the previous full precision root you cited.

In your earlier thread, Carl Love mentioned controlling accuracy and precision separately. I too had been working along those lines, to investigate. I did not find a root, but it's not clear whether that is because the ranges are so great, or the functions are so nonlinear, or perhaps there is no root for those chopped a1 and a2 values.

Anyway, here is one way to have the evaluations done at higher precision while still allowing fsolve to do its usual thing and compute an accuracy tolerance based on the inbound Digits value.

restart:

read "Cvector_with_a.m":

Cnew:=eval(C,[a1=1.514,a2=5.049]):

for i from 1 to 9 do
   Cproc[i]:=codegen[optimize](
                unapply(Cnew[i],[indets(Cnew[5],name)[]])):
end do:

save Cproc, "Cproc.m":  # you only ever want to do this once

And then,

restart:

read "Cproc.m":

for i from 1 to 9 do
   f[i]:=subs(counter=i,proc(cc,C1,C2,C3,C4,C5,C6,C7,HB)
      Digits:=100:
      Cproc[counter](cc,C1,C2,C3,C4,C5,C6,C7,HB);
   end proc);
end do:

Cproc[5](1.0,
         246.851, 9.597, 7.153, 56.14, 18.91, 2.656, 7.73,
         0.012);

                                              11
                              -0.9931016321 10

f[5](1.0,
     246.851, 9.597, 7.153, 56.14, 18.91, 2.656, 7.73,
     0.012);

-0.992357302391520534067512587085549339330928664664135982173091890590934784\

                                   11
    6556567851795358311687275614 10

and so on. Of course not all those decimal digits in the result from f[5](...) above are correct. But hopefully there are enough correct (more than from Cproc[5] called at default working precision, anyway) to allow a numerical solver to get a grip and converge if it can.

The procedures f[i] use the codegen[optimized] procedures's formulas -- by virtue of calling Cproc[i] -- but they each temporarily raise Digits to 100 before making those calls.

You can then try feeding those to a rootfinder (or global optimization scheme, if that's what you want). And you could replace the hard-coded 100 by some name, so that the higher working precision was under your control.

acer

You could try,

seq([eval(t,1),length(t)], t in anames(':-user'));

or make your own more complicated version. See ?length

acer

You can apply both those rules, together or separately. Below they are applied separately. Is that what you mean?

restart:

ee:=cos(3*x)*cos(x) + sin(3*x)*sin(x):

applyrule(cos(a::integer*x)*cos(b::integer*x)
          =1/2*(cos((a-b)*x)+cos((a+b)*x)),
          ee);

           1            1                           
           - cos(2 x) + - cos(4 x) + sin(3 x) sin(x)
           2            2                           

applyrule(sin(a::integer*x)*sin(b::integer*x)
          =1/2*(cos((a-b)*x)-cos((a+b)*x)),
          ee);

                             1            1         
           cos(3 x) cos(x) + - cos(2 x) - - cos(4 x)
                             2            2         

If both are applied, then the terms collapse, and produce just cos(2*x).

acer

A := Matrix([[a, b], [b, c]], shape = symmetric, attributes = [positive_definite]):

B := LinearAlgebra:-LUDecomposition(A, method = Cholesky);
Bhat := simplify(B) assuming real, positive;
map(normal, Bhat.Bhat^%T - A );

F := LinearAlgebra:-LUDecomposition(A, method = Cholesky, conjugate=false);
map(normal, F.F^%T - A );

If you go with conjugate=false then you might still want to simplify under the assumption of being positive, depending on example.

acer

In 16.01 on 64bit Linux on an Intel i5, after about 5 minutes,

SolveTools:-SemiAlgebraic({f1+f2 = 5, f3+f4 = 11/2,
0 <f1, 0 < f2, 0 < f3, 0 < f4,
6 < f2+f3, 6 < f2+f4,
7 < f1+f4, 8 < f1+f3},
[f1,f2,f3,f4]);

[]

which, with [] as output, means no solutions.

More quickly,

SolveTools:-SemiAlgebraic({ a + b = 1, a > 0, b > 1, c > 0 }, {a,b,c});

[]

However, the computation time grows quickly with problem size, for this approach.

The `solve` command is using an older, faster, and buggier algorithm here.

acer

Instead of using a table (or list) and complexplot you could instead use a hardware datatype Matrix and call plots:-pointplot.

First, create a Matrix M like so,

   M := Matrix( 10^6, 2, datatype=float[8] ):

Then populate the first column of M with the real values and the second column with the imaginary values. If possible, alter your original code to do this directly instead of producing a table (or list) of complex values and copying values over into M by using `Re` and `Im`.

Then you could plot it with something like,

   plots:-pointplot( M );

If done right then this should be able to display the graphic of the plot in a matter of seconds rather than hours.

Your statement that you don't want an approximation is unclear. All plots turn any exact data into floating-point numbers at some point before rendering.

acer

Are you looking for something like this, perhaps?

F:=x->parse(sprintf("%.7f",x)):

F(66666.123456789);
                                66666.1234568

F(0.123456789);
                                  0.1234568

F(0.12345);
                                  0.1234500

expr:=(lambda1-2.9881355172134833*lambda2
       +16.044692204509636*lambda1^2
       -114.31727386075951*lambda1*lambda2
       +202.36051939112247*lambda2^2):

subsindets(expr,float,F);
                                                           2
           lambda1 - 2.9881355 lambda2 + 16.0446922 lambda1 

                                                                 2
              - 114.3172739 lambda1 lambda2 + 202.3605194 lambda2 

acer

If I have a hyperlink in a Maple worksheet which is of type URL and with target,

file:///C:/Users/Acer/Documents/22-19217.pdf

then clicking on it opens a new page in my web-browser. In my Maple 16's Tools->Options->General the web-browser is set as "System Default", which happens to be Chrome on my Windows 7. And my Chrome browser is configured to open PDF files inlined.

This is just simple use of the "file" URI scheme. Is that the kind of thing that you are trying to do?

Or were you hoping to have Maple itself directly open some other (dedicated) PDF reader, bypassing web-browsers, and using some mechanism other than URL?

acer

First 256 257 258 259 260 261 262 Last Page 258 of 338