acer

32373 Reputation

29 Badges

19 years, 333 days
Ontario, Canada

Social Networks and Content at Maplesoft.com

MaplePrimes Activity


These are answers submitted by acer

Here's one way. I tried to focus on making it easier to understand, instead of using fancy terse syntax or making it super efficient.

restart;

mylist := [x^4 + (-4*m - 7)*x^3 + (m + 4)*x^2 + (3*m - 5)*x, x^4 + (-4*m - 7)*x^3 + (m + 5)*x^2 + (5*m - 7)*x, x^4 + (-4*m - 7)*x^3 + (m + 5)*x^2 + (7*m - 5)*x, x^4 + (-4*m - 7)*x^3 + (2*m + 5)*x^2 + (3*m - 5)*x, x^4 + (-4*m - 7)*x^3 + (3*m + 1)*x^2 + (7*m - 10)*x, x^4 + (-4*m - 5)*x^3 + (2*m + 1)*x^2 + (7*m - 9)*x]:

L := map~(normal, mylist):

ans := [seq([ L[i], diff(L[i],x),
              solve([ eval(diff(L[i],x),x=2),
                      eval(diff(L[i],x,x),x=2)>0 ], m) ],
            i = 1 .. nops(L))]:

print~(ans):

[x^4-(4*m+7)*x^3+(m+4)*x^2+(3*m-5)*x, 4*x^3-3*(4*m+7)*x^2+2*(m+4)*x+3*m-5, {m = -1}]

[x^4-(4*m+7)*x^3+(m+5)*x^2+(5*m-7)*x, 4*x^3-3*(4*m+7)*x^2+2*(m+5)*x+5*m-7, {m = -1}]

[x^4-(4*m+7)*x^3+(m+5)*x^2+(7*m-5)*x, 4*x^3-3*(4*m+7)*x^2+2*(m+5)*x+7*m-5, {m = -1}]

[x^4-(4*m+7)*x^3+(2*m+5)*x^2+(3*m-5)*x, 4*x^3-3*(4*m+7)*x^2+2*(2*m+5)*x+3*m-5, {m = -1}]

[x^4-(4*m+7)*x^3+(3*m+1)*x^2+(7*m-10)*x, 4*x^3-3*(4*m+7)*x^2+2*(3*m+1)*x+7*m-10, {m = -2}]

[x^4-(4*m+5)*x^3+(2*m+1)*x^2+(7*m-9)*x, 4*x^3-3*(4*m+5)*x^2+2*(2*m+1)*x+7*m-9, {m = -1}]

# For each `ee` in `ans`,
#   diff(ee[1],x)=ee[2] is 0 at x=2,
# and
#   diff(ee[1],x,x)=diff(ee[2],x) is positive at x=2

seq(eval([ee[2],diff(ee[2],x)],[x=2,ee[3][]]), ee=ans);

[0, 18], [0, 20], [0, 20], [0, 18], [0, 50], [0, 34]

Download minhthien_s1.mw

[edit] This is also not efficient. But, since you asked, it uses the prime notation for differentiation in 2D Input, and applies at argument 2 (instead of evaluating at point x=2).
minhthien_s2.mw

I did not look for any pattern in the entries of mylist (as a possible shortcut to the solutions, say). I don't know whether you have any rules for forming those entries.

Here's one way to get a more accurate plot. See attachment.

This is an interesting plot. (Is it meaningful, in terms of the as-yet unstated background theory here?)

Reza_Tchno_v3_ac.mw

Have you investigated the system symbolically, possibly using exact rationals? That might possible be able to demonstrate the feasible bounds for alpha, or piecewise explicit formulae, etc.

Do you know whether there is only ever a single solution that satisfies all the constraints?

Do you know for sure that it wouldn't be possible to churn out a relatively simple piecewise solution (in terms of alpha, say), that might be super fast to compute/plot?

Are you interested in finding the critical alpha values?

Have you considered trying to use purely floating-point optimization (eg. alpha->1 as constant objective, and the rest as mixed equality/inequating constraints, thus leveraging Optimization:-Minimize to compute/find a feasible point), for better speed?

Here's one way to handle your given polynomial example.

expr := x^3 + (-m - 1)*x^2 + (-4*m - 3)*x;

x^3+(-m-1)*x^2+(-4*m-3)*x

map(s->sign(s)*normal(s*sign(s)), expr);

x^3-(m+1)*x^2-(4*m+3)*x

Download sgn_ex.mw

ps. You gave the target form as x^3 - (m-1) x^2 - (4m + 3) which isn't equivalent to your input expression. And so I have guessed as to what you wanted to accomplish with the coefficients.

restart

ts := [csc = (proc (u) options operator, arrow; 1/sin(u) end proc), sec = (proc (u) options operator, arrow; 1/cos(u) end proc), cot = (proc (u) options operator, arrow; 1/tan(u) end proc)]


Expression to apply an identity to:

orig := JacobiSN(sin((1/2)*`ϕ__0`)*t, csc((1/2)*`ϕ__0`))

JacobiSN(sin((1/2)*varphi__0)*t, csc((1/2)*varphi__0))

(1)

We don't need to solve anything at this stage. We may simply assert that z,k equal the
arguments of JacobiSN in orig.
We can then quite reasonably utilize the general form Jacobi(z,k) later below.

C := eval({k = op(2, orig), z = op(1, orig)}, ts)

{k = 1/sin((1/2)*varphi__0), z = sin((1/2)*varphi__0)*t}

(2)

Next, without having to cook up 1/k in Jacobi(z,1/k) ourselves
(which would make a somewhat ad hoc approach), which happens
to arise from the LHS of the 5th identity.

This obtains a formula with general Jacobi(z,k) on the LHS.

temp := convert(FunctionAdvisor(identities, JacobiSN(__A, __B))[5], `global`); solve([z = op(1, lhs(temp)), k = op(2, lhs(temp))], [__A, __B])[1]; eqn := subs(%, temp)

JacobiSN(z, k) = JacobiSN(z*k, 1/k)/k

(3)

Easy way,

orig = subsindets(orig, specfunc(JacobiSN), proc (u) options operator, arrow; eval(rhs(eqn), C) end proc)

JacobiSN(sin((1/2)*varphi__0)*t, csc((1/2)*varphi__0)) = sin((1/2)*varphi__0)*JacobiSN(t, sin((1/2)*varphi__0))

(4)

Alternative using applyrule:
Using a dummy name like __d prevents runaway recursion under applyrule.

R := subsindets(lhs(eqn), And(name, Not(constant)), proc (nm) options operator, arrow; nm::anything end proc) = subs(JacobiSN = __d, rhs(eqn))

JacobiSN(z::anything, k::anything) = __d(z*k, 1/k)/k

(5)

orig = subs(__d = JacobiSN, C, applyrule(R, lhs(eqn)))

JacobiSN(sin((1/2)*varphi__0)*t, csc((1/2)*varphi__0)) = sin((1/2)*varphi__0)*JacobiSN(t, sin((1/2)*varphi__0))

(6)

restart


Now  a new example: Converting InverseJacobinAM to InverseJacobiSN

orig := InverseJacobiAM((1/2)*`ϕ__0`, sqrt(2)/sqrt(1-cos(`ϕ__0`)))

InverseJacobiAM((1/2)*varphi__0, 2^(1/2)/(1-cos(varphi__0))^(1/2))

(7)

eqn := convert(FunctionAdvisor(identities, InverseJacobiSN(z, k))[3], `global`)

InverseJacobiSN(z, k) = InverseJacobiAM(arcsin(z), k)

(8)

Use eval(eqn,1) for 1-level evaluation (to prevent its RHS from evaluating to an EllipticF call).

C := solve(`~`[`=`]([op(orig)], [op(rhs(eval(eqn, 1)))]), [k, z])[1]

[k = 2^(1/2)/(1-cos(varphi__0))^(1/2), z = sin((1/2)*varphi__0)]

(9)

orig = subs(C, lhs(eqn))

InverseJacobiAM((1/2)*varphi__0, 2^(1/2)/(1-cos(varphi__0))^(1/2)) = InverseJacobiSN(sin((1/2)*varphi__0), 2^(1/2)/(1-cos(varphi__0))^(1/2))

(10)

NULL

Download Applying_identities_from_FunctionAdivisor_ac.mw

restart;

with(InertForm):

 

The following automatic simplification occurs even before evaluation (which is why it cannot be prevented by single right-tick unevaluation quotes).

 

'3*(a-1)';

3*a-3

expr := 3*(a-1);

3*a-3

 

We can separate the numeric factor, and then use %* for inert multiplication.

 

f := content(expr) %* primpart(expr);

`%*`(3, a-1)

Display(f, inert=false);

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

value(f);

3*a-3

Download auto_simp_ex.mw

note: A little bit more on automatic simplification can be had in Section 3.4 "Unevaluated Expressions" of Chapter 3 of the Programming Guide, in the subsection "Evaluation and Automatic Simplification".

This question pops up now and then. This old Reply contains a couple of other approaches, including using MakeInert instead of %*.

Another alternative is to wrap one of the factors in a call to ``(...). For example,
   content(expr) %* ``(primpart(expr))
as shown also in that very same linked old Reply of mine. I don't generally prefer that approach over using InertForm. I also dislike how its results can sometimes get display with extra brackets in output, even though that doesn't occur for this particular example. I generally prefer how InertForm works, over the ``(...) approach.

With so many floating-point evaluations of those special functions the internal memory management system (gc) may collect some of the computed results from evalf's remember table.

But your expressions have two pairs of repeated special-function calls, ie. to WeierstrassZeta and to WeierstrassPPrime.

By forcing the computed values to be remembered (even when gc occurs) I get about another factor of two speed up here.

restart;

plots:-setoptions(size=[350,350]);

g := evalf((Beta(1/4, 1/4)/2)^4):

WZ := proc(a,b,c) option remember;
  if not [a,b,c]::list(complexcons) then return 'procname'(args);
  else evalf(WeierstrassZeta(a,b,c)); end if;
end proc:

WPP := proc(a,b,c) option remember;
  if not [a,b,c]::list(complexcons) then return 'procname'(args);
  else evalf(WeierstrassPPrime(a,b,c)); end if;
end proc:

AA := Pi*r*cos(theta) - Re(WZ(r*exp(theta*I), g, 0) + Pi*WPP(r*exp(theta*I), g, 0)/g);
BB := Pi*r*sin(theta) + Im(WZ(r*exp(theta*I), g, 0) - Pi*WPP(r*exp(theta*I), g, 0)/g);
C := sqrt(6*Pi/g)*Re(WeierstrassP(r*exp(theta*I), g, 0));

Pi*r*cos(theta)-Re(WZ(r*exp(I*theta), 189.0727201, 0)+0.1661579022e-1*WPP(r*exp(I*theta), 189.0727201, 0))

Pi*r*sin(theta)+Im(WZ(r*exp(I*theta), 189.0727201, 0)-0.1661579022e-1*WPP(r*exp(I*theta), 189.0727201, 0))

.3157447408*Re(WeierstrassP(r*exp(I*theta), 189.0727201, 0))

CodeTools:-Usage(plot3d([AA, BB, C], r = 1/5 .. 4/5, theta = -Pi .. Pi,
                        view = [-8 .. 8, -8 .. 8, -8 .. 8], shading = zhue,
                        style = surface, grid = [33, 129]));

memory used=3.61GiB, alloc change=142.00MiB, cpu time=33.23s, real time=30.42s, gc time=4.48s

 

Download WZ.mw

An alternative way to speed it up would be to parallelize the population of the float[8] Array inside the resulting MESH plotting structure. [see below]

Here is the most that I am able to recover.

Basismat_2_noter_acc.mw

I suspect that it is missing just two inputs, which I'd guess are 2D Math Matrices that occurred here:
1) just after "Løsning er da:"
2) just before "Bytter række 2 og 3:"

I suggest that you try to use Maple 2023 instead of your current Maple 2022. This kind of worksheet corruption appears (so far) to be much less common in 2023.

You can supply an option to the call to the plot command, for specifying the relative location of the first axis.

And you can change the size option below, to get dimensions more to your liking.

restart

k := 0.1e-1:

Pinit := 100:

V := 25:

rho := 1:

eq := diff(P(t), t) = -k*P(t)+rho*Cp*(int(V*sin(alpha), t));

diff(P(t), t) = -0.1e-1*P(t)+12.5*sin(5)*t

solution := dsolve({eq, P(0) = Pinit}, P(t));

P(t) = 1250*sin(5)*t-125000*sin(5)+exp(-(1/100)*t)*(100+125000*sin(5))

plot(rhs(solution), t = 0 .. 600, axis[1] = [location = low], size = [500, 400])

Download time-dependent_aerodynamic_pressure_ac.mw

See the Help page for topic  plot,options  for more details on 2D plotting choices.

The following reuable procedure P produces an expression (which is intended as pretty-printing in Maple's GUI similarly as does the the unevaluated ContinuedFraction function call).

Are you able to utilize the following expression assigned to res (either by copying its pretty-printed 2D output as special, rich-text, etc, or as generated MathML?)

restart

f := NumberTheory:-ContinuedFraction((3+sqrt(7))*(1/2)); lprint(%)

_m139948790579680

NumberTheory:-ContinuedFraction(3/2+1/2*7^(1/2))

P := proc (ee, N::posint) local f, H, i, T; H := Term(NumberTheory:-ContinuedFraction(ee), 0 .. N-1); T := `..`; for i from N by -1 to 2 do if (H[i])::list then T := H[i, 1]/(H[i, 2]+T) else T := 1/(H[i]+T) end if end do; H[1]+T end proc

res := P((3 + sqrt(7))/2, 10);

2+1/(1+1/(4+1/(1+1/(1+1/(1+1/(4+1/(1+1/(1+1/(1+`..`)))))))))

lprint(MathML:-ExportPresentation(res));

"<math xmlns='http://www.w3.org/1998/Math/MathML'><mrow><mn>2</mn><mo>+</mo><mf\
rac><mn>1</mn><mrow><mn>1</mn><mo>+</mo><mfrac><mn>1</mn><mrow><mn>4</mn><mo>+<\
/mo><mfrac><mn>1</mn><mrow><mn>1</mn><mo>+</mo><mfrac><mn>1</mn><mrow><mn>1</mn\
><mo>+</mo><mfrac><mn>1</mn><mrow><mn>1</mn><mo>+</mo><mfrac><mn>1</mn><mrow><m\
n>4</mn><mo>+</mo><mfrac><mn>1</mn><mrow><mn>1</mn><mo>+</mo><mfrac><mn>1</mn><\
mrow><mn>1</mn><mo>+</mo><mfrac><mn>1</mn><mrow><mn>1</mn><mo>+</mo><mi>..</mi>\
</mrow></mfrac></mrow></mfrac></mrow></mfrac></mrow></mfrac></mrow></mfrac></mr\
ow></mfrac></mrow></mfrac></mrow></mfrac></mrow></mfrac></mrow></math>"

#res := P(cos(x), 12);

#ContinuedFraction(cos(x), x);

#MathML:-ExportPresentation(res);

 

Download ContFrac_MathML.mw

Your worksheet has some coding mistakes.

1) The definition of prof1 has some function calls where multiplcation was likely intended. Eg,
   (...)(ps - s)
is is a function call, and not a product of two brackets terms. For multiplication you'd need either an explicit multiplication symbol, or (in 2D Input) a space, between the adjacent closing&opening brackets )( .  This occurs a few times.

2) The attempted multiple assignments like,
   [p,q] := solve(...)
is invalid. This occurs a few times.

It's often useful to test your procedure with some actual float arguments, before trying to plot it. That would reveal to you that WhyNot3 wasn't working at all.

I've corrected those aspects, in this attachment.  HmaxProc_ac.mw

Here's another approach.

If you'd prefer the result as just a list of two numbers then you could leave off the [a,b]=~ bit at the start.

L := [a,b]=~applyrule(a::numeric*f^b::numeric=[a,-b],X1);

          L := [a = 1.234, b = 6.789]

A couple of other examples,

[a,b]=~applyrule(a::numeric*f^b::numeric=[a,-b],-1/(f^6.789));
                      [a = -1, b = 6.789]

[a,b]=~applyrule(a::numeric*f^b::numeric=[a,-b],f^6.789);
                      [a = 1, b = -6.789]

Here's your original example, round-tripped,

X := a/f^b;

a/f^b

X1 := 1.234/(f^6.789);

1.234/f^6.789

L := [a,b]=~applyrule(a::numeric*f^b::numeric=[a,-b],X1);

[a = 1.234, b = 6.789]

eval(X, L);

1.234/f^6.789

Download rt.mw

ps. The order of terms in the result is forced by us, here, since the order of terms in the returned list is specified in the code.

pps. You could also tweak the specified types, eg. only looking for negative b, etc. You could also wrap it up in a reusable procedure. You could also use indets (and a structured type) to yank out multiple such subexpressions from a larger expression. Hard to say more because we don't yet know exactly how you intend to utilize this.

@Samuel Hesselberg Any new document will have its own Startup code region (which is empty, at first). The Startup code is not automatically inherited from any other already opened sheet.

Your original .mw file is a document with its own, separate Startup code.

So if you copy&paste the main body of that older document then you'll also have to provide access to the code in the older document's Startup region as well. The embedded components you're pasting in rely upon and call that code.

You can add it to the Startup code region of your newly opened document, or you can add in somewhere else in your new document's body (any document block, or even a code-edit-region).

First, find the Startup code in your original. Then figure out how you want your new document to be able to also use it. The easiest way is also to copy the older document's Startup code to the new document's Startup code region.

I am not sure about a good approach for Maple 2015.

In Maple 2018 and later you could use plots:-textplot with its rotation option.

An example is given below.

Note that I did not have the lines extend down to y=0, since your image didn't either. I also didn't make it handle the lines if f(x) becomes negative. I figure that if you want such a more general example or different behaviour you would provide it -- and all your requirements -- explicitly and clearly.

restart;

with(plots):

f := x -> exp(cos(x)+sin(x));

proc (x) options operator, arrow; exp(cos(x)+sin(x)) end proc

(a,b) := 0,3:

n := 12:

h := (b-a)/n:

L := [seq([x,evalf(f(x))],x=a..b,h)]:

plots:-display(
  plot(f,a..b),
  pointplot(L[2..-2],symbolsize=10,symbol=solidcircle),
  seq(textplot([L[i][],
                typeset(" ",eval('Typesetting:-Typeset'('f'(evalf[3](L[i][1]))))
                            =evalf[5](L[i][2])),
                rotation=Pi/2],align=above,font=[DEFAULT,10]),i=2..nops(L)-1),
  seq(plottools:-line([L[i][1],min(L[..,2])],L[i],thickness=0),i=1..nops(L)-1),
  tickmarks=[L[..,1],[$ceil(min(L[..,2]))..floor(max(L[..,2]))]],
  font=[DEFAULT,9], size=[700,450]);

Download rotated_textplot_ex2.mw

Are you asking about codegen:-optimize? You haven't usefully provided an example of what exactly you're trying to accomplish.

If it's replacing generated names in results from codegen:-optimize then it's not clear to me whether you're dealing with Maple expressions of procedures. Also, while codegen:-optimize does have a documented mechanism to utilize a given list of preferred names, it's somewhat awkward IMO.

Perhaps something like the following might serve:

restart;

 

rep := (L::list(`=`),s::symbol)
  -> subs((nm->nm=parse(s||(String(nm)[2..])))~(lhs~(L)),
          L)[]:

 

repproc := (P::procedure,s::symbol)
  -> subs((nm->nm=parse(s||(String(nm)[2..])))~([op(2,eval(P))]),
          eval(P)):


We use rep to replace the temporary names in the optimized computation sequence.

f := x*ln(x)+2*x^2*ln(x)+3*x^3*ln(x);

x*ln(x)+2*x^2*ln(x)+3*x^3*ln(x)

temp := codegen:-optimize(f);

t1 = ln(x), t3 = x^2, t9 = 3*t1*t3*x+2*t1*t3+t1*x

rep([temp], b);

b1 = ln(x), b3 = x^2, b9 = 3*b1*b3*x+2*b1*b3+b1*x

 

We use repproc to replace the locals in the optimized procedure.

 

p := proc(x) x*ln(x)+2*x^2*ln(x)+3*x^3*ln(x); end proc;

proc (x) x*ln(x)+2*x^2*ln(x)+3*x^3*ln(x) end proc

temp2 := codegen:-optimize(p);

proc (x) local t1, t3; t1 := ln(x); t3 := x^2; 3*t1*t3*x+2*t3*t1+t1*x end proc

repproc(temp2, b);

proc (x) local b1, b3; b1 := ln(x); b3 := x^2; 3*b1*b3*x+2*b3*b1+b1*x end proc

Download cg_rep.mw

You have p(lambda__1 + p) which is a function call (to p), and not a multiplication. Perhaps you intended that as p*(lambda__1 + p) instead.

restart

psi = Pi^2*n^2+p; q := diff(theta(y), y, y) = psi^2*theta(y)

bc := theta(0) = 0, theta(1) = `&lambda;__1`*(cos(n*Pi)-1)/(p*(`&lambda;__1`+p)*n)

q := dsolve({bc, q})

theta(y) = lambda__1*(cos(n*Pi)-1)*exp(psi*y)/(n*p*(exp(psi)*p+exp(psi)*lambda__1-p*exp(-psi)-lambda__1*exp(-psi)))-lambda__1*(cos(n*Pi)-1)*exp(-psi*y)/(n*p*(exp(psi)*p+exp(psi)*lambda__1-p*exp(-psi)-lambda__1*exp(-psi)))

simplify(convert(q, trigh))

theta(y) = lambda__1*(cos(n*Pi)-1)*sinh(psi*y)/(p*(lambda__1+p)*n*sinh(psi))

Download nima_ac.mw

First 38 39 40 41 42 43 44 Last Page 40 of 336