acer

32405 Reputation

29 Badges

19 years, 350 days
Ontario, Canada

Social Networks and Content at Maplesoft.com

MaplePrimes Activity


These are answers submitted by acer

You could assign the string (full and explicit name of the .mla file, with path) to a Maple name. Then you could pass that very same name to both commands, and in doing so avoid mismatch mistakes.

For example,

filenm := "C:\\Home\\Maple projects\\Test_package\\test.mla":
LibraryTools:-Create(filenm);
LibraryTools:-Save('testSavelib', filenm);

The above correction would also work in your original code, and avoid your mismatch between the your two lines related to the extra whitespace (that the second line had, while the first did not).

Regardless of the nature of your particular problem, I think that using LibraryTools:-Save is generally better (less error prone) than using savelib.

Why aren't you using fsolve for this?

CodeTools:-Usage(fsolve(1 + 0.15 = (1 + r__x)^1.28858));
memory used=1.33MiB, alloc change=0 bytes, cpu time=16.00ms, real time=16.00ms, gc time=0ns

              0.1145625357

Is there a numeric value for parameter a?

From your grainy image it looks like you may have accidentally only made the equation,

   a = 1

instead of doing the assignment,

   a := 1

It's more useful to upload and attach your actual worksheet (green up-arrow in the Mapleprimes editor) than it is to show only an image of your code.

Are you hoping to utilize all of those equations rC1,rC2,rf1,rA,rh to substitute for all of C1,C2,f1,A, and h?

restart

with(plots):

``

psi := I*(C1*exp(A)-C2*exp(-A))*exp(-I*t*(1/2))/f1;

I*(C1*exp(A)-C2*exp(-A))*exp(-((1/2)*I)*t)/f1

rf1 := f1 = sqrt(h^2-1);

f1 = (h^2-1)^(1/2)

rh := h = f^2+1;

h = f^2+1

rA := `assuming`([A = simplify(subs(rh, sqrt(h^2-1)*(x+I*h*t)))], [f > 0]);

A = f*(f^2+2)^(1/2)*(I*f^2*t+I*t+x)

rC1 := `assuming`([C1 = simplify(subs(rf1, rh, sqrt(h-f1)))], [f > 0]);

C1 = (f^2+1-f*(f^2+2)^(1/2))^(1/2)

rC2 := `assuming`([C2 = simplify(subs(rf1, rh, sqrt(h+f1)))], [f > 0]);

C2 = (f^2+1+f*(f^2+2)^(1/2))^(1/2)

psi_f := eval[recurse](psi, [rC1, rC2, rf1, rA, rh]);

I*((f^2+1-f*(f^2+2)^(1/2))^(1/2)*exp(f*(f^2+2)^(1/2)*(I*f^2*t+I*t+x))-(f^2+1+f*(f^2+2)^(1/2))^(1/2)*exp(-f*(f^2+2)^(1/2)*(I*f^2*t+I*t+x)))*exp(-((1/2)*I)*t)/((f^2+1)^2-1)^(1/2)

simplify(limit(psi_f, f = 0, right));

exp(-((1/2)*I)*t)*((2*I)*x-I-2*t)

simplify(limit(psi_f, f = 0, left));

-exp(-((1/2)*I)*t)*((2*I)*x-I-2*t)

 

NULL

Download limit1.mw

At lower (and default) working precision this example can encounter enough numeric round-off error that the symbolic solve command finds spurious roots if you perform non-inert symbolic solving and integration on the expression containing floating-point coefficients.

I suggest:
1) don't use RealDomain (it's not doing what you think, anyway, in your example).
2) use fsolve for numeric root-finding instead of solve here.
3) Don't introduce floating-point values instead of exact rationals if you can practically avoid doing so, especially if you intend on attempting symbolic integration or solving. (Symbolic integration/solving of an expression containing floats can sometimes do worse than trying to solve it numerically, or trying to solve an exact problem.)

If you left out RealDomain, and used exact rationals, then the spurious solutions from solve (due to floating-point round-off error) could be avoided by raising Digits. But it'd be better and faster to simply use a numeric root-finder for this problem which likely doesn't even have an explicit symbolic solution.

restart

f := proc (x) options operator, arrow; 6-sqrt(-x^2+8*x+9) end proc

eq := 10 = int(f(x), x = 0 .. L)

10 = -6-(25/2)*arcsin(4/5)+6*L-(1/2)*(-L^2+8*L+9)^(1/2)*L+2*(-L^2+8*L+9)^(1/2)-(25/2)*arcsin((1/5)*L-4/5)

fsolve(eq, L)

6.810993081

fsolve(eq, L = -10 .. 10, maxsols = 5)

6.810993081

Download Working_ac.mw

 

If you really want to see more of the specifics of how it went wrong...

restart

f := proc (x) options operator, arrow; 6-sqrt(-x^2+8*x+9) end proc

proc (x) options operator, arrow; 6-sqrt(-x^2+8*x+9) end proc

eq := 10 = int(f(x), x = 0 .. L)

10 = -6-(25/2)*arcsin(4/5)+6*L-(1/2)*(-L^2+8*L+9)^(1/2)*L+2*(-L^2+8*L+9)^(1/2)-(25/2)*arcsin((1/5)*L-4/5)

oof := 10 = int(evalf(f(x)), x = 0 .. L)

10 = -17.59119023+6.*L-.5000000000*(-1.*L^2+8.*L+9.)^(1/2)*L+2.*(-1.*L^2+8.*L+9.)^(1/2)-12.50000000*arcsin(.2000000000*L-.8000000000)

feq := evalf(eq)

10. = -17.59119022+6.*L-.5000000000*(-1.*L^2+8.*L+9.)^(1/2)*L+2.*(-1.*L^2+8.*L+9.)^(1/2)-12.50000000*arcsin(.2000000000*L-.8000000000)

solve(eq, L)

5*sin(RootOf(25*(cos(_Z)^2)^(1/2)*sin(_Z)-60*sin(_Z)+25*_Z+25*arcsin(4/5)-16))+4

RealDomain:-solve(feq, L)

8.855609656, -.4725737212, 6.810993078

solve(feq, L)

8.855609656, -.4725737212, 6.810993078, -.7764912837+15.69863487*I, 8.138450513-15.36358890*I, -.7764912837-15.69863487*I, 8.138450513+15.36358890*I, 2.341723694-6.481846912*I, 2.341723694+6.481846912*I, 13.34312523-18.81939387*I, -5.700438002+19.09985409*I, 13.34312523+18.81939387*I, -5.700438002-19.09985409*I, 16.71699261-17.45692870*I, -9.003420055+17.75423373*I, 16.71699261+17.45692870*I, -9.003420055-17.75423373*I

If we approximate eq at higher working precision then terms like arcsin(4/5) don't incur as must loss of accuracy (which propagates forward as greater inaccuracy, leading to spurious solutions).
Digits := 30; solve(evalf(eq))

6.81099308113251113676562942220

NULL

Download Not_working_ac.mw

You can call fsolve to compute a root between each of those.

But don't try and do it by appending (to a list, or Vector). That's a poor way to program, and habitually leads to unnecessary inefficiency in programs as they scale.

Instead, you could simply use seq (or map) to form the list. Or you could use Vector with a constructor operator.

restart

gamma1 := .1093:

f := q-(1-alpha)*tan(q)-c*tan(q)

q-0.4646295e-3*tan(q)-0.2298931352e-4*tan(q)/(0.3800000000e-3*q^2-0.2300000000e-4)

n := 10:

Asymptotes

fOdd := proc (j) options operator, arrow; (1/2)*(2*j+1)*Pi end proc; OddAsymptotes := Vector[row](n, fOdd)

fOdd := proc (j) options operator, arrow; (j+1/2)*Pi end proc

OddAsymptotes := Vector[row](10, {(1) = (3/2)*Pi, (2) = (5/2)*Pi, (3) = (7/2)*Pi, (4) = (9/2)*Pi, (5) = (11/2)*Pi, (6) = (13/2)*Pi, (7) = (15/2)*Pi, (8) = (17/2)*Pi, (9) = (19/2)*Pi, (10) = (21/2)*Pi})

[seq(fsolve(f, q = OddAsymptotes[i] .. OddAsymptotes[i+1]), i = 1 .. n-1)];

[7.853797468, 10.99548650, 14.13711266, 17.27872097, 20.42032239, 23.56192056, 26.70351698, 29.84511237, 32.98670709]

Vector[row](n-1, proc (i) options operator, arrow; fsolve(f, q = OddAsymptotes[i] .. OddAsymptotes[i+1]) end proc)

Vector[row](9, {(1) = 7.853797468, (2) = 10.99548650, (3) = 14.13711266, (4) = 17.27872097, (5) = 20.42032239, (6) = 23.56192056, (7) = 26.70351698, (8) = 29.84511237, (9) = 32.98670709})

 

``

Download Asymptotes_acc.mw

Isn't this a repeat of one of your previous Questions?

For larger (symmetric) examples you will get a result more quickly by computing the eigenvalues of the Matrix (cast suitably to float[8] datatype, symmetric indexing function, etc).

You wrote in a followup comment, "In general, the matrices I encounter are real symmetric", so below I'll focus on code for that case.

Here's code for doing it via Eigenvalues, and an example where the SingularValues computation takes about 3.2 times as much wall-clock time, on my Linux machine. (That slowdown ratio could typically grow, as the Matrix size increases, but it can depend on the nature of the Graph data.)

SR := proc(G::Graph)
  local M := GraphTheory:-AdjacencyMatrix(G);
  max(abs~(LinearAlgebra:-Eigenvalues(
    rtable(`if`([rtable_indfns(M)]=[':-symmetric'],
           ':-symmetric',NULL), M,
           'datatype'=':-hfloat','storage'=':-rectangular',
           'order'=':-Fortran_order'))));
end proc:

M22:= GraphTheory:-SpecialGraphs:-M22Graph();

   M22 := Graph 1: an undirected unweighted graph with
          77 vertices and 616 edge(s)

CodeTools:-Usage(SR(M22), iterations=100);

memory used=270.28KiB, alloc change=38.59MiB, cpu time=1.84ms,
real time=1.86ms, gc time=110.19us

            16.0000000000000

G := GraphTheory:-SpecialGraphs:-CirculantGraph(2345,[1,2]);

    G := Graph 2: an undirected unweighted graph with
         2345 vertices and 4690  edge(s)

CodeTools:-Usage(SR(G));

memory used=242.85MiB, alloc change=107.52MiB, cpu time=5.75s,
real time=2.73s, gc time=312.03ms

            4.00000000000000

SpectralRadius_ac.mw

If you know certain additional qualities of the Graph up front (eg. Bipartite) then you might be able to attain even better speed by computing just the "first" selected eigenvalue. I had an old post on computing only selected eigenvectors/values. But it's not helpful if you have to call it twice to compute first and last (largest-positive and smallest-negative). I'm not sure whether this is workable.

These are all using external calling to LAPACK for dense matrix operations. ie. not sparse stuff.

Are you asking about how to get the exploration to appear in a new document, in a separate GUI tab? If so, then it can be done by using the Explore command explicitly, instead using rigght-click context-menu.

For example,

restart;

Explore(plot(sin(a*x)*x^2 + cos(b*x),
             x = -2*Pi .. 2*Pi, view = -30 .. 30),
        a = 1.0 .. 10.0, b = 1.0 .. 10.0,
        animate, newsheet, showbanner = false);

 

If you mean something else then you really ought to explain it more clearly. What do you mean by "save", and "extension file"?

That particular accent you've chosen may be getting parsed (in 2D Input) as conjugate(a). So you might be seeing a 2D Input version of this:

conjugate(a):=3;

3

conjugate(a)+2;

conjugate(a)+2

Download conjugate_sigh.mw

The same kind of thing seems to happen in the 2D Input parsing if I utilize the Layout palette (for "Over") and the Punctuation palette for inserting the "macr" over-bar. That is, I get a call to conjugate as the parsed meaning.

You could get the expected behaviour using the "Over" item from the Layout palette and then some other kind of over-bar from the Punctuation palette, eg. the "mdash", "minus", or "ndash" items. As over-bars those don't seem to make the names get parsed as calls to conjugate, instead parsing as distinct (typeset) names that behave more as you'd expected.

Please check whether you have shown the differential equation, and the value of lambda, correctly. (For example, what do you expect for 4*f(8) = x*f(eta) when x=4 and eta=8?)

As I parse the DE that you've shown, the value of 4*f(8) is only about 3.2 when lambda=0.5. But that doesn't correspond to the red color for the values near 4.0 in your colorbar. So either I've misunderstood the DE, or something else is not right.

Below I use lambda=0.005.

restart;

kernelopts(version);

`Maple 2021.2, X86 64 LINUX, Nov 23 2021, Build ID 1576349`

Typesetting:-Settings(typesetprime=true):
Typesetting:-Suppress(f(x));
interface(typesetting=extended):

de := diff(f(x),x,x,x) + (f(x)*diff(f(x),x,x) - diff(f(x),x)^2)
      - lambda*diff(f(x),x) = 0;

diff(diff(diff(f(x), x), x), x)+f(x)*(diff(diff(f(x), x), x))-(diff(f(x), x))^2-lambda*(diff(f(x), x)) = 0

lambda := 0.005:

desol := dsolve({de, f(0)=0, D(f)(0)=1, D(f)(10)=0}, f(x),
                numeric, output=listprocedure):

fsol := eval(f(x), desol):

4 * fsol(8);

HFloat(3.987786123641983)

plots:-display(
  plots:-densityplot(x*fsol(y), x=0..4, y=0..8,
                     style=surface, restricttoranges,
                     colorscheme=["zcoloring",
                                  z->2/3*(1-z/4)]),
  plots:-contourplot(x*fsol(y), x=0..4, y=0..8,
                     contours=[seq(0..4,0.5)], thickness=0),
  seq(plottools:-line([0,0],[0,0],thickness=10,
                      color=ColorTools:-Color("HSV",[2/3*(4-i)/(4-0),1,1]),
                      legend=sprintf("%.3f",i)),
      i=4..0,-0.5),
  labels=[x,eta], labeldirections=[horizontal,vertical],
  legendstyle=[location=right], axes=box,
  axis[1]=[tickmarks=[seq(0.0..4.0,0.5)]],
  axesfont=[Times,14], size=[500,400]
);

Download de_col.mw

 

[edit] Looking at the OP's previous postings it seems as if he might be using the significantly older Maple 17, in which case some of the above functionality is not available. Here is a revision, for that older version.

My query about the correctness of the DE and value of lambda still applies.

restart;

kerneopts(version);

kerneopts(version)

de := diff(f(x),x,x,x) + (f(x)*diff(f(x),x,x) - diff(f(x),x)^2)
      - lambda*diff(f(x),x) = 0;

diff(diff(diff(f(x), x), x), x)+f(x)*(diff(diff(f(x), x), x))-(diff(f(x), x))^2-lambda*(diff(f(x), x)) = 0

lambda := 0.005:

desol := dsolve({de, f(0)=0, D(f)(0)=1, D(f)(10)=0}, f(x),
                numeric, output=listprocedure):

fsol := eval(f(x), desol):

4 * fsol(8);

HFloat(3.9877861236419827)

plots:-display(
  plots:-densityplot(piecewise( x<-0.75, 0, x<-0.5, 1, x<0, 2/3,
                                2/3*(1 - x*fsol(y)/4) ),
                     x=-1..4, y=0..8,
                     style=patchnogrid,
                     colorstyle=HUE),
  plots:-contourplot(x*fsol(y), x=0..4, y=0..8,
                     contours=[seq(0..4,0.5)], thickness=0),
  seq(plottools:-line([0,0],[0,0],thickness=10,
                      color=ColorTools:-Color("HSV",[2/3*(4-i)/(4-0),1,1]),
                      legend=sprintf("%.3f",i)),
      i=4..0,-0.5),
  labels=[x,eta], labeldirections=[horizontal,vertical],
  legendstyle=[location=right], axes=box, view=[0..4,0..8],
  axis[1]=[tickmarks=[seq(0.0..4.0,0.5)]],
  axesfont=[Times,14]
);

Download de_col17.mw

You are substituting 4 for instances of 3.

subs(3=4,3);

             4

subs(3=4,9);

             9

x:=3;

             3

subs(x=4, 3);

             4

subs(x=4, 9);

             9

Check whether this (or any other Answer) provides as many roots as expected. (Eg, compare against the plot, or compare methods against each other, etc.)

restart:

with(Student:-Calculus1):

a[1] := .1093: k[3] := 7.5*10^(-12): k[2] := 3.8*10^(-12):
d := 0.2e-3: eta[1] := 0.240e-1: alpha[2] := -.1104:
alpha[3] := -0.1104e-2: eta[2] := .1361: xi := 1.219*10^(-6):
alpha := 1-alpha[3]^2/(a[1]*eta[1]):theta[0] := 0.5e-1:
Hc := (Pi/d)*sqrt(k[2]/xi): H := 5.5*Hc: lambda := a[1]/(xi*H^2):

 

sols := Roots((H/(Hc))^2 -(4*q^2)/Pi^2*((tan(q)- q/(1-alpha))/(tan(q)-q)),
              q=0..100, numeric, maxsols=200);

[1.579423766, 4.712622100, 7.853994064, 10.99555811, 14.13714635, 17.27873942, 20.42033357, 23.56192783, 26.70352198, 29.84511595, 32.98670974, 36.12830339, 39.26989691, 42.41149032, 45.55308364, 48.69467689, 51.83627007, 54.97786320, 58.11945627, 61.26104931, 64.40264231, 67.54423529, 70.68582823, 73.82742115, 76.96901405, 80.11060693, 83.25219980, 86.39379265, 89.53538549, 92.67697831, 95.81857112, 98.96016393]

nops(sols);

32

for i from 1 to nops(sols) do
 p[i] := alpha*lambda/(1-(4*sols[i])^2*(Hc/H)^2/Pi^2);
end do;

4.446263879

-20.25154621

-1.670508907

-.7029799724

-.3966600833

-.2567907775

-.1804394088

-.1339680365

-.1035031442

-0.8241813528e-1

-0.6720612510e-1

-0.5586412774e-1

-0.4717829446e-1

-0.4037722719e-1

-0.3495120384e-1

-0.3055223066e-1

-0.2693603757e-1

-0.2392701066e-1

-0.2139623682e-1

-0.1924733487e-1

-0.1740706641e-1

-0.1581896595e-1

-0.1443892890e-1

-0.1323209774e-1

-0.1217062731e-1

-0.1123205608e-1

-0.1039810281e-1

-0.9653765777e-2

-0.8986640837e-2

-0.8386399464e-2

-0.7844385207e-2

-0.7353298899e-2

 

Download rf_examp.mw

Take any odd number X greater than 1 (ie. X is 3 or higher).

Subtract 2 from it. So,

   X = 2 + (X-2)

The result, X-2, is also odd. And 2 is even.

So X = 2 + (X-2) which is the sum of an even number and an odd number.

In programming, a common easy way to generate odd numbers is to add 1 to a multiple of 2. (Similarly, a common way to generate even numbers is to take multiples of 2.) For example,

seq( 2*n, n=1..5 );

          2, 4, 6, 8, 10

seq( 2*n + 1, n=1..5 );

          3, 5, 7, 9, 11

We can use that trick to generate the following output, using the print command to show them on separate lines. (Recall, 2*n+1 is the odd number...)

seq( print( 2*n+1 =
            InertForm:-Display( 2 %+ (2*n+1-2), 'inert'=false ) ),
     n=1..5 ):

0, "%1 is not a command in the %2 package", _Hold, Typesetting

0, "%1 is not a command in the %2 package", _Hold, Typesetting

0, "%1 is not a command in the %2 package", _Hold, Typesetting

0, "%1 is not a command in the %2 package", _Hold, Typesetting

0, "%1 is not a command in the %2 package", _Hold, Typesetting

Download odd_stuff.mw

Something I have done in the past is redefine userinfo.

I once redefined it so as to redirect its output to TextArea Embedded Components. But in your case perhaps it might suffice to alter it to always act as if called from the Command-Line Interface (CLI, tty), ie. plaintext output.

I may have overlooked some corner-case(s), in which case it might require a tweak or two.

restart;

 

unprotect(userinfo):
userinfo:=proc(lev::nonnegint,
               nms::{name,set(name)})
  local ee,L,snms;
  if nms::name then snms:={nms}
  else snms:=nms; end if;
  if nargs<3 then return NULL; end if;
  L:=[args[`if`(args[3]=':-NoName',4,3)..-1]];
  if nops(L)=0 then return NULL; end if;
  if ormap(nm->assigned(infolevel[nm])
               and infolevel[nm]>lev,
           snms) then
    if args[3]<>':-NoName' then
      printf("%s",sprintf("%s: ",
                          debugopts(':-callstack')[5]));
    end if;
    printf("%s",convert(L[1],string));
    for ee in L[2..] do
      printf(" %s",convert(ee,string));
    end do;
    printf("\n");
  end if;
  NULL:
end proc:
protect(userinfo):

 

infolevel[:-dsolve]:=5:

dsolve(diff(y(x),x)=sin(x),y(x));

Methods for first order ODEs:
--- Trying classification methods ---
trying a quadrature
<- quadrature successful

y(x) = -cos(x)+_C1

writeto("output_of_dsolve.txt");  #send all output to file
sol:=dsolve(diff(y(x),x)=sin(x),y(x));
close("output_of_dsolve.txt"):
writeto(terminal); #to send output back to terminal

 

Download userinfo_redef.mw

You wrote, "Is it possible to recover that _R is in fact A, that _R0 is in fact B and so on?".

But question includes the assertion that "_R is in fact A", which may not be true. If you had also executed,
   P := A;
at the top level then both P and A would be assigned with value _R. So why would _R "in fact" be A, as opposed to being P? The claim, as stated, would not be true.

Having said that, you can cobble together (dodgy) schemes to find out what assigned names might evaluated to _R.

restart:

with(Statistics):

A := RandomVariable(Uniform(0, 1));

_R

B := RandomVariable(Uniform(0, 2));

_R0

C := RandomVariable(Uniform(1, 2));

_R1

F := (A+2*B)/C;

(_R+2*_R0)/_R1

recover := proc(ee::algebraic) local S;
  S := {anames(user)};
  subsindets(ee, name,
             proc(nm) local c;
               c := select(v->eval(v)=nm,S);
               if nops(c)=0 then nm; else c[1]; end if;
            end proc);
end proc:

 

recover(F);

(A+2*B)/C

Download ooph.mw

I am not a fan of this kind of thing, which seems quite backwards to me. I think it's far more sensible to utilize lists of equations and then evaluate formulas using such equations. (Ie, don't make assignments that get in your way.) For example,

restart:

with(Statistics):

a := RandomVariable(Uniform(0, 1)):

b := RandomVariable(Uniform(0, 2)):

c := RandomVariable(Uniform(1, 2)):

F := (A+2*B)/C;

(A+2*B)/C

eval(F, [A=a, B=b, C=c]);

(_R+2*_R0)/_R1

F;

(A+2*B)/C

Download form_forwards.mw

First 72 73 74 75 76 77 78 Last Page 74 of 336