Axel Vogt

5936 Reputation

20 Badges

20 years, 252 days
Munich, Bavaria, Germany

MaplePrimes Activity


These are answers submitted by Axel Vogt

It is not quite clear for me, what you mean and do, but it seems you think of T as 't is a function of x'.

However that is already part of a correct definition for your (1) and (2).

?

why not save or export as Maple code and read it in?

or deleting all output and save as text file?

if you are sure that's unique (which is not clear), then you can try "select(has, ...)" and work with that or try the "pattern" command

It may be that your expression is already in a quite compact form and best would be to work on your eq1 and eq2. While it is not understandable why combining it into one monster is desirable, having 3 parts it is a little more readable (if it is at all).

If you only want it to have something shorter to type for Latex: there is not need at all, just say latex(A2).

  length(sol); # your solution
                                 5827


  simplify(sol, trig):
  combine(%, exp) assuming l::real:
  simplify(%, size) assuming l::real;
  T:=Tryhard(%);
  length(T);

                                 839

This is not quite perfect, but the obvious rest for sinh should be easy.

For 'Tryhard' you may use www.mapleprimes.com/blog/axelvogt/simplifyingthroughcodegenoptimizetryhard

PS: I avoid to use the letter l, as it is hard to distinguish from 1, the number,
so typos are not easy to find (but that is a matter of taste)
you may try [map(CodeGeneration[Matlab],J)]; 
I played with it (no, not with the mathematical part, just with the
Numerics, even if the Math would be more interesting - but too hard
for me (have not found Schanuel's conjecture for trigonometrics, I
guess that might be a way).

Since sin(x)/cos(x) = tan(x) is better to solve sin(x) = cos(x)*x and
using asympotics a choice is x = n*Pi + 1/2*Pi - epsilon, n an integer.

  sin(x) - cos(x)*x:
  'eval(%, x=n*Pi + Pi/2 - epsilon)';
  expand(%): 
  simplify(%) assuming n::posint:
  collect(%, sin):
  EQ:=simplify(%, size);

    EQ := (-1)^(1+n)*(((n+1/2)*Pi-epsilon)*sin(epsilon)-cos(epsilon))

The asymptotic series for solve(EQ, epsilon) is 1/(Pi*n) +- ... and
this way even for high values one finds a solution:

  nTst:=1+'`^`(10,10)'; # even or odd does not matter

  eval(EQ, n=nTst);
  epsilon0:=fsolve(%, epsilon=1/Pi/nTst, fulldigits);
  #epsilon0:= convert(%, rational);

Check it (I work with Digits = 14):

  eval(sin(x) - cos(x)*x, x=nTst*Pi + Pi/2 - epsilon0); evalf(%);

                                  0.


Doing the same with 'RootOf' instead of 'fsolve' gives a formal solution,
which (of course) Maple agrees to be correct:

  eval(EQ, n=nTst);
  epsilon0:=RootOf(%, epsilon);
  eval(sin(x) - cos(x)*x, x=nTst*Pi + Pi/2 - epsilon0); simplify(%);

                                  0

and if one wants to know it, then

  'evalf[100](epsilon0)': '%'=%;
                                             -18
     evalf[100](epsilon0) = 0.3183 ... 983 10


To compute or verify it directly for tan(x)=x is not my thing ...

I work with The Reasonable Interface AKA as the classical one, on Windows

Marking the section, ctrl + C and ctrl + V does is, where then I choose format = formatted in the editor here from the menu bar. This also pastes the usual '> ' for inputs (and if I do not want that, then I just remove it before ina word processor)

If marking the output as well at the time, then I usually get what I want (may be removing some outputs).

That's all (for me).

Salvo, at a blackboard it would be better to follow you, but there is none.

My 'feeling' says that should written a bit different, but I do not really
recognize what you intend. If it is to apply the integrator to a sum or so
then would choose a different way.

Perhaps you can write it down not as code, but as task with some explanations:
"expression(s) = ... now want Int(..) and sum it" or so (extrema = bounds
for the integrals(s) ?).


May be it helps if you see, that concerning 2 cases for the integrand is
not needed: if g=1, then 'x -> g(x)*K(x,a)' = 'x -> K(x,a)', since g(x)=1.
Just try it by defining gg:= x -> 1 and use it (gg:=1 also would do it).

Salvo, That fine you cleared that and hope you can nail down the rest.

PS. The routine has a typo: at the beginning Pi/2 is computed as
pi2 := evalhf(1/2*Pi) and that has to be pi2 := evalf(1/2*Pi) to
have that with the correct precision.

is 'if eval(A, a=0) = A then ...' enough for your needs?
As far as I understand you only have a handling problem Maple's syntax and
you want some of the following as additional factor for K(x,a):

  [1, ln(f), (ln(f))^2, 1/(1+x^2)]; 
  eval(%, f=t0*x);

                                         2    
                   [1, ln(t0 x), ln(t0 x) , 1/(1+x^2)]

Just define them as functions, then it is easy to write it down:

  myInt2_DE:=proc(g, a,b,c, DigitsForWorking)
  local x, f0;
  Digits:=DigitsForWorking;
  f0:=215;
  f0^(1+a)*intDE('x -> g(x)*K(x,a)',  f0/c,f0/b, Digits)
  end proc;

  gTst:= unapply(ln(f0*x), x);
  #gTst:= x -> 1/(1+x^2);
  myInt2_DE(gTst, 1,20,2500, 32);

                  11746.667173055677376120240588752


A variant of the integration routine (similar to what is possible in Maple)
may be to use high Digits, but setting the desired relative error 'eps' 
much lower. In this way one can care for possible peaks or ugly situations
without enforcing to compute with a very small relative, which actually
may not be needed later (i.e. one can try to save computational time).
But such 'tuning' should be done after all is seen to be ok.

Salvo, no problem, over all the years I get used to that, but only almost :-)

Salvo,

double exponential integration is not my method, it is due to Mori and known to work well even for higher precision. All I did was translating/adapting existing code to Maple, as stated.

Roughly it computes Int(f(x), x = -1 .. 1) by transforming it through x = tanh( Pi/2 * sinh(t) ) to the full real line and using the trapez rule there.

PS: Alex is not quite the same as Axel ...

I lost a bit the overview for the variants you have in mind, but hope the
attached sheet does it for your sub-task Int(t^a/S(t), t = b ... c).

A bit strange, that Maple's on-board-method double exponential integration
seems to be on strike for that. So I modified what I once had translated
for my needs, hope it helps.

For the jittering between levels: it is a 'zooming' or scaling problem and
does not stand for errors - it more or less says that the 'resulution' is
not fine enough (and one can be trapped by graphics as well). May be you
also consider the advice 'make it dimensionless' (but as a former Math guy
I am not used really to that recipe).

The last thing i tried (not today) was towards approximating the integrand
by a Hermite expansion. But your stuff as pdf (i.e. with a=0) is quite far
away from a normal distribution (even on that finite range and certainly
not over the Reals). But as 'zooming' is your problem one would need quite
a good one (even if symbolic solutions exist). So I did not: time for a beer.

Download 102_double_exponential_integration_with_varaible_digits.mws
View file details

Edited 29 Jan 2010

The routine has a typo: at the beginning Pi/2 is computed as
pi2 := evalhf(1/2*Pi) and that has to be pi2 := evalf(1/2*Pi) to
have that with the correct precision.

First 56 57 58 59 60 61 62 Last Page 58 of 93