Axel Vogt

5936 Reputation

20 Badges

20 years, 250 days
Munich, Bavaria, Germany

MaplePrimes Activity


These are answers submitted by Axel Vogt

Is this similar to your older question at http://www.mapleprimes.com/questions/201767-Fsolve-Cannot-Solve-This-Equations-And but now with log(-exp(...)) ?

I am not sure whether 3 1/2 is a good idea in Maple (but I always do use multiplication sign).
However the following gives something close
3 + ``(1/2); 
                              3 +  (1/2)

expand(%);

                                 7/2

There are various ways. The inner integral is Int(f1, y=0 .. h(x)).
Compute an anti-derivative and name ist V. To apply the fundamental
theorem V has to be continuous along the curve h(x). Substitude y =
ymax and plot real and imaginary part

  V;
  eval(%, y=ymax);
  plot([Re(%),Im(%)], x=0 .. Pi, color=[red,blue]);

You 'see' one will split in x=Pi/2 (and moreover the imaginary part
will drop out). Finally I get 3.13316121793696 (15 decimals)

PS: What is the result given by Matlab?
evalf(Int(u(t)^2, t = 0..1)); # upper case
    .448698034001023
Note that your 1st integral needs 0 < k (take t=0, k= -1, lambda= -1,
then it is not finite).

Changing by k*s = x, abbreviating k=lambda/a one finally has the task

task:=Int(exp(-exp(-x)+(a-1)*x),x = -infinity .. t); # a new t

Changing by x = -ln(-ln(y)) this is

task1 := Int((-ln(y))^(-a),y = 0 .. exp(-exp(-t)))

The exp term is between 0 and 1 (plot it, for example), so y is between
0 and 1. Hence z = ln(y) is between -oo and 0 (smells like a pdf, no?)

We have (-z)^(-a) = z^(-a)*(-1)^(-a)*exp(2*I*Pi*a) and for z=ln(y):

task1 = faktor * task2, task2 := Int(ln(y)^(-a),y = 0 .. exp(-exp(-t))),
faktor := (-1)^(-a)*exp(2*I*Pi*a)

In that form one can cheat asking Mathematica (or perhaps Gradsteyn?),
Anti := GAMMA(-a+1,-ln(y))*ln(y)^(-a)/((-ln(y))^(-a)) is an antiderivative

eval(Anti, y= exp(-exp(-t))) - MultiSeries:-limit(Anti, y=0, right);
my:=simplify(%) assuming t::real;

is an intermediate solution (no check for continuity) and then

faktor*my;
simplify(%) assuming a::real;
evalc(%);
My:=simplify(%);

                     My := GAMMA(-a + 1, exp(-t))
 
is a solution, I think.

Check 1: in integers

task = My;
'[seq(%, a= - 3 .. 3)]'; value(%): expand(%): simplify(%); map(is, %);

              [true, true, true, true, true, true, true]

 
Check 2: in ugly floats

Digits:=15;
task = My;
subs(t=1, %);
'[seq(%, a= - 3 .. 3, 1 + evalf(1/Pi))]'; evalf(%); #map(is, %);

  [5.99658218952224 = 5.99658218952226,

        1.50284712032528 = 1.50284712032527,

        0.737383845311871 = 0.737383845311871,

        0.750709217778417 = 0.750709217778418,

        1.29447312393282 = 1.29447312393282]


So 'task' is an integral for the incomplete Gamma function and knowing
that there should be a more direct way.

From that your initial 1st task will follow. The 2nd is just a product, no?

PS: Where does your series of questions come from?
The inner integral is of the form Int(f(x), x=0.. phi(y)) and using
IntegrationTools[Change](%, x=z*phi(y), z) this becomes

phi(y)*Int(f(z*phi(y)),z = 0 .. 1)


Now combine with the original task one can write a double integral
Int(K, [z=0..1, y=0..Pi], epsilon=1e-8) and evalf gives a good result:

      3.13316050306181 + .647491e-17*I

I did it a bit differently:

Maple finds an antiderivative w.r.t. z. Call it val, "starring" at it
result this seems to be continous in for those y, so I get a formal
solution for the inner integral Int(K, z=0 .. 1).

val:=int(K, z);
S:=eval(%,z=1) - eval(%,z=0)

Now Int(Re(S), y=0..Pi, method = _d01ajc); evalf(%);

      3.13316050516358

= -.666758443186807e-2 + (11/2*Ei(1,11/2)-exp(-11/2)+1)*Pi

For Im(S) like the above it is almost vanishing.

Plotting K's imaginary part a point-symmetry is seen in y=Pi/2, and
plotting left side + right side gives me a zero plot (ok, no proof),
so I would only consider the real part.

I used Digits:=15

The most easy way ist Compiler:-Compile, it comes with Maple.

Else you may search this site for "DLL", this should bring up most of stuff (may be it does not dicuss issues of 32 vs 64 bit DLL). And of course you can can look into the examples coming with Maple (+ coding manual), subdirectory "samples" of your installation.

I would suggest to use Compiler:-Compile - except you have _good_ reason to take another route.

I think it is roughly  exp( 0.3 E16331116 ), could that be?

I made an error, sorry. I get Carl's result as .4666594270e-3: converting to rationals the inner integral can be done anlytically (assuming 1 < eta) and then the overall integral can be computed with relativ error of 1E-10 quickly using method = _d01ajc

It "reduces" to 0 =  -2^(-1/4+1/4*c^2+c)+2^(-1/4-c+1/4*c^2)-c, c=cos(2*x)
by using 'combine'

Maple can not show "only for c=0", neither can (online-) MMA, where I mean
c (or x) to be complex.

From that it follows x = 1/4 *(2 *Pi *n+Pi), using Maple, as MMA asserts.
cos(x)^(n-1)*sin(x)+n*(Int(cos(x)^(n-2), x))-(Int(cos(x)^(n-2), x)); 
combine(%, Int);
simplify(%, size);

           /
          |        (n - 2)                    (n - 1)
          |  cos(x)        (n - 1) dx + cos(x)        sin(x)
          |
         /

Your recursion formula is at Wikipedia as well (with ugly proof), http://en.wikipedia.org/wiki/Integration_by_reduction_formulae#Examples

Using Maple you can "verfiy" it:

F(n) = cos(x)^(n-1)*sin(x)/n +(n-1)/n*F(n-2);
eval(%, F= 'n -> Int(cos(t)^n, t=0 .. x)' ); # substitute the integral
'diff(%, x)';
simplify(%); # so up to constants it is ok
 
eval(%%%, x=0);

Now one can use 'rsolve' to solve the recursion, giving a noisy result

Have not tried cos(x)^n = Sum in lower terms seriously (that sum is what is given by the hypergeometric 2F1 and the correcting factors)

CAmax:=(-k2+(k2*k1)^(1/2))*CA0/(k1-k2);
CAmaxgiven:= CA0*(k2/k1)^(1/2)/(1+(k2/k1)^(1/2));
CAmax - CAmaxgiven:
simplify(%) assuming 0<k1,0<k2;

                                  0

Edited: or even more simple:

CAmax - CAmaxgiven;
subs(CA0=1, %); subs(k1=1,%); # justify that !
normal(%);

                                  0

I do not give much for a factor of 2. For real speed one may want to have specific SW for specific tasks. But the price may be high (setting it up, learning). So I would ask: for what? I am satisfied, that Mpl covers that, if I use it as "working horse". It would be more interesting to compare speed beyond double precision (kicking out Matlab). And even then: personally I will spend most of time to feed the problem. Else: use specific SW. Just my 0.02 €.

I would prefer to map a function (also see Preben Alsholm below)

  f:= x -> `if`(1=x, 0, x);
  map(f,A);
                     f := x -> `if`(1 = x, 0, x)


                            [0    0    6]
                            [           ]
                            [5    3    2]
                            [           ]
                            [0    0    0]

A mathematical answer:

This is an implicit definition, valid in Math, if you *know* there is
an explicit, unique solution for f(x) = something (i.e. that f has a
[partial] left inverse).

That is not a good idea for a CAS (and perhaps related to formal proofs)
it is even fals for f = sqrt (having 2 solutions).

One way may be:

  omega:= 'omega':                        # clear variable
  omega^2 = sqrt(w0^2*(1+((z-zf)/z0)^2)); # set up taks
  [solve(%, omega)];                      # *solve* it (if possible)
 
  omega:= %[1];                           # select (!) the very solution


NB: note that this might make no sense in some context if you leave the
Reals and work in other domains, https://en.wikipedia.org/wiki/Modular_square_root
First 15 16 17 18 19 20 21 Last Page 17 of 93