dharr

Dr. David Harrington

7109 Reputation

21 Badges

20 years, 148 days
University of Victoria
Professor or university staff
Victoria, British Columbia, Canada

Social Networks and Content at Maplesoft.com

Maple Application Center
I am a retired professor of chemistry at the University of Victoria, BC, Canada. My research areas are electrochemistry and surface science. I have been a user of Maple since about 1990.

MaplePrimes Activity


These are replies submitted by dharr

@MaPal93 I wasn't suggesting the different results were pathological, rather that there was some minor error somewhere. A pathological case might occur for example if you want lim(a(x)*b(x),x=infinity). You evaluate lim(a(x),x= infinity and find it is infinity. Or you evaluate limit(b(x),x=infinity) and find it is zero. You cannot conclude anything from this, because we could have a=x^2 and b=1/x^2 and the true limit is 1. But if both limits are finite and non-zero, you should be OK.

I redid the calculation and simplification steps, checking numerically at each step. The result is different from above. To find out why, you could apply the same numerical checks to the other procedures to track done what is probably some small error in one step. (Edit: I see you figured this out.) But W4 is achieved without any limiting process, and is verifiably correct. If you rerun with gamma = 1e8, you see that W4 and the limiting values agree. So I am confident this version is correct.

limits.mw

@MaPal93 I wanted to check the limit for the RootOf, since there was some error with an earlier example that I submitted an SCR for, but at least for lambda it looks OK. You could pursue the sort of numerical tests I do here to resolve the discrepancy.

same3.mw

@MaPal93 I can imagine some pathological cases, but in general I would expect that the results should be the same. If I use limit(X[i], gamma = infinity) in generic_n-form_of_infinity_limit2_MaPal.mw I don't get the expression for X[i] in  n_from_2_to_6_-_not_coincide.mw so perhaps something is incorrectly entered, or perhaps I am not understanding what you are doing.

@MaPal93 I rewrote with as few summations as possible and was able to get the following simple result:

generic_n-form_of_infinity_limit2.mw

The _a indices arise when it would otherwise lead to a summation over i inside another summation over i, which is confusing since they are only dummy indices and not related. So Maple makes sure they have different names.

@MaPal93 I don't have any specific suggestions other than try to figure out the patterns.

restart;

coeff_term := (sqrt(3)*sigma__v*w/(4*sqrt(1 + 2*rho)*sigma__dr) + sqrt(3)*sigma__v*w*rho/(2*sqrt(1 + 2*rho)*sigma__dr))*delta__r^2 + sigma__v*delta[2]*delta[1]*rho/(9*sigma__d*sqrt(1 + 2*rho)) + sigma__v*delta[1]^2/(36*sigma__d*sqrt(1 + 2*rho)) + sigma__v*delta[2]^2/(36*sigma__d*sqrt(1 + 2*rho)) + sigma__v*delta[3]^2/(36*sigma__d*sqrt(1 + 2*rho)) + sigma__v*delta[1]*delta[3]*rho/(9*sigma__d*sqrt(1 + 2*rho)) + sigma__v*delta[2]*delta[3]*rho/(9*sigma__d*sqrt(1 + 2*rho)) + sigma__v*delta[1]^2*rho/(18*sigma__d*sqrt(1 + 2*rho)) + sigma__v*delta[2]^2*rho/(18*sigma__d*sqrt(1 + 2*rho)) + sigma__v*delta[3]^2*rho/(18*sigma__d*sqrt(1 + 2*rho)) + sigma__v*delta[2]*delta[1]/(18*sigma__d*sqrt(1 + 2*rho)) + sigma__v*delta[1]*delta[3]/(18*sigma__d*sqrt(1 + 2*rho)) + sigma__v*delta[2]*delta[3]/(18*sigma__d*sqrt(1 + 2*rho));

((1/4)*3^(1/2)*sigma__v*w/((1+2*rho)^(1/2)*sigma__dr)+(1/2)*3^(1/2)*sigma__v*w*rho/((1+2*rho)^(1/2)*sigma__dr))*delta__r^2+(1/9)*sigma__v*delta[2]*delta[1]*rho/(sigma__d*(1+2*rho)^(1/2))+(1/36)*sigma__v*delta[1]^2/(sigma__d*(1+2*rho)^(1/2))+(1/36)*sigma__v*delta[2]^2/(sigma__d*(1+2*rho)^(1/2))+(1/36)*sigma__v*delta[3]^2/(sigma__d*(1+2*rho)^(1/2))+(1/9)*sigma__v*delta[1]*delta[3]*rho/(sigma__d*(1+2*rho)^(1/2))+(1/9)*sigma__v*delta[2]*delta[3]*rho/(sigma__d*(1+2*rho)^(1/2))+(1/18)*sigma__v*delta[1]^2*rho/(sigma__d*(1+2*rho)^(1/2))+(1/18)*sigma__v*delta[2]^2*rho/(sigma__d*(1+2*rho)^(1/2))+(1/18)*sigma__v*delta[3]^2*rho/(sigma__d*(1+2*rho)^(1/2))+(1/18)*sigma__v*delta[2]*delta[1]/(sigma__d*(1+2*rho)^(1/2))+(1/18)*sigma__v*delta[1]*delta[3]/(sigma__d*(1+2*rho)^(1/2))+(1/18)*sigma__v*delta[2]*delta[3]/(sigma__d*(1+2*rho)^(1/2))

simplify(coeff_term);

(1/2)*sigma__v*(3^(1/2)*w*delta__r^2*sigma__d+(1/9)*sigma__dr*(delta[3]+delta[1]+delta[2])^2)*(1/2+rho)/((1+2*rho)^(1/2)*sigma__dr*sigma__d)

 

NULL

Download patterns.mw

@MaPal93 Actually, I realized that if you use collect(W, gamma) (first method) you see W = (stuff1)*gamma+(stuff2), so the limit will be +infinity if stuff1 > 0 and -infinity if stuff1 < 0. So it seems correct that is infinite. (But now you have new X[i] etc.)

In general if you see something like signum(n*rho - rho + 1) you can just make some assumption like n*rho - rho + 1>0 to make it go away, assuming that such an assumption is valid.

@salim-barzani So if the input is 1/G(xi)*($$$$$$$$$$$$+diff(G(xi),xi)*(##########) what is the output you want?

I don't understand what you are trying to do, but this might be useful.

collect.mw

I don't think there is any simple way. Since you want several matrices, GenerateMatrix won't work directly. You can use dcoeffs to find the pieces and then assemble them into the matrices.

@janhardo The DLMF (and the corresponding hardback book version) is produced by NIST (National Institute of Standards and Technology) in the U,S., and is the successor to Abramowitz and Stegun's, "Handbook of Mathematical Functions"; it is considered an authoritative source, perhaps "the" authoritative source.

These differential equations are standard forms in the sense that they are often quoted, e.g., one speaks of "Bessel's equation"  and there is a common way to write it, though of course there is not universal consistency in these things. 

But in terms of classifications in a text as I think you want, there is probably not much of a standard. I personally like Zwillinger's "Handbook of Differential Equations". Of course these classifications exist because of different solving methods, but since the objective is often a special function, they are related to the other standard forms.

Maple' odeadvisor(verg) gives: [[_2nd_order, _missing_x], [_2nd_order, _reducible, _mu_x_y1]]

Since many of the second-order differential equations and their "special" functions are related through transformations to the Sturm-Liouville eqation, perhaps it is the standard form to rule them all :-).

@MichaelVio Here's an example.

restart

with(PDEtools):

Differential equation (=0)

det := diff(E(t),t) - A*E(t)/(exp(E(t)/k/T)-1);

diff(E(t), t)-A*E(t)/(exp(E(t)/(k*T))-1)

dchange({t=1/x},det);

-(diff(E(x), x))*x^2-A*E(x)/(exp(E(x)/(k*T))-1)

dchange({t=1/x},det/x^2,expand);

-(diff(E(x), x))-A*E(x)/(x^2*(exp(E(x)/(k*T))-1))

NULL

Download dchange.mw

@salim-barzani I'm not sure I fully understand the order you want to do things or what you want to store, but perhaps this:

restart

with(PDEtools); with(DEtools)

undeclare(prime)

`There is no more prime differentiation variable; all derivatives will be displayed as indexed functions`

c := combinat:-powerset({A, B}); Cases := [seq({seq(x = 0, `in`(x, s))}, `in`(s, c))]

{{}, {A}, {B}, {A, B}}

[{}, {A = 0}, {B = 0}, {A = 0, B = 0}]

Fode := (-delta*eta^2+alpha*eta)*(diff(diff(U(xi), xi), xi))-U(xi)*(2*eta*gamma*theta*(delta*eta-alpha)*U(xi)^2+eta^2*delta*k^2+(-alpha*k^2-2*delta*k)*eta+2*k*alpha+delta)

(-delta*eta^2+alpha*eta)*(diff(diff(U(xi), xi), xi))-U(xi)*(2*eta*gamma*theta*(delta*eta-alpha)*U(xi)^2+eta^2*delta*k^2+(-alpha*k^2-2*delta*k)*eta+2*k*alpha+delta)

case1 := {alpha = delta*(A^2*eta^2-2*eta^2*k^2-4*B*eta^2+4*eta*k-2)/(A^2*eta-2*eta*k^2-4*B*eta+4*k), eta = eta, a[0] = -A/(2*gamma*theta*RootOf(_Z^2*gamma*theta+1)), a[1] = RootOf(_Z^2*gamma*theta+1)}

{alpha = delta*(A^2*eta^2-2*eta^2*k^2-4*B*eta^2+4*eta*k-2)/(A^2*eta-2*eta*k^2-4*B*eta+4*k), eta = eta, a[0] = -(1/2)*A/(gamma*theta*RootOf(_Z^2*gamma*theta+1)), a[1] = RootOf(_Z^2*gamma*theta+1)}

K := U(xi) = a[0]+a[1]*G(xi)

U(xi) = a[0]+a[1]*G(xi)

S := diff(G(xi), xi) = G(xi)^2+A*G(xi)+B

diff(G(xi), xi) = G(xi)^2+A*G(xi)+B

"for i,case in Cases do    sol:=dsolve(eval(S,case));    case1AB:=eval(case1,case):     Usol:=eval(eval(K,{sol} union case1AB),case);   simplify(odetest(Usol,eval(Fode,case1AB)));  end do;  "

G(xi) = -(1/2)*A-(1/2)*tanh((1/2)*(A^2-4*B)^(1/2)*(c__1+xi))*(A^2-4*B)^(1/2)

{alpha = delta*(A^2*eta^2-2*eta^2*k^2-4*B*eta^2+4*eta*k-2)/(A^2*eta-2*eta*k^2-4*B*eta+4*k), eta = eta, a[0] = -(1/2)*A/(gamma*theta*RootOf(gamma*theta*_Z^2+1)), a[1] = RootOf(gamma*theta*_Z^2+1)}

U(xi) = -(1/2)*A/(gamma*theta*RootOf(gamma*theta*_Z^2+1))+RootOf(gamma*theta*_Z^2+1)*(-(1/2)*A-(1/2)*tanh((1/2)*(A^2-4*B)^(1/2)*(c__1+xi))*(A^2-4*B)^(1/2))

0

G(xi) = tan(B^(1/2)*(c__1+xi))*B^(1/2)

{alpha = delta*(-2*eta^2*k^2-4*B*eta^2+4*eta*k-2)/(-2*eta*k^2-4*B*eta+4*k), eta = eta, a[0] = 0, a[1] = RootOf(gamma*theta*_Z^2+1)}

U(xi) = RootOf(gamma*theta*_Z^2+1)*tan(B^(1/2)*(c__1+xi))*B^(1/2)

0

G(xi) = A/(-1+exp(-A*xi)*c__1*A)

{alpha = delta*(A^2*eta^2-2*eta^2*k^2+4*eta*k-2)/(A^2*eta-2*eta*k^2+4*k), eta = eta, a[0] = -(1/2)*A/(gamma*theta*RootOf(gamma*theta*_Z^2+1)), a[1] = RootOf(gamma*theta*_Z^2+1)}

U(xi) = -(1/2)*A/(gamma*theta*RootOf(gamma*theta*_Z^2+1))+RootOf(gamma*theta*_Z^2+1)*A/(-1+exp(-A*xi)*c__1*A)

0

G(xi) = 1/(c__1-xi)

{alpha = delta*(-2*eta^2*k^2+4*eta*k-2)/(-2*eta*k^2+4*k), eta = eta, a[0] = 0, a[1] = RootOf(gamma*theta*_Z^2+1)}

U(xi) = RootOf(gamma*theta*_Z^2+1)/(c__1-xi)

0

NULL

NULL

Download find_all_case_solution_of_ode_.mw

For Vectors, use a:=<2,1>, not <<2,1>>, which is a 2x1 Matrix. Then you can use "." for the dot product.

@Carl Love I think shape = symmetrix or shape = hermitian change the algorithm for float matrices to produce only real eigenvalues and eigenvectors but have no effect in the symbolic case.

@rzarouf You're welcome. You found a Maple bug, and I have submitted an SCR (bug report), so it is fixed in future versions.

3 4 5 6 7 8 9 Last Page 5 of 72