Axel Vogt

5936 Reputation

20 Badges

20 years, 253 days
Munich, Bavaria, Germany

MaplePrimes Activity


These are replies submitted by Axel Vogt

1. The task: 
I did not follow your 'physical' thinking, just use IR^2 ~ C and take polar coordinates. Here it goes:
#        IR^2              ->  IR
  f:=   (x,y)              ->  exp(-x^2-y^2);

Compute it using Maple's routines:

  Int(Int(f(x,y),x=-infinity..infinity),y=-infinity...infinity);
  value(%);
                infinity    infinity
               /           /
              |           |                2    2
              |           |          exp(-x  - y ) dx dy
              |           |
             /           /
               -infinity   -infinity


                                  Pi

Now do a change of variables using the formula above:

#        IR * ]-Pi,Pi[ = W ->  V = IR^2
  Phi:= (r,theta)          ->  [r*cos(theta),r*sin(theta)];

  'f( op(Phi(r,theta)) )*detJ(Phi)(r,theta)':    # the change of variables
  'Int(Int(%, r=0..infinity), theta=-Pi .. Pi)'; # according integral
  %;
  simplify(%); value(%);

               f(op(Phi(r, theta))) detJ(Phi)(r, theta)


     Pi    infinity
    /     /
   |     |
   |     |          f(op(Phi(r, theta))) detJ(Phi)(r, theta) dr d
   |     |
  /     /
    -Pi   0

        theta


     Pi    infinity
    /     /
   |     |                2           2    2           2             2                 2
   |     |          exp(-r  cos(theta)  - r  sin(theta) ) (cos(theta)  r + r sin(theta) ) dr dtheta
   |     |
  /     /
    -Pi   0



                  infinity                  Pi
                 /                         /
                |                2        |
                |          exp(-r ) r dr  |    1 dtheta
                |                         |
               /                         /
                 0                         -Pi


                                  Pi


Note that one needs op(phi(...)) to switch from the 'vector notation' to sequences
which is what f expects as input.

--------------------------------------------------------------------

2. I am aware that the practical problem is to find appropriate transforms & coordinates.
But that fits into Alejandro Jakubi's answers above :-)

Perhaps I do not quite understand the question, but Rudin "Real and Complex Analysis"
has it (never liked that great book as a student):

Let be f a C^1 function on V, V bounded and open IR^k and phi: W -> V a differentiable
homeomorphism (or what notion one prefers), psi its inverse.

  inv := f -> (f@@(-1));
  psi:=inv(phi);
                                      (-1)
                            psi := phi



  Int(f(x),x=`V` ..``) = Int(f(phi(xi))*detJacobi(phi)(xi),xi =`W` ..``); `W`= psi(`V`) ;

    /              /
   |              |
   |   f(x) dx =  |   f(phi(xi)) detJacobi(phi)(xi) dxi,
   |              |
  /              /
    V              W

                (-1)
        W = (phi    )(V)


where detJacobi is the determinant of the Jacobian.

For the 1-dimensional case that equals D, the differential, and it reads as

  Int(f(x),x=`V` ..``) = Int(f(phi(xi))*detJacobi(phi)(xi),xi =`W` ..``):
  eval(%, detJacobi=D);

              /              /
             |              |
             |   f(x) dx =  |   f(phi(xi)) D(phi)(xi) dxi
             |              |
            /              /
              V              W

where now V and W can be thought as (an union of) intervals.

Now using intervals we write the integrals in the common way 

  I=[`a...b`], `J`= psi(`I`) ;
  Int(f(x),x=a ..b) = Int(f(x),x=`I` ..``), `...`;

                                         (-1)
                    I = [a...b], J = (phi    )(I)


                      b
                     /              /
                    |              |
                    |   f(x) dx =  |   f(x) dx, ...
                    |              |
                   /              /
                     a              I


and get

  Int(f(x),x=a ..b) = Int(f(phi(xi))*diff(phi(xi),xi),xi = psi(a) .. psi(b));

                           (-1)
        b              (phi    )(b)
       /              /
      |              |                         / d         \
      |   f(x) dx =  |              f(phi(xi)) |--- phi(xi)| dxi
      |              |                         \dxi        /
     /              /
       a                  (-1)
                      (phi    )(a)

which is the same as using "Int(f(x),x=a ..b); Change(%, x=phi(xi));" getting
the ugly 'RootOf' notation.


For the Jacobian one would need to know the number of variables to set up a
general procedure like the 'D' operator, but I only have a lame way (which
works more or less for procedures, but not for functions like 'exp', I would
need to re-write it as proc(x) x -> [exp(x)] end proc as vector version).

  detJ := proc( h::procedure )
  local n,X,i,Y, result;
 
  [op(eval(h))]; 
  remove(has,%,{arrow,operator}); 
  n:=nops(%); # number of variables of h
  
  if(n < 1) then
    return NULL;
  elif n=1 then 
  #  return D(h)
  else
  end if;
 
  X:=seq( cat( x, i ), i = 1..n); i:='i':
  Y:= seq(h(X)[i], i = 1 .. n);  i:='i': 

  result:=Student[MultivariateCalculus][Jacobian]([Y], [X], output = matrix); 
  if(1 < n) then result:=LinearAlgebra[Determinant](result) end if;

  unapply(result, X);
  end proc;

One needs the partial derivatives to commute to use Maple's Jacobian, yes?

I have not tried it for product situations.

The usual way is ... to buy it ... or get it from a library.

You may accept that not every thing is for free ... of course only on a legal basis ...

The usual way is ... to buy it ... or get it from a library.

You may accept that not every thing is for free ... of course only on a legal basis ...

Danke - das freut mich :-)

Gruß
Axel

just use the "Standard interface", than you can see it

for this either use thw Windows context menu to open using ...

or make a copy, rename it choosing mw as extension (not mws), than it should open with
the Standard interface - run it again and you should see the desired picture (at least for me
it then works)

just use the "Standard interface", than you can see it

for this either use thw Windows context menu to open using ...

or make a copy, rename it choosing mw as extension (not mws), than it should open with
the Standard interface - run it again and you should see the desired picture (at least for me
it then works)

same problem with that approach: ok in standard and funny in classical interface, but ok using Maple 9.5

same problem with that approach: ok in standard and funny in classical interface, but ok using Maple 9.5

Even if it was my intention to give a somwwhat detailed error report I put together
what was done above do give the result in Maple.

  Int( arccos(x)/(x - a), x = 0 .. 1 ) =
= Int(-t*sin(t)/(-cos(t)+a),t = 0 .. 1/2*Pi) =
= Int(ln(-cos(t)+a),t = 0 .. 1/2*Pi) - 1/2*Pi*ln(a) =
= Int(ln(-x+a)/(1-x^2)^(1/2),x = 0 .. 1) - 1/2*Pi*ln(a) = # let Maple evaluate it
= ln( (1+sqrt(1-1/(a^2)))/2 )*Pi/2 - 1/a*hypergeom([1/2, 1, 1],[3/2, 3/2],1/(a^2))

Let me call E the function given by the last expression.

Then it makes sense to use (complex) partial fraction for the original problem
using convert(1 / ( 1+x^4),fullparfrac) and this then writes as

  'Sum(-1/4*alpha*E(alpha),alpha = RootOf(1+Z^4))';
  allvalues(%): simplify(%,size); 

  evalf(%);

                     -----
                      \
                       )      (-alpha * E(alpha) / 4)
                      /
                     -----
              alpha = RootOf(1+Z^4)

  1/2 hypergeom([1/2, 1, 1], [3/2, 3/2], I) + 1/2 hypergeom([1/2, 1, 1], [3/2, 3/2], -I)

                         0.9225481756 + 0. I

which is the numerical value for  Int( arccos(x) / ( 1+x^4) , x=0 .. 1).


One can write the result as Re( hypergeom([1/2, 1, 1], [3/2, 3/2], I) ), but still
in terms of hypergeometrics 3F2, which I do not like too much - but have nothing
better.

May be one can do Int(ln(-x+a)/(1-x^2)^(1/2),x = 0 .. 1) through contour integrals,
but I am too long out of practicing that.

admin: "Hello, here is your service desk, may I help you?"

user:  "jlkajghdlazböjvöldavörlaeknrunröaelknoinätorksm ...!!!"

admin: "So the pages are saved as text files?"

user:  "Damned, YES!"

admin: "May it be you save them as text?"

user:  "ARGGHH - I am not an idiot!"

admin: "Sure ... I mean: sure, you are not. Could we try?"

user:  "If YOU think that could help you to solve YOUR task ..."

.
.
.

admin: "What can you read in the 2nd box?"

user:  "Hm ..."

admin: "Pardon?"

user:  " ... just try to reproduce the error ..."

admin: "No problem ... I have remote access if you allow that"

user:  "Will see later, have a VERY important meeting NOW ..."

admin: "Sure ... if the problem persist please contact me again."

Using Firefox I just use the menu to save the file ( = current page) as html - as I always did & it worked ... but if nothing changed at this site for weeks the reason must be my installation ... will check it ... may be an upgrade to the latest Firefox can solve it.

Thanks for the answer, Will, will let you know

I especially like this one, might be a prologue for the homework section:

 

asnwer times 2

It reminds me a lot how I learned Math :-)

k, from which I took as motivation:

groups.google.de/group/comp.soft-sys.math.maple/browse_thread/thread/c71f510677102dca/
where MMA is said to give an answer, hypergeom([1/4, 1/2, 1, 1], [3/4, 5/4, 5/4], -1) for mine,
just by symmetry, but a 4F3 hypergeometric is not quite nice, even if it is not a series in
special functions (ok, it may be written as such I think).

Herman Rubin suggests to use u=x and dv=sin(x)/(1+(cos(x))^4)dx with integration by parts to
solve Int(x*sin(x)/(1+(cos(x))^4),x=0..Pi/2).

This - more or less - with (real) partial fraction leads to similar questions/problems as we
looked at the thread http://www.mapleprimes.com/forum/mapleintegrationerror0, where the task
was integrating log ( polynom( cos(t) ) ) over 0 ... Pi.

I think the latter always cam be done in Maple, if the polynom has conjugated complex roots
(treating them by hand) and that should be the case for the problem here.

However that seems not to work for degree=1, i.e. for Int(ln(cos(t)-a),t = 0 .. 1/2*Pi) and
'my' eval(%, a = (1 - I ) / sqrt(2)).

The translation is as follows:
  
  Int( arccos(x)/(x - a), x = 0 .. 1 ) =                 # x=arccos(t) 
= Int( t*sin(t)/(cos(t)-a), t = 0 .. 1/2*Pi ) =          # using Parts(%,t)
= Int( ln(cos(t)-a),t = 0 .. 1/2*Pi ) - 1/2*Pi*ln(-a)

For our a = (1 - I ) / sqrt(2)) Maple 11.2 returns the last integral unevaluated, numerical
it is -0.41481325995246+2.5875517785926*I.

However one *can* do it symbolical as below, where it is nicer to have it up to Pi (not Pi/2)
where M11 also would return unevaluated:

  Int(ln((cos(t)-a)),t=0..Pi);
  ``=PDEtools[dchange](t=arccos(x),%,[x]);
  ``=value(rhs(%)):
  simplify(%); 
  eval(%, a = (1 - I ) / sqrt(2)); evalf(%);

                          Pi
                         /
                        |
                        |    ln(cos(t) - a) dt
                        |
                       /
                         0

                              1
                             /
                            |    ln(x - a)
                         =  |   ----------- dx
                            |         2 1/2
                           /    (1 - x )
                             -1

                   /                        /      2\1/2 \
                   |                        |-1 + a |    |
              = Pi |ln(-a) - ln(2) + ln(1 + |-------|   )|
                   |                        |   2   |    |
                   \                        \  a    /    /


                 = 0.22348749526260 + 6.7313498260590 I

That is the same value as a pure numerical quadrature gives.

Doing it for Pi/2 would show hypergeometrics 3F2, that just simplified for Pi as upper bound.

Now as we have it for 1/linear one can get the inital problem (have not carried it out) using
partial fraction decomposition.

i would try Tsat_p@8 as function name, not Tsat_p  and try to find a "Dependency walker" on the web to look at your DLL (which should show you th exact name of the function in the DLL)

First 179 180 181 182 183 184 185 Last Page 181 of 209