5^(1/11)*2^(1/22)*exp(5/44*I*Pi); # wrong answer

                   (1/11)  (1/22)    /5      \
                  5       2       exp|-- I Pi|
                                     \44     /

evalf(%);

                  1.119288902 + 0.4174731763 I

5^(1/11)*2^(1/22)*(cos(5*Pi/44)+I*cos(17*Pi/44)); # wrong answer

           (1/11)  (1/22) /   /5    \        /17   \\
          5       2       |cos|-- Pi| + I cos|-- Pi||
                          \   \44   /        \44   //

evalf(%);

                  1.119288902 + 0.4174731760 I

5^(1/11)*2^(1/22)*exp(-3/44*I*Pi);

                  (1/11)  (1/22)    / 3      \
                 5       2       exp|--- I Pi|
                                    \ 44     /

evalf(%);

                  1.167308782 - 0.2539324823 I

5^(1/11)*2^(1/22)*exp(5/44*I*Pi); # wrong answer

                   (1/11)  (1/22)    /5      \
                  5       2       exp|-- I Pi|
                                     \44     /

evalf(%);

                  1.119288902 + 0.4174731763 I

5^(1/11)*2^(1/22)*(cos(5*Pi/44)+I*cos(17*Pi/44)); # wrong answer

           (1/11)  (1/22) /   /5    \        /17   \\
          5       2       |cos|-- Pi| + I cos|-- Pi||
                          \   \44   /        \44   //

evalf(%);

                  1.119288902 + 0.4174731760 I

5^(1/11)*2^(1/22)*exp(-3/44*I*Pi);

                  (1/11)  (1/22)    / 3      \
                 5       2       exp|--- I Pi|
                                    \ 44     /

evalf(%);

                  1.167308782 - 0.2539324823 I

@mgu 

Maple 15.01

restart;

assume(n::posint);

dsolve({diff(y(x),x,x)+diff(y(x),x)/x-y(x)*(n^2)/(x^2)=-4*(n^2)/(x^2),
        y(1)=0,y(0)=y(0)}) assuming y(0)::finite;

                                       n    
                            y(x) = -4 x  + 4

@mgu 

Maple 15.01

restart;

assume(n::posint);

dsolve({diff(y(x),x,x)+diff(y(x),x)/x-y(x)*(n^2)/(x^2)=-4*(n^2)/(x^2),
        y(1)=0,y(0)=y(0)}) assuming y(0)::finite;

                                       n    
                            y(x) = -4 x  + 4

Could you post lprint'd (plaintext) code for your three inequalities, and the ranges, in a Comment to this Question?

@mgu 

Please reread the Answer and followup comments more carefully.

I explained that the following would not produce the result (or, indeed, any result at all):

restart;
assume(K<infinity);
dsolve({diff(y(x),x,x)+(diff(y(x),x))/x+8,y(1)=0,y(0)=K});

And my Answer showed that the following would produce your desired result (in Maple 15.01) for your form:

restart:
dsolve({x*ln(x)*diff(y(x),x)-y(x)=2*x^2-4*x^2*ln(x), y(0)=y(0)},
        y(x)) assuming y(0)<infinity;

And, with that same approach:

restart;
dsolve({diff(y(x),x,x)+(diff(y(x),x))/x+8,y(1)=0,y(0)=y(0)}) assuming y(0)<infinity;

                                   2    
                        y(x) = -2 x  + 2

@mgu 

Please reread the Answer and followup comments more carefully.

I explained that the following would not produce the result (or, indeed, any result at all):

restart;
assume(K<infinity);
dsolve({diff(y(x),x,x)+(diff(y(x),x))/x+8,y(1)=0,y(0)=K});

And my Answer showed that the following would produce your desired result (in Maple 15.01) for your form:

restart:
dsolve({x*ln(x)*diff(y(x),x)-y(x)=2*x^2-4*x^2*ln(x), y(0)=y(0)},
        y(x)) assuming y(0)<infinity;

And, with that same approach:

restart;
dsolve({diff(y(x),x,x)+(diff(y(x),x))/x+8,y(1)=0,y(0)=y(0)}) assuming y(0)<infinity;

                                   2    
                        y(x) = -2 x  + 2

@Markiyan Hirnyk I was surprised that it worked in the first place, since these both produced NULL.

dsolve({ODE, y(0)=K},y(x)) assuming K<infinity;

assume(K<infinity);
dsolve({ODE, y(0)=K},y(x));

@Markiyan Hirnyk I was surprised that it worked in the first place, since these both produced NULL.

dsolve({ODE, y(0)=K},y(x)) assuming K<infinity;

assume(K<infinity);
dsolve({ODE, y(0)=K},y(x));

The performance of these various methods vary (differently from each other) with both the number of repetitions and the number of elements of the list.

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^5)]:
str:=time[real]():
for j to 10 do nops(remove(c->0<c,thelist)):end do:
time[real]()-str;
                             0.651

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^5)]:
str:=time[real]():
for j to 10 do nops(select(`<`,thelist,0)):end do:
time[real]()-str;
                             0.446

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^5)]:
str:=time[real]():
for j to 10 do add(`if`(x<0,1,0),x=thelist):end do:
time[real]()-str;
                             0.480

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^5)]:
str:=time[real]():
for j to 10 do Threads:-Add(`if`(x<0,1,0),x=thelist):end do:
time[real]()-str;
                             0.255

That was in 64bit Maple 15.01 in Windows 7 on a quad-core i7.

Now, increasing the length

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^7)]:
str:=time[real]():
for j to 1 do nops(remove(c->0<c,thelist)):end do:
time[real]()-str;
                             7.216

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^7)]:
str:=time[real]():
for j to 1 do nops(select(`<`,thelist,0)):end do:
time[real]()-str;
                             5.276

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^7)]:
str:=time[real]():
for j to 1 do add(`if`(x<0,1,0),x=thelist):end do:
time[real]()-str;
                             5.024

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^7)]:
str:=time[real]():
for j to 1 do Threads:-Add(`if`(x<0,1,0),x=thelist):end do:
time[real]()-str;
                             4.819

Interestingly, Threads:-Add does more poorly as the number of repetitions gets high.

The performance of these various methods vary (differently from each other) with both the number of repetitions and the number of elements of the list.

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^5)]:
str:=time[real]():
for j to 10 do nops(remove(c->0<c,thelist)):end do:
time[real]()-str;
                             0.651

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^5)]:
str:=time[real]():
for j to 10 do nops(select(`<`,thelist,0)):end do:
time[real]()-str;
                             0.446

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^5)]:
str:=time[real]():
for j to 10 do add(`if`(x<0,1,0),x=thelist):end do:
time[real]()-str;
                             0.480

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^5)]:
str:=time[real]():
for j to 10 do Threads:-Add(`if`(x<0,1,0),x=thelist):end do:
time[real]()-str;
                             0.255

That was in 64bit Maple 15.01 in Windows 7 on a quad-core i7.

Now, increasing the length

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^7)]:
str:=time[real]():
for j to 1 do nops(remove(c->0<c,thelist)):end do:
time[real]()-str;
                             7.216

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^7)]:
str:=time[real]():
for j to 1 do nops(select(`<`,thelist,0)):end do:
time[real]()-str;
                             5.276

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^7)]:
str:=time[real]():
for j to 1 do add(`if`(x<0,1,0),x=thelist):end do:
time[real]()-str;
                             5.024

restart:
gen:=rand(-100..100):
thelist:=[seq(gen(),i=1..10^7)]:
str:=time[real]():
for j to 1 do Threads:-Add(`if`(x<0,1,0),x=thelist):end do:
time[real]()-str;
                             4.819

Interestingly, Threads:-Add does more poorly as the number of repetitions gets high.

Sorry, but I can't tell if you want to integrate independently over each component of some Vectors (in which case you could just `map` the integration command) or whether you are looking to make your multiple integrals more terse.

If the latter then note that both exact and numeric (nested) multiple integrals can be done with a shorthand syntax:

P:=Int(Int(Int(x*y^2-x*y*z, x=0..Pi), y=0..Pi/4), z=0..2*Pi);
evalf(P);
value(P);
P:=Int(x*y^2-x*y*z, x=0..Pi, y=0..Pi/4, z=0..2*Pi);
evalf(P);
value(P);
P:=Int((x,y,z)->x*y^2-x*y*z, [0..Pi, 0..Pi/4, 0..2*Pi]);
evalf(P);
value(P);

@dman I believe that you are wrong in stating that the one of the three values returned by my applying my appoach to the modified equation (6, for 4) is "incorrect". It seems that there are at least four roots, and that my and your methods have each found only three of them (with 2 in common, found by both).

However, by applying `simplify` to the RHS -the LHS of the log[10] transformed equation, a missing root is found (and another dropped). So, four distinct roots can be computed programmatically, by taking the union of both attempts.

restart:
eq:=tan(6*(1.1575*10^12/(1-.5)^1.2-x^2)^(1/2)*10^(-6)*(1/2))
    /tan(6*(5.555*10^11/(1-.5)^1.2-x^2)^(1/2)*10^(-6)*(1/2))
  = -4*x^2*(5.555*10^11/(1-.5)^1.2-x^2)^(1/2)*(1.1575*10^12
    /(1-.5)^1.2-x^2)^(1/2)/(1.1575*10^12/(1-.5)^1.2-2*x^2)^2:

neweq:=subs(x=10^(y),eq):
Digits:=15:

newans:=Student:-Calculus1:-Roots(neweq, y = log[10](1*10^5) .. log[10](9*10^6));

                [5.91575925330614, 6.21237849725823, 6.26387733777593]

newans2:=Student:-Calculus1:-Roots(simplify((rhs-lhs)(neweq)),
                                   y = log[10](1*10^5) .. log[10](9*10^6));

                [5.91575925330614, 6.05300347077205, 6.26387733777594]

all:=[op(evalf[14]({op(newans), op(newans2)}))];

          [5.9157592533061, 6.0530034707720, 6.2123784972582, 6.2638773377759]

Digits:=100:
allmod:=map(t->10^t,all):
Digits:=15:
newestmod:=evalf(allmod);

                      [                   5                     6  
                      [8.23681389099905 10 , 1.12980494375297 10 , 

                                        6                     6]
                     1.63071661915845 10 , 1.83601970432110 10 ]

Digits:=100:
check:=map(t->eval((rhs-lhs)(eq),x=t),allmod):
Digits:=15:
evalf(check);

           [                   -12                     -9                       -7    
           [2.10833371198086 10   , 9.22936651431375 10   I, 8.05724958238040 10   I, 

                                 -13]
             -1.09761159667640 10   ]

This is becoming more ad hoc.

@dman I believe that you are wrong in stating that the one of the three values returned by my applying my appoach to the modified equation (6, for 4) is "incorrect". It seems that there are at least four roots, and that my and your methods have each found only three of them (with 2 in common, found by both).

However, by applying `simplify` to the RHS -the LHS of the log[10] transformed equation, a missing root is found (and another dropped). So, four distinct roots can be computed programmatically, by taking the union of both attempts.

restart:
eq:=tan(6*(1.1575*10^12/(1-.5)^1.2-x^2)^(1/2)*10^(-6)*(1/2))
    /tan(6*(5.555*10^11/(1-.5)^1.2-x^2)^(1/2)*10^(-6)*(1/2))
  = -4*x^2*(5.555*10^11/(1-.5)^1.2-x^2)^(1/2)*(1.1575*10^12
    /(1-.5)^1.2-x^2)^(1/2)/(1.1575*10^12/(1-.5)^1.2-2*x^2)^2:

neweq:=subs(x=10^(y),eq):
Digits:=15:

newans:=Student:-Calculus1:-Roots(neweq, y = log[10](1*10^5) .. log[10](9*10^6));

                [5.91575925330614, 6.21237849725823, 6.26387733777593]

newans2:=Student:-Calculus1:-Roots(simplify((rhs-lhs)(neweq)),
                                   y = log[10](1*10^5) .. log[10](9*10^6));

                [5.91575925330614, 6.05300347077205, 6.26387733777594]

all:=[op(evalf[14]({op(newans), op(newans2)}))];

          [5.9157592533061, 6.0530034707720, 6.2123784972582, 6.2638773377759]

Digits:=100:
allmod:=map(t->10^t,all):
Digits:=15:
newestmod:=evalf(allmod);

                      [                   5                     6  
                      [8.23681389099905 10 , 1.12980494375297 10 , 

                                        6                     6]
                     1.63071661915845 10 , 1.83601970432110 10 ]

Digits:=100:
check:=map(t->eval((rhs-lhs)(eq),x=t),allmod):
Digits:=15:
evalf(check);

           [                   -12                     -9                       -7    
           [2.10833371198086 10   , 9.22936651431375 10   I, 8.05724958238040 10   I, 

                                 -13]
             -1.09761159667640 10   ]

This is becoming more ad hoc.

Preben's observation that those numeric integrals don't depend upon `t`, and can be computed just once up front (instead of repeatedly, for each step) is very good.

As a minor comment, that would also work with the original mode of using dsolve/numeric.

restart:
st:=time():
# Maple15\xsys11demo.mw, standard worksheet, BoltzQM, Shift & enter, ADD not SUM, 1-D math for procedure
with(StringTools):
FormatTime(%c);
with(DEtools):with(plots):
printlevel := 1:
mm := 3: nn := 3: a := 5.: e := 1.:  R := .1: h := 1.:

for n to nn do
for nd to nn do
for md to mm-1 do
integ[n,nd,md]:=evalf(Int( sin(Pi*n*x/a)*sin(Pi*nd*x/a)*sin(Pi*n*y/a)
                           *sin(Pi*nd*y/a)/sqrt((x-y)^2+4*R^2*sin(Pi*md/mm)^2),
                           y = 0. .. a, x = 0. .. a))
end do
end do
end do;

dproc3 := proc (nn, t, Y, YP) local n, nd, md;
  for n to nn do
    YP[n] := -((I*e^2*4)/(h*a^2))*Y[n]
             *add(add(Y[nd]*conjugate(Y[nd])*integ[n,nd,md],
                      md = 1 .. mm-1), nd = 1 .. nn);
  end do;
end proc:

ic3 := array([seq(1/sqrt(nn), n = 1 .. nn)]):
dvars3 := [seq(c[n](t), n = 1 .. nn)]:

p := dsolve(type=numeric, number = nn, procedure = dproc3, start = 0.,
            initial = ic3, procvars = dvars3, range=0..0.1, method=rkf45):

FormatTime(%c);
odeplot(p,[seq([t,abs(c[n](t))],n=1 .. nn)]);
time()-st;