Axel Vogt

5408 Reputation

20 Badges

17 years, 230 days
Munich, Germany

MaplePrimes Activity


These are answers submitted by Axel Vogt

I would proceed similar, but a bit different (except you only care for the 'length' - but why?).

Usually one expects such a fraction to be 'exact' for the precision one needs.

Thus I would do

v := 1433162862156581221897903/4340366287316864960103999;
oldDigits:= Digits:
  Digits:=10:        # fall down to desired precision
  evalf(v):
  vNew:=convert(%, rational);
Digits:=oldDigits:   # restore settings


                         1433162862156581221897903
                    v := -------------------------
                         4340366287316864960103999


                                    15251
                            vNew := -----
                                    46188

In this example vNew = v is correct up to at most the last decimal place
in a precision of 10 decimal places.

# symmetric in x and then integrand has an inverse
ann:= 0<x,0<q,1<y;

Int(sech(x)^q, x);
IntegrationTools:-Change(%, x=arccosh(y), y);
combine(%) assuming ann;
simplify(%) assuming ann;
value(%) assuming ann;
A:=eval(%, y = cosh(x)) assuming ann;

 

# result:

-1/q*cosh(x)^(-q)*hypergeom([1/2, 1/2*q],[1+1/2*q],1/(cosh(x)^2));

Your f is abs(g) and g is positive over the range 1 ... 9. So just write Int(g, p= 1 .. 9), which numerically evaluates quickly (and even faster if you use 'simplify' before feeding it to evalf).


 

#VectorCalculus:-int(sin(x^2)*cos(y^2), [x, y] = Circle(`<,>`(0, 0), 1), inert);

sin(x^2*cos(y)^2)*cos(x^2*sin(y)^2)*x; #int(%, y);
4*Int(%, y=0..Pi/2);
IntegrationTools:-Change(%, y = arcsin(t), t); value(%);
Int(%, x=0..1);
value(%);
evalf(%);

sin(x^2*cos(y)^2)*cos(x^2*sin(y)^2)*x

 

4*(Int(sin(x^2*cos(y)^2)*cos(x^2*sin(y)^2)*x, y = 0 .. (1/2)*Pi))

 

-4*x*(Int((sin(x^2*t^2)*cos(-x^2)+cos(x^2*t^2)*sin(-x^2))*cos(x^2*t^2)/(-t^2+1)^(1/2), t = 0 .. 1))

 

x*sin(x^2)*Pi

 

Int(x*sin(x^2)*Pi, x = 0 .. 1)

 

(1/2)*Pi-(1/2)*cos(1)*Pi

 

.722091449378415

(1)

 


 

Download MP-220842.mw

Shift it to 0, find an ODE for it, solve it and determine parameters:

task:=hypergeom([1/3, 2/3],[3/2],27/4*z^2*(1-z));
h:=eval(%, z=1+x);

# find differential equation
f(x)=h;
PDEtools[dpolyform](%, no_Fn):
ode:=op(1, %)[1];
dsolve(ode):
H:=rhs(%);

# check ode
odetest(f(x)=H,ode);
odetest(f(x)=h,ode);

# find parameters
'eval(h, x=1)'; convert(%, StandardFunctions): evala(%):
c1:=%=eval(H, x=1); #convert(%, StandardFunctions): evala(%);

'eval(h, x=2)'; convert(%, StandardFunctions): evala(%):
c2:=%=eval(H, x=2);

param:=solve({c1,c2});

# re-substitute
eval(H, param);
eval(%, x=z-1);

 

Occasionally I do the following

f := x -> piecewise(0 < x,-1,x <= 0,1): # double point = quiet
'f(x)': '%'=%;                          # a simple way to display such things
                             { -1        0 < x
                      f(x) = {
                             { 1         x <= 0

 

One can read it as a(i+1) = i^2 - a(i), which is a recursion.
Maple can solve this, i.e. the term a(k) is given "directly":

  restart;
  rsolve(a(i)+a(i+1)=i^2, a(k)); # see the help for the command

                        k
               a(0) (-1)  - 1 + (k + 1) (k/2 + 1) - 2 k

 


For specific values one has to provide a(0), of course.

In "standard interface" one can insert a table (through the menu bar) and each cell can have execution groups to display the plots.

If you dislike to show the commands the you can hide them (via menu/show) - but (for me) that only works for the entire worksheet.

The correct way ist

x^(5/3)-5*x^(2/3);
[Re(%), Im(%)];
plot(%, x=-3..7, color=[red, blue]);

Instead of using floating point numbers, like 0.3, one can use decimal fraction, like 3/10.

Then the DE can be solved symbolically. Plotting indicates a zero at r ~ 1.5*10^5  1.5*10^6 and that complex values might occure.

Analyzing and using my own zero finder (like bracketting) gives a zero between 1337252.92663784 and 1337252.92663785


Be aware that this depends on the exactness of your input parameters (i.e. how many decimal places are used) and how many Digits are used for computations.

I used Digits:=15

 

Int(1/u*t*exp(-t/u), t = 0 .. infinity);
value(%) assuming 0<u;

Numerically one splits in Pi/2 * k, k integer, because of needed 'smoothness' ~ subdivisions. One would need k ~ 10^6 for double precision, I think.

Manual support gives Sum(2*Si(Pi*(2*k+1))-Si(2*k*Pi)-Si((2*k+2)*Pi)+(Si(k*Pi)-Si((k+1)*Pi))*(-1)^k,k = 0 .. infinity) and that numerically gives the expected result.

Plot the anti derivative and remember the very theorem you want to use.

Int(sqrt(exp(2*I*t)-1),t);
J:=value(%);
[Re(%), Im(%)];
plot(%, t=0 .. Pi, color=[red, blue]);

You get it through the bivariate standard normal distribution N2 using
a transformation of variables:

restart;

f := proc (t1, t2, alpha1, beta1, alpha2, beta2, rho) options operator, arrow;
  (1/4)*(sqrt(beta1/t1)+(beta1/t1)^(3/2))*(sqrt(beta2/t2)+(beta2/t2)^(3/2))*
  exp(-((sqrt(t1/beta1)-sqrt(beta1/t1))^2/alpha1^2+
  (sqrt(t2/beta2)-sqrt(beta2/t2))^2/alpha2^2-2*rho*(sqrt(t1/beta1)-
  sqrt(beta1/t1))*(sqrt(t2/beta2)-sqrt(beta2/t2))/(alpha1*alpha2))/(2-2*rho^2))/
  (alpha1*beta1*alpha2*beta2*Pi*sqrt(1-rho^2)) end proc;

N2:= (x,y,rho) -> 1/sqrt(1-rho^2)/2/Pi*
  Int(Int(exp(-(xi^2-2*rho*xi*eta+eta^2)/(2*(1-rho^2))),
                    eta=-infinity..y), xi=-infinity..x);

trans := {x1 = (sqrt(t1/beta1)-sqrt(beta1/t1))/alpha1,
          x2 = (sqrt(t2/beta2)-sqrt(beta2/t2))/alpha2};

assume(-1 < rho and rho < 1,
       alpha1 > 0, beta1 > 0, alpha2 > 0, beta2 > 0,
       t1 > 0, t2 > 0);

'N2(x1,x2,rho)': 'eval(%, trans)':
'diff(%, t1)': 'diff(%, t2)':;
'1/2*f(t1, t2, alpha1, beta1, alpha2, beta2, rho) = %';
 
simplify(%):
is(%);

  1/2 f(t1, t2, alpha1, beta1, alpha2, beta2, rho) =

           2
          d
        ------- eval(N2(x1, x2, rho), trans)
        dt2 dt1


                                 true


The task can be seen as to find Beta(1/n,n-1/n) = n, a fixed point. Using the integral theorem for the Beta function it means that the integrals over t = 0 ... 1 for  t^(1/n-1)*(1-t)^(n-1/n-1) and t^(1/n-1) are equal.

This certainly is the case if it holds true for the integrands. For that 'solve' includes the desired formal solution.

5 6 7 8 9 10 11 Last Page 7 of 89