## 4009 Reputation

6 years, 176 days

## How to integrate this expression?...

Maple

Let

```f := beta/(2*sigma*GAMMA(1/beta))*exp(-(abs(x-mu)/sigma)^beta);
```

where mu::real, beta > 0, sigma > 0.
How can we obtain the expression of F?

`F := int(f, x=-infinity..s)`

Reverse problem (a priori simpler):
It's known that

`F := 1/2+signum(x-mu)/(2*GAMMA(1/beta))*(GAMMA(1/beta)-GAMMA(1/beta, abs((x-mu)/sigma)^beta))`

How can we check that diff(F, x)=f?
Even with the assumptions on beta, mu and sigma and additional assumptions x>0 or x<0, I can't verify that diff(F, x)=f.

## Could you help me find the value of P?...

Maple 2015

I want to estimate numerically the value of P (that is the probability that f exceeds the value 12 when x, y and z are uniformly distributed within the box [-Pi, Pi]3).

```restart
f := sin(x)+7*sin(y)^2+0.1*z^4*sin(x);
Omega := [x, y, z] =~ [(-Pi..Pi)\$3]:
2        4
sin(x) + 7 sin(y)  + 0.1 z  sin(x)

h := Heaviside(f-12):
P := Int(h, Omega);
```

Here is a simple Monte Carlo estimation of P.

```f_MC  := x -> sin(x[1])+7*sin(x[2])^2+0.1*x[3]^4*sin(x[1]);
h_MC  := x -> Heaviside(f_MC(x) - 12);
omega := -Pi, Pi;

P_MC := proc(N)
local Z := Statistics:-Sample(Uniform(omega), [N, 3]):
local F := Vector(N, n -> h_MC(Z[n])):
local K := add(F):
local P := evalf(add(F)/N):
local q := Statistics:-Quantile(Normal(0, 1), 0.005, numeric):
local e := -q*sqrt(P*(1-P)/N):
printf("A bilateral 99%% confidence interval that Prob(f > 12) is %1.3e +/- %1.2e\n", P, e);
end proc:

P_MC(10^6)
A bilateral 99% confidence interval that Prob(f > 12) is 1.643e-02 +/- 3.27e-04
```

I was hoping to get a value for P using  evalf/Int with some method.
But I didn't succeed with any of the methods for multiple integration:

• The following command returns 0 for n=1 and 2 and Int(h, Omega) if n >=3
```evalf( Int(h, Omega, 'epsilon=1e-n', 'method=_MonteCarlo') );
```

• And all my attempts with methods from the Cuba library have resulted in an unevaluated Int(h, Omega) integral.

Could you help me to estimate P using  evalf/Int ?
TIA

## Wrong result in integrating a Heaviside ...

Maple 2015

Hi,
Can anyone tell me what's going wrong with int(Heaviside(f-1/2), x=0..1, y=0..1); ?

```restart:

f := x*y:

JH := int(Heaviside(f-1/2), x=0..1, y=0..1);  #???
1
-
4

JP := int(piecewise(x*y>1/2, 1, 0), x=0..1, y=0..1);
1   1
- - - ln(2)
2   2

# integrate "1" over the domain where f > 1/2

S := solve(f=1/2, x):
JS := int(1, x=S..1, y=1/2..1);
1   1
- - - ln(2)
2   2
```

## Do you have any idea how to solve this p...

Maple

In a french magazine written by High Schools teachers I found this problem:

```let a, b, p, q four strictly positive integers such that a > b^2 and p > q+1;
find 4-tuples (a, b, p, q) such that

(a^2 - b^4) = p!/q!```

Given the source of this problem I suspect that there is a trick to answering this question.
After some hours spent, I have found no general method to solve it, only a few solutions (first one and second one are almost obvious), for instance

```rel := a^2 - b^4 = p!/q!:

eval(rel, [a= 5, b=1, q=1, p=4]);
eval(rel, [a=11, b=1, q=1, p=5]);
eval(rel, [a=71, b=1, q=1, p=7]);
eval(rel, [a= 2, b=1, q=2, p=3]);
eval(rel, [a=19, b=1, q=2, p=6]);
eval(rel, [a=21, b=3, q=2, p=6]);

```

Do you have any idea how to solve this problem?
Could it be handled by Maple (without a systematic exploration of a part of N^4)?

## How to find the maximum value of this fu...

Maple 2015

Hi,

I solve numerically an ODE system which depends on 25 parameters.
I want to know the maximum value of the time at which a specific event is triggered, when these parameters vary (independently the one of the other) within a 25 dimensional hyperbox.

To solve this maximization problem, which I assume is local1, I would like to use NLPSolve

The attached zip (I use Mac OSX) contains:

• a file that gives a full description of the problem (the original ".pages" file and its export as pdf),
• an m file wich contains the solution procedure plus a few other variables needed to solve the problem,
• the mw file which contains:
• the reading of the m file,
• the procedure (OBJ) which returns the trigerring time,
• many attemps to find its maximum value using NLPSolve.

None of my attempts at using NLPSolve gave me the expected answer

Could you help me to fix this?