sursumCorda

997 Reputation

13 Badges

1 years, 199 days

MaplePrimes Activity


These are answers submitted by sursumCorda

Here is a somewhat inefficient approach: 

e1 := a/sqrt(tan(x+c__1)^2+1):
(_ -> simplify(ifelse(type(_, atomic), _, map(thisproc, _)), trig))(e1);
 = 
                     

But the advantage is that there is no need to understand the structure of the expression beforehand.

There exist at least two approaches, 

e1 := a/sqrt(tan(x + c__1)**2 + 1):
convert(e1, exp, 'simplifier' = evala@rcurry(simplify, trig, exp)@normal);
 = 
                     convert(simplify(convert(a/sqrt(tan(x + c__1)^2 + 1), 'sincos'), sqrt, constant), 'sectan')

simplify(simplify(convert(e1, 'sincos'), constant, sqrt), trig);
 = 
                     convert(convert(a/sqrt(tan(x + c__1)^2 + 1), exp, 'simplifier' = evala@rcurry(simplify, trig, exp)@normal), sec)

even though neither of them looks elegant …. (By the way, I believe it should be a BUG that the output of simplify(e1, trig); remains unchanged.) 

If you set interface('typesetting' = "standard"):, the eqn will look somewhat shorter: 

"A& B& \\textrm{The fundamental matrix has } &C ":
StringTools:-RegSubs("(.*)\\\\textrm\\{(.*)\\}(.*)" = "\\1\\begin{minipage}{\\linewidth}\\textrm{\\2}\\end{minipage}\\3", `%`);
 = 
  "A& B& \begin{minipage}{\linewidth}\textrm{The fundamental 

     matrix has }\end{minipage} &C "


StringTools:-RegSubs("(.*)\\\\textrm\\{(.*)\\}(.*)" = "\\1\\begin{minipage}{\\linewidth}\\textrm{\\2}\\end{minipage}\\3", "\\textrm{ WILDCARD\n* }");
 = 
        "\begin{minipage}{\linewidth}\textrm{ WILDCARD

          * }\end{minipage}"


 Unfortunately, Maple's regular expression engine does not offer certain features, and therefore there exists a problem: 

StringTools:-RegMatch("(.*)\\\\textrm\\{(.*)\\}(.*)", " A&B&{\\textrm{The fundamental matrix has {}}}&C ", '_0', '_1', '_2', '_3'):
StringTools:-IsBalanced~([_ || (0 .. 3)], "{", "}");
 = 
                   [true, false, false, true]

In my view, a better approach is using Python's `regex` library (in Maple) instead. 

Although inefficient (much slower than the direct nested loops), this is still a good idea: 

s := (n::posint) -> 
   add(add(1 /~ 
      mul~(ListTools:-FlattenOnce(
        combinat:-permute~(
         combstruct['allstructs']('Partition'(l), 
          'size' = 3))))), l = 3 .. n):
Warning, (in s) `l` is implicitly declared local
seq(s(n), n = 1 .. 10);
 = 
                    5  17  49  967  801  4523  84095
           0, 0, 1, -, --, --, ---, ---, ----, -----
                    2  4   8   120  80   378   6048 

time[real](s(100));
 = 
                             6.046


Why not just 

[LinearAlgebra:-Column](Matrix(3,4,'symbol'=m),[1..-1]);
         [LinearAlgebra:-Column](Matrix(3,4,'symbol'=m),[1..-1])


whattype~(`%`);
[Vector[column], Vector[column], Vector[column], Vector[column]]


I'm not sure why, but there're at least 2 workarounds: 

RealDomain:-solve('identity'(eq,x),trial_solution_constants);
PDETools:-Solve(eq,trial_solution_constants,'independentof'=x);

I am unclear about why the latter command is so inefficient. (Actually, Mma carried out the same thing without any manual interference within 2 seconds!) 
Anyway, as regards your example, certain manual interventions appear necessary. 

restart;
B_poly := -(x + 4)*(x + 3)*(x + 2)*(x + 1)*(x - 1)*(x - 2)*(x - 3)*(x - 4):
g := B_poly/100000*(((B_poly - 669)/670)^136 - ((B_poly + 669)/670)^136):
f := -4347225/87808*x^8 - 17375/392*x^7 + 629491375/395136*x^6 + 266375/252*x^5 - 200677775/12544*x^4 - 3174625/504*x^3 + 11067842125/197568*x^2 - 53625/98*x - 126496075/4116:

Firstly, we can see that we only need test four (closed) intervals. 

rts := applyrule('RealRange(a::anything, b::anything) = a .. b', [RealDomain:-solve](B_poly >= 0, x));
 = 
              [-4 .. -3, -2 .. -1, 1 .. 2, 3 .. 4]

Next, we shall check whether g - f changes signs on these four intervals. (In theory, using sturm(sturmseq(g - f, x), x, a, b) should be enough; unfortunately, this built-in function seems to be less efficient as well (and so is fsolve(g - f, x, a .. b, 'maxsols' = 1)).) Here is a workaround. 

CodeTools:-Usage([seq](traperror(RootFinding:-RefineRoot(k, expand(g - f), x)), k = rts)); # treated as open intervals 
 = 
memory used=11.05MiB, alloc change=0 bytes, cpu time=515.00ms, real time=520.00ms, gc time=0ns

        ["given interval does not contain any roots", 

          "given interval does not contain any roots", 

          "given interval does not contain any roots", 

          "given interval does not contain any roots"]


Now we know g - f is either positive or negative in each of (-4, -3), (-2, -1), (1, 2), and (3, 4), so we just need check if g - f >= 0 at the endpoints and some interior point of each interval. However, since the inequality is not strict, in accordance with the continuity, it is enough to consider one interior point of each interval. 

andseq(eval(g - f >= 0, x = stats['describe', 'mean'](convert(k, list))), k = rts); # verify at the midpoint 
 = 
                             false

Accordingly, the region {B_poly ≥ 0, g - f ≥ 0} is empty. 

These computations can be finished in one second on my low-end computer. 

First and foremost, “:=” should be used instead: 

cont_real := Re(continuous_solution):
cont_imag := Im(continuous_solution):
disc_real := Re(discrete_solution):
disc_imag := Im(discrete_solution):

To display a list of 2-D points, “plots:-pointplot” should be used: 

p1 := plots:-pointplot(x_values, cont_real, 'color' = "blue", 'legend' = "Continuous - Real"):
p2 := plots:-pointplot(x_values, disc_real, 'color' = "green", 'legend' = "Discrete - Real"):
p3 := plots:-pointplot(x_values, cont_imag, 'color' = "red", 'legend' = "Continuous - Imag"):
p4 := plots:-pointplot(x_values, disc_imag, 'color' = "orange", 'legend' = "Discrete - Imag"):

However, since some data in `disc_real` and `disc_imag` are too large, it is better to draw "log10~(disc_real)" and "log10~(disc_imag)" instead.

How about

is_symbol_inside_func_only := proc(expr::anything, f::name, y::symbol, ` $`)::truefalse; 
	type(applyrule('conditional'(f(_x::anything), _depends(_x, y)) = `tools/gensym`(f), expr), freeof(y)) 
end: 
expr:=3*ln(1+y)+ln(3*y)*y+ln(y)+cos(7*y):
is_symbol_inside_func_only(expr,ln,y); #should return false

expr:=3*ln(1+y)+ln(3*y):
is_symbol_inside_func_only(expr,ln,y); #should return true

expr:=ln(y)+ln(3*y)+cos(y):
is_symbol_inside_func_only(expr,ln,y); #should return false

expr:=3+cos(y):
is_symbol_inside_func_only(expr,cos,y); #should return true

expr:=y+ln(y):
is_symbol_inside_func_only(expr,ln,y); #should return false
 = 
                             false

                              true

                             false

                              true

                             false


Sometimes it would be better if there exists a special cbrt function. However, here you can use surd(x, 3) directly
Edit. Note that 

showstat(RealDomain:-`^`, 3);

RealDomain:-`^` := proc()
       ...
   3   clean(`assuming`([proc( s, n )
           if type(n,'fraction') then
               surd(s,denom(n))^numer(n);
           else
               s^n;
           end if;
       end proc(_passed)],['real']))
end proc


So it is somewhat unnecessary to with(RealDomain, `^`):.

Alternatively, 

r1 := -1 <= x and x <= 0:
r2 := 0 <= x and x <= 1:
convert(solve(r1 or r2, {x}), `and`);
 = 
                       -1 <= x and x <= 1

 

Some help pages provide a rough indication of basic methods and algorithms used for a particular purpose inside. Sometimes setting infolevel[…] to a positive integer (or just using showstat(…)) is also conducive to seeing such information, but it may still be difficult to determine which formula was used simply by reading the source code.
For instance, according to 

showstat((combinat::Partitions)::NumberOfPartitions, 13 .. 14):

(combinat::Partitions):-NumberOfPartitions := proc(n, {method::identical(recursive,dp,hrr,auto) := 'auto'})
local pn;
       ...
  13     userinfo(1,':-numbpart',"Using Hardy-Ramanujan-Rademacher series algorithm");
  14     pn := combinat:-Partitions:-NumberOfPartitions:-numbpart1_hrr(n)
       ...
end proc

, Maple by default make use of the so-called Hardy–Ramanujan–Rademacher method for larger n as well. However, the artificial benchmark - partition function claims that

all implementations (except Maple?) use the numerical Hardy-Ramanujan-Rademacher formula.
My questions are then, is there a way to get Maple to solve this? 

Yes, there is a way: 
 

restart

Warning, inserted missing semicolon at end of statement

 

simplify(eval(student['changevar'](u = sqrt(sqrt(x^4+1)+x^2), int(sqrt(sqrt(x^4+1)+x^2)/((x+1)*sqrt(x^4+1)), x), u), u = sqrt(sqrt(x^4+1)+x^2)), 'size' = true, assume = negative)
NULL

(1/4)*2^(1/2)*(2*(2^(1/2)-1)^(1/2)*arctan(((x^4+1)^(1/2)+x^2)^(1/2)/(2^(1/2)-1)^(1/2))+(2^(1/2)-1)^(1/2)*arctan((1/2)*((2^(1/2)-1)*(x^4+1)^(1/2)+2^(1/2)*x^2-x^2+1)/((2^(1/2)-1)^(1/2)*((x^4+1)^(1/2)+x^2)^(1/2)*x))-2*(1+2^(1/2))^(1/2)*arctanh(((x^4+1)^(1/2)+x^2)^(1/2)/(1+2^(1/2))^(1/2))+(1+2^(1/2))^(1/2)*arctanh((1/2)*(2^(1/2)*x^2+(x^4+1)^(1/2)*2^(1/2)-1+x^2+(x^4+1)^(1/2))/((1+2^(1/2))^(1/2)*x*((x^4+1)^(1/2)+x^2)^(1/2))))

(1)

radnormal(diff((1/4)*2^(1/2)*(2*(2^(1/2)-1)^(1/2)*arctan(((x^4+1)^(1/2)+x^2)^(1/2)/(2^(1/2)-1)^(1/2))+(2^(1/2)-1)^(1/2)*arctan((1/2)*((2^(1/2)-1)*(x^4+1)^(1/2)+2^(1/2)*x^2-x^2+1)/((2^(1/2)-1)^(1/2)*((x^4+1)^(1/2)+x^2)^(1/2)*x))-2*(1+2^(1/2))^(1/2)*arctanh(((x^4+1)^(1/2)+x^2)^(1/2)/(1+2^(1/2))^(1/2))+(1+2^(1/2))^(1/2)*arctanh((1/2)*(2^(1/2)*x^2+(x^4+1)^(1/2)*2^(1/2)-1+x^2+(x^4+1)^(1/2))/((1+2^(1/2))^(1/2)*x*((x^4+1)^(1/2)+x^2)^(1/2)))), x)-sqrt(sqrt(x^4+1)+x^2)/((x+1)*sqrt(x^4+1)));

0

(2)

int((-7+x)/((-11+5*x)*sqrt(x^4-3*x^3-21*x^2+83*x-60)), x, 'method' = 'Trager')NULL

radnormal(diff((1/3)*sqrt(2/3)*arctanh(-3*sqrt(3/2)*(x-1)*(x-3)/sqrt((x+5)*(x-1)*(x-3)*(x-4))), x)-(-7+x)/((-11+5*x)*sqrt(x^4-3*x^3-21*x^2+83*x-60)))

-(1/18)*RootOf(_Z^2-6)*ln(-(29*RootOf(_Z^2-6)*x^2-106*RootOf(_Z^2-6)*x+41*RootOf(_Z^2-6)+36*(x^4-3*x^3-21*x^2+83*x-60)^(1/2))/(-11+5*x)^2)

(3)


 

Download Don't_Forget_Integrals!.mw 
Compare:

I wonder why by default Maple is unable to calculate the second antiderivative mentioned in the link given by the OP (Note that the substitution x = 3 - 1/(t - 1/8) can reduce some degree.) in terms of elementary functions. According to Integration Methods - method=pseudoelliptic, now Maple is capable of integrating it heuristically, but here Maple simply fails to do so (unless using method=Trager (yet the result is still less nice)).  (Besides, Mathematica is also rumored to discover human-friendly solutions to the so-called pseudo-elliptic integrals: https://github.com/stblake/algebraic_integration.)

restart;
_seed := 1234:
Warning, the use of _seed is deprecated.  Please consider using one of the alternatives listed on the _seed help page.
M := LinearAlgebra:-RandomMatrix(4, 'generator' = 0 .. 1., 'datatype' = float[4]);
st := time():
linalg:-exponential(M, -I*t):
time() - st;
                             0.141

st := time():
:-LinearAlgebra:-MatrixExponential(M, -I*t):
time() - st;
                             1.828

M := LinearAlgebra:-RandomMatrix(8, 'generator' = 0 .. 1., 'datatype' = float[4]);
st := time():
linalg:-exponential(M, -I*t):
time() - st;
                             4.219

st := time():
:-LinearAlgebra:-MatrixExponential(M, -I*t):
time() - st;
                             24.594

So the LinearAlgebra package does not always seem to be more powerful and efficient in doing linear algebra calculations?!

1 2 3 4 5 Page 1 of 5