Markiyan Hirnyk

## 7198 Reputation

11 years, 270 days

## By Quantile...

The Quantile command does the job:

```restart;
with(Statistics):
X := RandomVariable(Poisson(2.6)):
Quantile(X, .95);
5.
```

## By intersectplot...

See the first example in  ?intersectplot to this end.

## Syntax...

This works for me

```int((x^2+y^2)*piecewise(y-x-1/2 <= 0 and y+2*x-2 <= 0 and y+(1/2)*x-1/2 >= 0, 1, 0),
[x = 0 .. 1, y = 0 .. 1]);
7/32```

## Complex-valued function with branch poin...

Let us consider Int(exp(LambertW(1/(-1+x))*(-1+x)), x)+1. First, the term 1 is of no importance. Second, let us consider the int(exp(LambertW(1/(-1+x))*(-1+x)), x = 1 .. t) command instead of an inert Int one. Third,

`FunctionAdvisor(branch_cuts, exp(LambertW(1/(-1+x))*(-1+x)), plot = 2.);`

[exp(LambertW(1/(-1+x))*(-1+x)), And(-exp(1)+1 < x, x < 1)]

says the integrand is a complex-valued function with branch points. Because of this,

```a := Re(int(exp(LambertW(1/(-1+x))*(-1+x)), x = 1 .. t)):
plot(a, t = 0 .. 2, numpoints = 5, thickness = 3);
```

```int(exp(LambertW(1/(-1+x))*(-1+x)), x = 1 .. 2.1);
1.670860297
int(exp(LambertW(1/(-1+x))*(-1+x)), x = 1 .. -2.1+1.*I);
-7.784875149 - 3.643815590 I

```

int_of_Lambert.mw

## By DirectSearch:-SolveEquations...

This can be done as follows (The DirectSearch package should be downloaded from http://www.maplesoft.com/applications/view.aspx?SID=101333 and installed in your Maple >=12.).

```sol4 := DirectSearch:-SolveEquations([sol1, sol2], {T = 0 .. .4}, AllSolutions, solutions = 4);
sol4 := Matrix(1, 4, [[0.1040897473e-19,
Vector[column](2, [-0.6873203012e-10, -0.7539816154e-10]),
[T = .321611763404719, W = 29.4643511745074], 1013]])```

The plot produced by

`plots:-implicitplot([sol1, sol2], T = 0 .. 5, W = 0 .. 50, gridrefine = 2, color = [blue, red]);`

system.mw

## DiscreteUniform...

This can be done in such a way:

```with(Statistics):
X := RandomVariable(DiscreteUniform(i-x, i+x)):
Mean(X);
i
Variance(X);

(1/12)*(i-x)^2+(1/12)*(i+x)^2+(1/3)*x-(1/6*(i-x))*(i+x)
i := 6:x := 11:
Sample(X,6);
Vector[row](6, [13., 15., -3., 16., 9., -3.])```

## Two-dimensional case...

This is not so simple. Here is one way. We start from plotting

`restart; plots:-implicitplot(t^8+s^6+s*t^5 = 1, s = -10 .. 10, t = -10 .. 10, gridrefine = 2);`

We see two almost vertical pieces. This means that a small change of s implies a relatively big change of t.

Now we find the range of t by

```Optimization:-Maximize(t, {t^8+s^6+s*t^5 = 1})[1];
1.08827028359556088
```

Then we solve

```sol := solve(t^8+s^6+s*t^5 = 1);
sol := {s = RootOf(t^8+_Z^6+_Z*t^5-1), t = t}
allvalues(sol);
{s = RootOf(t^8+_Z^6+_Z*t^5-1, index = 1), t = t},
{s = RootOf(t^8+_Z^6+_Z*t^5-1, index = 2), t = t},
{s = RootOf(t^8+_Z^6+_Z*t^5-1, index = 3), t = t},
{s = RootOf(t^8+_Z^6+_Z*t^5-1, index = 4), t = t},
{s = RootOf(t^8+_Z^6+_Z*t^5-1, index = 5), t = t},
{s = RootOf(t^8+_Z^6+_Z*t^5-1, index = 6), t = t}
```

and build the parametric plots, substituting the above in

```x := t*cos(t^2+s); y := s^3-sin(t);
Digits := 15; [seq(plot([eval(x, allvalues(sol)[j]), eval(y, allvalues(sol)[j]),
t = -1.08827028359556088 .. 1.08827028359556088], thickness = 2), j = 1 .. 6)];
plots:-display(%);
```

The result is not satisfactory for me: we obtain two warnings

Warning, unable to evaluate the function to numeric values in the region; see the plotting command's help page to ensure the calling sequence is correct
Warning, unable to evaluate the function to numeric values in the region; see the plotting command's help page to ensure the calling sequence is correct
and a somewhat srange plot

which includes a superfluous interval and two breaks. I think the first one is caused by Indexed RootOf (especially by s = RootOf(t^8+_Z^6+_Z*t^5-1, index = 1) ) and the second one is caused by roundoff errors which produce complex numbers (Hope Maple developers will shed light on it.).

In oder to overcome the first deficiency, I do the following (The problem with sol[1] and t = -1 were found through trial-and-error method.).

```a := plot([eval(x, allvalues(sol)[1]), eval(y, allvalues(sol)[1]),
t = -1.08827028359556088 .. -1], thickness = 2, numpoints = 500):
b := plot([eval(x, allvalues(sol)[1]), eval(y, allvalues(sol)[1]),
t = -.99 .. 1.08827028359556088], thickness = 2, numpoints = 500):
plots:-display([seq(plot([eval(x, allvalues(sol)[j]), eval(y, allvalues(sol)[j]),
t = -1.08827028359556088 .. 1.08827028359556088], thickness = 2, numpoints = 500),
j = 2 .. 6), a, b]);
```

and

Warning, unable to evaluate the function to numeric values in the region; see the plotting command's help page to ensure the calling sequence is correct
Warning, unable to evaluate the function to numeric values in the region; see the plotting command's help page to ensure the calling sequence is correct
parametric_plot.mw

## All solutions...

can be obtained in such a way.

```expr := 1/(u+v)^2+4*u*v-1:
normal(expr);
(4*u^3*v+8*u^2*v^2+4*u*v^3-u^2-2*u*v-v^2+1)/(u+v)^2
```

Let us plot it in the positive quadrant:

`plots:-implicitplot(numer(normal(expr)) = 0, u = 0 .. 8, v = 0 .. 8, gridrefine = 2, thickness = 2);`

Now

```A := RealDomain:-solve(numer(normal(expr)) = 0, v);
(1/12)*(16*u^4-8*(64*u^6+48*u^4+12*sqrt(-192*u^6-144*u^4+288*u^2-3)*u-204*u^2+1)^(1/3)*u^2+(64*u^6+48*u^4+12*sqrt(-192*u^6-144*u^4+288*u^2-3)*u-204*u^2+1)^(2/3)+8*u^2+(64*u^6+48*u^4+12*sqrt(-192*u^6-144*u^4+288*u^2-3)*u-204*u^2+1)^(1/3)+1)/(u*(64*u^6+48*u^4+12*sqrt(-192*u^6-144*u^4+288*u^2-3)*u-204*u^2+1)^(1/3))                                                                                              ```

For example,

```eval(A, u = 2);
(1/24)*(289-31*(4049+24*sqrt(-13443))^(1/3)+(4049+24*sqrt(-13443))^(2/3))/(4049+24*sqrt(-13443))^(1/3)
evalf(%);

0.09656230030 + 2.101279539 10 ^(-10)   I
```

Let us chek it by plotting

`plot(A, u = 1 .. 8);`

Because expr is symmetric in u and v, the rest is clear.

## By solve...

One can solve this system analitically by the command of Maple 2016.2

`solve({-5 <= 4*y1+2*y2+3*y3, -2 <= 3*y1+5*y2+2*y3, -1 <= y1+2*y2+y3, 0 <= y2, y1 <= 0}, {y1, y2, y3});`

{y2 = 0, y3 = -(3/2)*y1-1, y1 < -2}, {y2 = 0, y1 < -2, -(3/2)*y1-1 < y3}, {y1 = -2, y2 = 0, 2 <= y3}, {y2 = 0, y3 = -(3/2)*y1-1, -2 < y1, y1 < 0}, {y2 = 0, -2 < y1, y1 < 0, -(3/2)*y1-1 < y3}, {y1 = 0, y2 = 0, -1 <= y3}, {y3 = -1-(3/2)*y1-(5/2)*y2, 0 < y2, y1 < -2+4*y2, y2 < 4/11}, {0 < y2, y1 < -2+4*y2, y2 < 4/11, -1-(3/2)*y1-(5/2)*y2 < y3}, {y1 = -2+4*y2, y3 = 2-(17/2)*y2, 0 < y2, y2 < 4/11}, {y1 = -2+4*y2, 0 < y2, y2 < 4/11, 2-(17/2)*y2 < y3}, {y3 = -1-(3/2)*y1-(5/2)*y2, 0 < y2, y1 < -y2, y2 < 4/11, -2+4*y2 < y1}, {0 < y2, y1 < -y2, y2 < 4/11, -2+4*y2 < y1, -1-(3/2)*y1-(5/2)*y2 < y3}, {y1 = -y2, y3 = -y2-1, 0 < y2, y2 < 4/11}, {y1 = -y2, 0 < y2, y2 < 4/11, -y2-1 < y3}, {y3 = -1-y1-2*y2, 0 < y2, y1 < 0, y2 < 4/11, -y2 < y1}, {0 < y2, y1 < 0, y2 < 4/11, -y2 < y1, -1-y1-2*y2 < y3}, {y1 = 0, y3 = -2*y2-1, 0 < y2, y2 < 4/11}, {y1 = 0, 0 < y2, y2 < 4/11, -2*y2-1 < y3}, {y2 = 4/11, y3 = -21/11-(3/2)*y1, y1 < -6/11}, {y2 = 4/11, y1 < -6/11, -21/11-(3/2)*y1 < y3}, {y1 = -6/11, y2 = 4/11, -12/11 <= y3}, {y2 = 4/11, y3 = -21/11-(3/2)*y1, -6/11 < y1, y1 < -4/11}, {y2 = 4/11, -6/11 < y1, y1 < -4/11, -21/11-(3/2)*y1 < y3}, {y1 = -4/11, y2 = 4/11, -15/11 <= y3}, {y2 = 4/11, y3 = -19/11-y1, -4/11 < y1, y1 < 0}, {y2 = 4/11, -4/11 < y1, y1 < 0, -19/11-y1 < y3}, {y1 = 0, y2 = 4/11, -19/11 <= y3}, {y3 = -1-(3/2)*y1-(5/2)*y2, 4/11 < y2, y1 < -2+4*y2, y2 < 2/5}, {4/11 < y2, y1 < -2+4*y2, y2 < 2/5, -1-(3/2)*y1-(5/2)*y2 < y3}, {y1 = -2+4*y2, y3 = 2-(17/2)*y2, 4/11 < y2, y2 < 2/5}, {y1 = -2+4*y2, 4/11 < y2, y2 < 2/5, 2-(17/2)*y2 < y3}, {y3 = -1-(3/2)*y1-(5/2)*y2, 4/11 < y2, y1 < -y2, y2 < 2/5, -2+4*y2 < y1}, {4/11 < y2, y1 < -y2, y2 < 2/5, -2+4*y2 < y1, -1-(3/2)*y1-(5/2)*y2 < y3}, {y1 = -y2, y3 = -y2-1, 4/11 < y2, y2 < 2/5}, {y1 = -y2, 4/11 < y2, y2 < 2/5, -y2-1 < y3}, {y3 = -1-y1-2*y2, 4/11 < y2, y1 < 4-11*y2, y2 < 2/5, -y2 < y1}, {4/11 < y2, y1 < 4-11*y2, y2 < 2/5, -y2 < y1, -1-y1-2*y2 < y3}, {y1 = 4-11*y2, y3 = -5+9*y2, 4/11 < y2, y2 < 2/5}, {y1 = 4-11*y2, 4/11 < y2, y2 < 2/5, -5+9*y2 < y3}, {y3 = -1-y1-2*y2, 4/11 < y2, y1 < 0, y2 < 2/5, 4-11*y2 < y1}, {4/11 < y2, y1 < 0, y2 < 2/5, 4-11*y2 < y1, -1-y1-2*y2 < y3}, {y1 = 0, y3 = -2*y2-1, 4/11 < y2, y2 < 2/5}, {y1 = 0, 4/11 < y2, y2 < 2/5, -2*y2-1 < y3}, {y2 = 2/5, y3 = -2-(3/2)*y1, y1 < -2/5}, {y2 = 2/5, y1 < -2/5, -2-(3/2)*y1 < y3}, {y1 = -2/5, y2 = 2/5, -7/5 <= y3}, {y2 = 2/5, y3 = -9/5-y1, -2/5 < y1, y1 < 0}, {y2 = 2/5, -2/5 < y1, y1 < 0, -9/5-y1 < y3}, {y1 = 0, y2 = 2/5, -9/5 <= y3}, {y3 = -1-(3/2)*y1-(5/2)*y2, 2/5 < y2, y1 < 4-11*y2, y2 < 1/2}, {2/5 < y2, y1 < 4-11*y2, y2 < 1/2, -1-(3/2)*y1-(5/2)*y2 < y3}, {y1 = 4-11*y2, y3 = 14*y2-7, 2/5 < y2, y2 < 1/2}, {y1 = 4-11*y2, 2/5 < y2, y2 < 1/2, 14*y2-7 < y3}, {y3 = -5/3-(4/3)*y1-(2/3)*y2, 2/5 < y2, y1 < -y2, y2 < 1/2, 4-11*y2 < y1}, {2/5 < y2, y1 < -y2, y2 < 1/2, 4-11*y2 < y1, -5/3-(4/3)*y1-(2/3)*y2 < y3}, {y1 = -y2, y3 = -5/3+(2/3)*y2, 2/5 < y2, y2 < 1/2}, {y1 = -y2, 2/5 < y2, y2 < 1/2, -5/3+(2/3)*y2 < y3}, {y3 = -5/3-(4/3)*y1-(2/3)*y2, 2/5 < y2, y1 < -2+4*y2, y2 < 1/2, -y2 < y1}, {2/5 < y2, y1 < -2+4*y2, y2 < 1/2, -y2 < y1, -5/3-(4/3)*y1-(2/3)*y2 < y3}, {y1 = -2+4*y2, y3 = -6*y2+1, 2/5 < y2, y2 < 1/2}, {y1 = -2+4*y2, 2/5 < y2, y2 < 1/2, -6*y2+1 < y3}, {y3 = -1-y1-2*y2, 2/5 < y2, y1 < 0, y2 < 1/2, -2+4*y2 < y1}, {2/5 < y2, y1 < 0, y2 < 1/2, -2+4*y2 < y1, -1-y1-2*y2 < y3}, {y1 = 0, y3 = -2*y2-1, 2/5 < y2, y2 < 1/2}, {y1 = 0, 2/5 < y2, y2 < 1/2, -2*y2-1 < y3}, {y2 = 1/2, y3 = -9/4-(3/2)*y1, y1 < -3/2}, {y2 = 1/2, y1 < -3/2, -9/4-(3/2)*y1 < y3}, {y1 = -3/2, y2 = 1/2, 0 <= y3}, {y2 = 1/2, y3 = -2-(4/3)*y1, -3/2 < y1, y1 < -1/2}, {y2 = 1/2, -3/2 < y1, y1 < -1/2, -2-(4/3)*y1 < y3}, {y1 = -1/2, y2 = 1/2, -4/3 <= y3}, {y2 = 1/2, y3 = -2-(4/3)*y1, -1/2 < y1, y1 < 0}, {y2 = 1/2, -1/2 < y1, y1 < 0, -2-(4/3)*y1 < y3}, {y1 = 0, y2 = 1/2, -2 <= y3}, {y3 = -1-(3/2)*y1-(5/2)*y2, 1/2 < y2, y1 < 4-11*y2}, {1/2 < y2, y1 < 4-11*y2, -1-(3/2)*y1-(5/2)*y2 < y3}, {y1 = 4-11*y2, y3 = 14*y2-7, 1/2 < y2}, {y1 = 4-11*y2, 1/2 < y2, 14*y2-7 < y3}, {y3 = -5/3-(4/3)*y1-(2/3)*y2, 1/2 < y2, y1 < -y2, 4-11*y2 < y1}, {1/2 < y2, y1 < -y2, 4-11*y2 < y1, -5/3-(4/3)*y1-(2/3)*y2 < y3}, {y1 = -y2, y3 = -5/3+(2/3)*y2, 1/2 < y2}, {y1 = -y2, 1/2 < y2, -5/3+(2/3)*y2 < y3}, {y3 = -5/3-(4/3)*y1-(2/3)*y2, 1/2 < y2, y1 < 0, -y2 < y1}, {1/2 < y2, y1 < 0, -y2 < y1, -5/3-(4/3)*y1-(2/3)*y2 < y3}, {y1 = 0, y3 = -(2/3)*y2-5/3, 1/2 < y2}, {y1 = 0, 1/2 < y2, -(2/3)*y2-5/3 < y3}

This approach works in higher dimensions too.

## By listtorec and rsolve and sum...

This can be done as follows.

`p := convert(a-(1/2)*beta*a^2*y^2+(1/24)*beta^2*a^3*y^4-(1/720)*beta^3*a^4*y^6+(1/40320)*beta^4*a^5*y^8-(1/3628800)*beta^5*a^6*y^10+(1/479001600)*beta^6*a^7*y^12-(1/87178291200)*beta^7*a^8*y^14+O(y^16), polynom);`

a-(1/2)*beta*a^2*y^2+(1/24)*beta^2*a^3*y^4-(1/720)*beta^3*a^4*y^6+(1/40320)*beta^4*a^5*y^8-(1/3628800)*beta^5*a^6*y^10+(1/479001600)*beta^6*a^7*y^12-(1/87178291200)*beta^7*a^8*y^14

`L := [seq(coeff(p, y, n), n = 0 .. 14)];`

[a, 0, -(1/2)*beta*a^2, 0, (1/24)*beta^2*a^3, 0, -(1/720)*beta^3*a^4, 0, (1/40320)*beta^4*a^5, 0, -(1/3628800)*beta^5*a^6, 0, (1/479001600)*beta^6*a^7, 0, -(1/87178291200)*beta^7*a^8]

```with(gfun):
rec:=listtorec(L, u(n));
```

[{beta*a*u(n+1)+(n^2+5*n+6)*u(n+3), u(0) = a, u(1) = 0, u(2) = -(1/2)*beta*a^2}, ogf]

`rsolve(rec[1], u);`

(1/2)*a*((-sqrt(-beta*a))^n+(-beta*a)^((1/2)*n))/GAMMA(n+1)

`sum(%*y^n, n = 0 .. infinity);`

a*cos(y*sqrt(beta*a))

rec.mw

## No solution...

```eq1 := x[1]-x[2] = 0:
eq2 := -x[1]+2*x[2]-x[3] = 0:
eq3 := -x[2]+2*x[3]-x[4] = 0:
eq4 := -x[3]+x[4]-t = 0:
solve({eq1, eq2, eq3, eq4}, [x[1], x[2], x[3], x[4]]);
[]```

Up to ?solve

• If the solve command does not find any solutions, then if the second argument is a name or set of names, then the empty sequence (NULL) is returned; if the second argument is a list, then the empty list is returned. This means that there are no solutions, or the solve command cannot find the solutions. In the second case, a warning is issued, and the global variable _SolutionsMayBeLost is set to true.

this indicates no solution.

## Another way...

is as follows.

``` restart; Digits := 15:
f := .9/abs(t-.4)^(1/3)+.1/abs(t-.6)^(1/2):
G1 := -.9445379894:
A := dsolve({a(0) = 0, (D(a))(t) = f*exp(t)}, numeric);
B := dsolve({b(0) = 0, (D(b))(t) = f*exp(-t)}, numeric);
rhs(B(.2)[2]);
.270547407517139
U1 := proc (x) options operator, arrow; -(1/2)*exp(-x)*(rhs((eval(A(t), t = x))[2])+G1)-(1/2)*exp(x)*(rhs((eval(B(t), t = x))[2])+G1) end proc:
U := proc (x) options operator, arrow; -(1/2)*exp(x)*(rhs((eval(B(t), t = x))[2])+G1)+(1/2)*exp(-x)*(rhs((eval(A(t), t = x))[2])+G1) end proc:
F := proc (s) options operator, arrow; U(s)^2+U1(s)^2-(2*.9/abs(s-.4)^(1/3)+2*.1/abs(s-.6)^(1/2))*U(s) end proc:
F(.2);
HFloat(-0.08307649472006268)
evalf(Int(F, 0 .. 1, method = _d01ajc), 8);
-0.37533515
evalf(Int(F, 0 .. 1, method = _d01ajc), 7);
-0.3753352
```

In the above the Joe Riel's idea is used.

another_way.mw

## Megahit...

Similar questions were asked and answered many times. Procedures are plotted specifically. Try

`plot(myPi_1, 0 .. 500, numpoints = 500, thickness = 2, color = black);`

## Diverges...

The integral

`int(U1(x)^2+U(x)^2-2*f(x)*U(x), x=0..1)`

diverges because of the singularity at x=3/5. Here are my arguments.

First, let us switch to exact calculations by

```f := convert(.9/abs(x-.4)^(1/3)+.1/abs(x-.6)^(1/2), rational):
G1 := convert(-.9445379894, rational):```

Second, the latest Maple versions (Maple 2015 and Maple 2016) handle  integrals of  moduli:

```a := int(f*exp(t), t = 0 .. x);
(1/10)*(exp(x)*abs(x-2/5)^(1/3)+9*exp(x)*sqrt(abs(x-3/5))-abs(x-2/5)^(1/3)-9*sqrt(abs(x-3/5)))/(abs(x-2/5)^(1/3)*sqrt(abs(x-3/5)))
b := int(f*exp(-t), t = 0 .. x);
(1/10)*(exp(x)*abs(x-2/5)^(1/3)+9*exp(x)*sqrt(abs(x-3/5))-abs(x-2/5)^(1/3)-9*sqrt(abs(x-3/5)))*exp(-x)/(abs(x-2/5)^(1/3)*sqrt(abs(x-3/5)))```

Now

```U1 := -(1/2)*exp(-x)*(a+G1)-(1/2)*exp(x)*(b+G1):
U := -(1/2)*exp(x)*(b+G1)+(1/2)*exp(-x)*(a+G1):
limit((U^2-2*U*f+U1^2)*abs(x-3/5), x = 3/5);
(1/200)*exp(6/5)-1/100+(1/200)*exp(-6/5)
evalf(%);
0.008106555680
```

Therefore, the integrand behaves as constant/|x-3/5| around x=3/5. This implies the divergence.

question.docx

 1 2 3 4 5 6 7 Last Page 1 of 139
﻿