Actually the int in the logarithm was the issue... Int doesnt crash...

Ah sry u mean coz its inside the log?

@Carl Love sry about the thread safe aspect.. I was just thinking maybe it is trying to write sth while still calculating...That is actually what i meant but probably I'm wrong...

for the int part...since I'm using evalf(int()) anway is Int really needed???

@tomleslie Ok: of the order 1-10 means 1.0 - 10.0 not 1e-10

So for z=1 do u have sth like 5*10^53

coz that would be wrong I suppose or?

@tomleslie 1. if you set z=1 instead of z=1e-3 then I'm at 5*10^53 in my plot using ur modified version:

sol:-plot(z = 1, numpoints = 500, color = blue);

2. Sure I know that. What I meant was this:

Y(eta,z)=Y0(z)*(1-eta^2)

Y0(z) is decreasing as i move along the z-axis...the radial part is here for simplicity written as (1-eta^2). Qualitatively that is was I expect somehow just that the amplitude decreases as I move along the z-axis.

This is just sth like a heat equation. At the beginning the temperature is constant e.g. 10*T0 everywhere and with time (dirichlet boundary conditions T(boundary)=T0) it will decay in a visual melting fashion towards T0 like 1-r^2...

Just that time here is z...and there are no oscillations...

3. But that is just the problem: this for fixed eta>0 and eta<1 should go to zero as Y(eta,z)=sum(c_n*besselj(0,alpha_n*eta)*exp(-lambda_n*z),n=1..infinity);

shows for constant vz.

@tomleslie Hey thx for the reply. The issue is that the solution is totaly not what I expect. In particular I dont expect sth oscilatory which I also saw, but the numbers where somewhere in the range 10^50

That's why I said divergent. I do expect the result to be 1 at z=0 and then slightly decreasing as u move along the z-axis...since I said it should vanish at eta=1, qualitatively it should look sth like

Y(eta,z)=Y0(z)*(1-eta^2) where Y0(z) is decreasing for increasing z.

But for sure not oscilating...

How large are the oscilating numbers? 10^50 too?

edit: even though I put z=1e-3 in the initial example it should be more of the order 1-10 and thats where it becomes too large to be believable...

edit2: You can just check the result for vz=const. The solution can be given analytically and the result looks completely different: i.e. it is Y(eta,z)=sum(c_n*besselj(0,alpha_n*eta)*exp(-lambda_n*z),n=1..infinity);
lambda_n is a function of the roots of the bessel function alpha_n and >0

This is what I expect from the solution here too... to decay exponentially

As a sidenote: I'm not quite sure about the fourth initial conditions, (D[2](Y))(eta, 0) = 0].

I know the profile at z=0 but since the equation is second order in z i obviously need a second one and I thought this one would be reasonable. Y(eta,infinity)=0 is not feasible I guess...

Anyway in general this shouldnt make the solution diverge, right?

Ok I found a solution sry

though I still dont understand what mo,mi,mn,mtext...

means

mo() somehow stretches something mi() shortens it mtext something like text but how u actually use it I dont know :-/

@Kitonum  THX

@Kitonum yeah thanks. I did have a procedure but I wasn't sure whether there might be some already built in function since it's quite elementary. BTW I used convert([...],`*`) since `*`([...]) doesn't work for me somehow...

is it the same or internally a bit different?

@Preben Alsholm Hey Thanks. Do you also no an answer to my other question?

the help says D and diff are similar. So in particular

D(f)(x) = diff(f(x),x) but why does not D(x)(x) equal diff(x,x) ??

Obviously in the first case maple thinks of x being some function of x but this is not the case...

Maybe some related question:

Why does maple not recognize in

((x*D[1])@@2)(F)(x,M)=x*((D[1](x))(x, M)*(D[1](F))(x, M)+x*(D[1, 1](F))(x, M))

the (D[1])(x)(x, M) as 1 ??

What is maple thinking?

@Carl Love ok thank you. I have found an answer to the deleted question!

@Carl Love Hey THANKS!

Do you know how it's possible to plot for the infinite product??

e.g.

S:=QPochhammer(q,q^5,infinity)
plot(S,q=-1..1)

Warning, unable to evaluate the function to numeric values in the region; see the plotting command's help page to ensure the calling sequence is correct

In general is there a way to work with expressions like

QPochhammer(-q,-q^5) since naturally I don't see a way around the minus sign in the second argument.

@Rouben Rostamian  Unfortunately I do not have such an old version but thanks ;)