ecterrab

14160 Reputation

24 Badges

19 years, 161 days

MaplePrimes Activity


These are answers submitted by ecterrab

@Carl Love 

Carl, there is only one definition of integrating factor. "Suppose the equation is of the form ode and has an integrating factor mu: then mu*ode is a total derivative, i.e.: mu ode = d/dx F(y(x), x)". This is also the first sentence you read in the help page ?DEtools[intfactor]

You see you can rewrite mu ode as (mu*f)*(ode/f), and so as soon as you change the form of the equation (now it is ode/f), it also changes the integrating factor (now it is mu*f). That explains the results you see in Maple that you posted above in this thread.

The wikipedia and mathworld articles are weak in that all their examples have the highest derivative isolated (equivalent to multiplying ode by some f such that the coefficient of the highest derivative is equal to 1), which is what seems to have confused you. Note that in the equation posted by "nm" the coefficient of the highest derivative is t^2 - t. 

The fact that you can rewrite the integrating factor of a given ODE in infinitely many manners (I mean: without changing ODE) is easy to understand: suppose you have an expression G(y(x), x) that on the solutions y(x) of your ODE becomes a constant. Then d/dx G(y(x),x) = 0 on the solutions y(x) of ODE. But then if mu ODE is a total derivative, so is mu * G(x, y(x)) * ODE because its d/dx is equal to d/dx G * (mu ODE) + G * d/dx(mu ODE) = 0 + G * d/dx(mu ODE) = 0.

In other words: both mu and mu * G(x, y(x)) are integrating factors of ODE. You know, given y(x) (the solution of ODE) we can write infinitely many functions G(x, y(x)) = constant.

In summary: For any ODE you can write integrating factors in infinitely many diffeerent ways. The output of Maple's intfactor is correct. The integrating factor posted by 'nm' is not correct (it is missing a factor). 

Edgardo S. Cheb-Terrab
Physics, Maplesoft

Hi Andriy

Recall there is a difference between `*` and `.`, the latter is used for the scalar product of Bras and Kets and operators, while you can use `*` to represent the operation without computing it. Note also that phi is commutative and you use it as such when performing the product, but afterwards you replace it by a noncommutative object a Ket. That is not a good strategy: the order of the operands gets swapped in this case, normalizing the product taking advantage that 'phi is commutative'.

In summary: before multiplying, enter Setup(op = phi) or Setup(noncommutativeprefix = phi) to state that phi is not commutative, then perform H . phi as you do in your worksheet, then enter eval(H, [phi = Ket(psi, 1, 0, 1, 1, 0, 0), f = 1, `*` = `.`]) (so as you did, but including now the equation `*` = `.`), in order to introduce a Ket replacing phi and perform the scalar product. That resolves your 1st problem.

Regarding your second question: conjugate is present in your H, and split into Re and Im within simplify/conjugate, because in this example that leads to a partial simplification. Although at first sight I would say that partial simplification is just of no value in this case and complicates the expression. Until that is reviewed, I would suggest you replacing conjugate by something else like Conjugate, then perform the scalar product (that is where simplify/conjugate is called), and at the end remove Conjugate reintroducing conjugate. If I tweak simplify/conjugate before that I will let you know but this week I'm rather busy with other stuff.

Edgardo S. Cheb-Terrab
Physics, Maplesoft

Hi

What is what you are finding strange in intfactor's result? It is correct: try mutest(intfactor(ode), ode) and it returns zero. You can test this result manually too. Regarding the integrating factor you show, simplify(exp(int( (t^2+2*t-2)/(t^2-t),t))), it is not correct. Multiplying the ODE by it you do not get an exact equation. Try testing it with DEtools[mutest].

Independent of that, note that the integrating factor can depend on as many variables as the differential order plus 1 (all the derivatives up to the ODE order minus one, plus the unknown y(x) and the independent variable x) and as such you can write an integrating factor for an ODE in infinitely many different manners. Some of that is explained in ?DEtools[intfactor] and in ?DEtools[redode] that computes the ODE reducible by a given integrating factor; see also PDEtools[DeterminingPDE] to compute the PDE system satisfied by it.

Edgardo S. Cheb-Terrab
Physics, Maplesoft

@Andriy 

To avoid having to sort products many times, use the optional argument 'evaluateexpression'. The 'too many levels of recursion' was happening only when assigning to geometrical coordinates (z is part of the cartesian coordinates). When working with Physics:-Vectors, you cannot assign to these coordinates. It is fixed now: if you do not have Vectors loaded, you can assign any of the geometrical coordinates without problem. I just posted a new update in the Maplesoft webpage "Maple Physics: Research & Development" with only this fix.

From your post it remains: 1) why it is not sorting the second A in ABAC - I will revise this in the next iteration, and 2) why do whe have A . A and A * A.

The answer to 2) is debatable. As you say, A * B looks as good as A . B to represent their product, there is indeed no difference whatsoever in meaning between these two constructions, and we could have implemented everything with a single product operator. I still preferred to distinguish between non commutative generic products and scalar products of Kets or Kets and operators. This distinction is useful at least in two ways: using `*` you can express your quantum computations and then replacing by `.`, as in eval(expression, `*` = `.`), you can execute them, also having A . B, you can insert a Projector (see ?Physics:-Projector) in between and have the operation performed exactly as expected. The inconvenience of using `*` and `.` is that for a computer these two constructions are different, and this 'computer language difference' forces us to take care of it within the code or interactively. I will think if there is a way to diminish the inconvenience to a minimum.

Edgardo S. Cheb-Terrab
Physics, Maplesoft

I gave a look at your worksheet: a DE system for x,y,z, this system is exactly solvable, then you solve for x,y,z another system of three functions f, h, k equal to 0. No text. I is not clear to me what you are asking, nor why are you setting f, h and k equal to zero or what would be the relation between that and the system sys where none of f, h or k are equal to zero. I guess that is why you received no feedback. In order to receive more concrete feedback I suggest you repost explaining your problem more clearly. 

Independent of that, there is in Maple a nice package for plotting Poincare sections of dynamical systems - check ?DEtools[poincare] and to see this at work see ?Poincare,examples.

Edgardo S. Cheb-Terrab
Physics, Maplesoft

In Maple 17, also some releases before that, you can use the Physics:-Vectors package. For examples of its use in Vector Analysis problems like the ones you are mentioning, see for instance the 1st section of the help page ?Physics[examples].

Edgardo S. Cheb-Terrab
Physics, Maplesoft

Hi

This issue got fixed in the latest update of Physics (Sep/26), available for download as usual in the "Maple Physics: Research & Development" updates page. 

Edgardo S. Cheb-Terrab
Physics, Maplesoft

Hi

The old 'tensor' package is not maintained anymore and is being deprecated: tensors are now implemented in the Physics package, in a more natural way, including covariant derivatives, see ?Physics[D_]. So my suggestion would be for you to try that - if it doesn't work then could you please post the details, e.g. upload a worksheet with them, so that we can help you more concretely.

Edgardo S. Cheb-Terrab
Physics, Maplesoft 

To solve an equation in terms of  Bessel functions, you can indicate the solving method dsolve(ode, [Bessel]) (see ?dsolve[details]). In most cases that would be equivalent to using Frobenius series then rewriting the series in terms of Bessel functions. You can do that with any of the three equations you show. The one solved in terms of sin and cos you can express that solution in terms of Bessels using convert(solution, Bessel, include = all) (see ?convert[to_special_function]).

Regarding Frobenius series solution, you can formulate such a solution around a regular singularity of the quation (check ?DEtools[singularities]). Given the location of such a singularity, you can use ?DEtools[formal_sol] to see the first terms around that location, or also dsolve with its series option, or when a closed form summation solution can be computed you can directly use ?dsolve[formal_solution] as Preben mentioned, optionally indicating the location of the singularity.

Edgardo S. Cheb-Terrab
Physics, Maplesoft

Hi Andriy, Carl

I am unable to reproduce the problem you mention. Here I am loading Maple with no initialization whatsoever

restart; anames()

`debugger/no_output`, TestTools, renumber, TRY, buildTRY

(1)

I am running Maple 17.02

version()

 User Interface: 872941

         Kernel: 872941
        Library: 872941

 

872941

(2)

The Physics package loaded in 17.02 has datestamp from September 5

Physics:-Version()

"/Library/Frameworks/Maple.framework/Versions/17/lib/update.mla", `2013, September 5, 15:45 hours`

(3)

And I cannot reproduce your problem:

with(Physics)

a := ``

``

(4)

``(3)

``(3)

(5)

Now loading the latest update of the package distributed inthe "Maple Physics: Research & Development" updates page.

restart; anames()

`debugger/no_output`, TestTools, renumber, TRY, buildTRY

(6)

For that, add the path to the udpated library

libname := "/Users/ecterrab/Maple/lib", libname;

"/Users/ecterrab/Maple/lib", "/Library/Frameworks/Maple.framework/Versions/17/lib", "/Users/ecterrab/maple/toolbox/emacs/lib", "/Library/Frameworks/Maple.framework/Versions/17/toolbox/NAG/lib", "."

(7)

The updated version of Physics has datestamp from yesterday

Physics:-Version();

"/Users/ecterrab/Maple/lib/Physics.mla", `2013, September 20, 12:37 hours`

(8)

with(Physics):

a := ``;

``

(9)

``(3);

``(3)

(10)

So everything works correctly, as expected. I'd need more info from you in order to understand what the problem you have is.


Download unabletoreproduce.mw

Edgardo S. Cheb-Terrab
Physics, Maplesoft

@Andriy 

Note that cc is not a Maple command but just written to address your particular problem: to get all the coefficients in products. You are now passing something that is not a product, the term '+ap[1]'. Just tune your cc accordingly. Here is a cc that meets this new requirement of your problem:

cc := proc (EE)
local ee, P, U;

ee := UnnestProductsInExpression(EE, ':-productsinoutput = P');
P := {seq(P[j] = U[j], j = 1 .. nops(P))};
ee := subs(P, ee);
coeffs(ee, map(rhs, P) union indets(ee, Library:-PhysicsType:-QuantumOperator))
end:

Note also that Library:-PhysicsType is a package, it has many other types that may be of use for you, for instance ExtendedQuantumOperator.

Regarding your second question: look at the operands of phi (as in: op(phi)) and you will understand why is that you are getting only a + b, and how to solve it; one way of doing that is select(type, phi, commutative).

Edgardo S. Cheb-Terrab
Physics, Maplesoft

Hi Andryi
Regarding the problem you posted and getting all the coefficients in one go, a new update of Physics is being posted today, with it:

Physics:-Version();

(1)

for i to 4 do ap[i] := Creation(psi, i, notation = explicit); am[i] := Annihilation(psi, i, notation = explicit) end do:

ApAm1 := [seq(ap[j], j = 1 .. 4), seq(am[j], j = 1 .. 4)]:

z := Typesetting:-delayDotProduct(Typesetting:-delayDotProduct(Typesetting:-delayDotProduct(Physics:-`*`(D1, ap[1]), am[2]), ap[4]), am[3])+Typesetting:-delayDotProduct(Typesetting:-delayDotProduct(Typesetting:-delayDotProduct(Physics:-`*`(A, am[2]), ap[2]), am[1]), ap[1])+Typesetting:-delayDotProduct(Typesetting:-delayDotProduct(Typesetting:-delayDotProduct(Physics:-`*`(B, ap[1]), am[1]), ap[2]), am[2])+Typesetting:-delayDotProduct(Typesetting:-delayDotProduct(Typesetting:-delayDotProduct(Physics:-`*`(B, ap[3]), am[3]), ap[4]), am[4])+Typesetting:-delayDotProduct(ap[1], Typesetting:-delayDotProduct(Typesetting:-delayDotProduct(Physics:-`*`(D2, am[2]), ap[4]), am[3]))+Typesetting:-delayDotProduct(Physics:-`*`(3, ap[1]), am[1])+Typesetting:-delayDotProduct(Physics:-`*`(R, ap[1]), am[1])+Typesetting:-delayDotProduct(ap[2], am[2])+Typesetting:-delayDotProduct(Physics:-`*`(Physics:-`*`(K[1], 1/(C+Physics:-`*`(5, D))), ap[3]), am[4])+Typesetting:-delayDotProduct(Physics:-`*`(Physics:-`*`(L[5], 1/(C+Physics:-`*`(5, D))), ap[3]), am[4]):

(-A+R+3)*Physics:-`.`(`a+`[psi[1]], `a-`[psi[1]])-B*Physics:-`.`(`a+`[psi[1]], Physics:-`.`(`a+`[psi[2]], Physics:-`.`(`a-`[psi[1]], `a-`[psi[2]])))-(D1+D2)*Physics:-`.`(`a+`[psi[1]], Physics:-`.`(`a+`[psi[4]], Physics:-`.`(`a-`[psi[2]], `a-`[psi[3]])))+A*Physics:-`.`(`a+`[psi[1]], Physics:-`.`(`a-`[psi[1]], Physics:-`.`(`a+`[psi[2]], `a-`[psi[2]])))-(A-1)*Physics:-`.`(`a+`[psi[2]], `a-`[psi[2]])+(L[5]+K[1])*Physics:-`.`(`a+`[psi[3]], `a-`[psi[4]])/(C+5*D)-B*Physics:-`.`(`a+`[psi[3]], Physics:-`.`(`a+`[psi[4]], Physics:-`.`(`a-`[psi[3]], `a-`[psi[4]])))+A

(2)

You asked about getting all the coefficients of products of noncommutative operands in one single command call, without having to tell what these products are, etc. Here is a cc procedure that accomplishes that:

cc := proc (EE) local ee, P, U; ee := UnnestProductsInExpression(EE, ':-productsinoutput = P'); P := [seq(P[j] = U[j], j = 1 .. nops(P))]; ee := subs(P, ee); coeffs(ee, map(rhs, P)) end proc:

cc(z)

A, -A+R+3, -A+1, (L[5]+K[1])/(C+5*D), -B, -D1-D2, A, -B

(3)


Download AllCoefficients.mw

Edgardo S. Cheb-Terrab
Physics, Maplesoft

@Carl Love 

Note that in the same page ?assuming[details] you are quoting, after these three items you quote, you also read:

Notes: The assuming command does not place assumptions on integration or summation dummy variables in definite integrals and sums, nor in limit or product dummy variables, because all these variables already have their domain restricted by the integration, summation or product range or by the method used to compute a limit. To obtain the simplification of the expression being summed, integrated or subject to a product taking into account the restriction on the values of the dummy variable implicit in the integration/summation/product range, use the simplify command -- see the Examples section of this help page.

The assuming command does not scan Maple programs regarding the potential presence of assumed variables. To compute taking into account assumptions on variables found in the bodies of programs, use assume -- see the Examples section.

These two paragraphs complete the description. Taking them into account, assuming does work according to its design.

Regarding the question by Sergio, the second paragraph quoted above is the answer. There are various reasons for this design, mainly that if you scan procedures, they typically call other procedures that call other ones etc. and you cannot go scanning all the tree, and "I only scan this many levels" would be of no use in various cases one could imagine, and you have the functionality anyway if you use assume alone or in combination with assuming (with or without its option additionally). There is one example in the Examples section illustrating this aspect of the design explicitly.

Edgardo S. Cheb-Terrab
Physics, Maplesoft

A new Physics:-Library:-Add command is included in today's update of Physics, implementing the ideas discussed in this post, so working around the model of evaluation of arguments of sum that sometimes generates unexpected results as discussed for instance in the Mapleprimes posts of August  annihilation/creation operators , Problem-With-A-Function-Within-Sum, and Sum-In-For-Loop. The next upated by the end of the week will include a `print/Add` to have the expected capital Sigma display an related copy & paste.

Edgardo S. Cheb-Terrab
Physics, Maplesoft

@Andriy 


You posted 3 replies. Let me go through them in sequence in a single reply here, and in what follows I am using the next update of Physics, to be posted tomorrow:

restart; with(Physics); Physics:-Version()

"/Users/ecterrab/Maple/lib/Physics.mla", `2013, September 16, 16:33 hours`

(1)

____________________________________

 

In your 'Great but not perfect', the problem is not in SortProducts but in AntiCommutator. You know, for Commutators you can always expand them when one of the operands is a product, say as in,

Setup(quantumoperators = {A, B, C, F, G})

Commutator(A*B, C)

Physics:-`*`(A, Physics:-Commutator(B, C))+Physics:-`*`(Physics:-Commutator(A, C), B)

(2)

But this formula does not exist for anticommutators, so anticommutators of products are not expanded

AntiCommutator(A*B, C)

Physics:-AntiCommutator(Physics:-`*`(A, B), C)

(3)

On the other hand - say in this example - if the anticommutator between B and C and also between A and C are known, then one can always compute this anticommutator as an expansion.


Example:

Setup(%AntiCommutator(A, C) = F, %AntiCommutator(B, C) = G)

[algebrarules = {%AntiCommutator(A, C) = F, %AntiCommutator(B, C) = G}]

(4)

The following anticommutator returns uncomputed with the version of Physics you have, but from the knowlege of the anticommutatores (rules) above it can always be computed:

AntiCommutator(A*B, C)

2*Physics:-`*`(C, A, B)+Physics:-`*`(A, G)-Physics:-`*`(F, B)

(5)

In steps:

%AntiCommutator(A*B, C); % = expand(%)

%AntiCommutator(Physics:-`*`(A, B), C) = Physics:-`*`(A, B, C)+Physics:-`*`(C, A, B)

(6)

Take the first term on the right-hand side and move C to the left in two steps: first anticommuting with C with B then with A.

Step 1

%AntiCommutator(B, C); % = expand(%)

%AntiCommutator(B, C) = Physics:-`*`(B, C)+Physics:-`*`(C, B)

(7)

isolate(%, B*C)

Physics:-`*`(B, C) = %AntiCommutator(B, C)-Physics:-`*`(C, B)

(8)

value(%)

Physics:-`*`(B, C) = G-Physics:-`*`(C, B)

(9)

A*lhs(%) = map2(`*`, A, rhs(%))

Physics:-`*`(A, B, C) = Physics:-`*`(A, G)-Physics:-`*`(A, C, B)

(10)

subs(%, %AntiCommutator(Physics[`*`](A, B), C) = Physics[`*`](A, B, C)+Physics[`*`](C, A, B))

%AntiCommutator(Physics:-`*`(A, B), C) = Physics:-`*`(A, G)-Physics:-`*`(A, C, B)+Physics:-`*`(C, A, B)

(11)

Step 2

%AntiCommutator(A, C); % = expand(%)

%AntiCommutator(A, C) = Physics:-`*`(A, C)+Physics:-`*`(C, A)

(12)

isolate(%, A*C)

Physics:-`*`(A, C) = %AntiCommutator(A, C)-Physics:-`*`(C, A)

(13)

value(%)

Physics:-`*`(A, C) = F-Physics:-`*`(C, A)

(14)

lhs(%)*B = map(`*`, rhs(%), B)

Physics:-`*`(A, C, B) = Physics:-`*`(F, B)-Physics:-`*`(C, A, B)

(15)

subs(%, %AntiCommutator(Physics[`*`](A, B), C) = Physics[`*`](A, G)-Physics[`*`](A, C, B)+Physics[`*`](C, A, B))

%AntiCommutator(Physics:-`*`(A, B), C) = 2*Physics:-`*`(C, A, B)+Physics:-`*`(A, G)-Physics:-`*`(F, B)

(16)

This mechanism to compute the anticommutator of products when the anticommutator for each of the operands of the product is known was not implemented, and so SortProducts was returning with an AntiCommutator uncomputed. After implementing this, for the example you posted in "Great but not perfect" we now have

Setup(anticommutativeprefix = psi)

[anticommutativeprefix = {_lambda, psi}]

(17)

for i to 4 do ap[i] := Creation(psi, i, notation = explicit); am[i] := Annihilation(psi, i, notation = explicit) end do

ee := ap[2].am[2].am[1].ap[1]

Physics:-`.`(`a+`[psi[2]], Physics:-`.`(`a-`[psi[2]], Physics:-`.`(`a-`[psi[1]], `a+`[psi[1]])))

(18)

ApAm1 := [seq(ap[j], j = 1 .. 4), seq(am[j], j = 1 .. 4)]

[`a+`[psi[1]], `a+`[psi[2]], `a+`[psi[3]], `a+`[psi[4]], `a-`[psi[1]], `a-`[psi[2]], `a-`[psi[3]], `a-`[psi[4]]]

(19)

No uncomputed anticommutators:

Library:-SortProducts(ee, ApAm1, useanticommutator)

Physics:-`.`(`a+`[psi[1]], `a+`[psi[2]], `a-`[psi[1]], `a-`[psi[2]])+Physics:-`.`(`a+`[psi[2]], `a-`[psi[2]])

(20)

____________________________________

 

In your next reply "Another issue appeared", you point out that simplify/size does not collect the terms the way you want and is expected.

Setup(clear, op = C, quiet)

[quantumoperators = {A, B, F, G}]

(21)

z2 := 3*ap[1].am[1]+R*ap[1].am[1]+K*(ap[2].am[2])/(C+D)

3*Physics:-`.`(`a+`[psi[1]], `a-`[psi[1]])+R*Physics:-`.`(`a+`[psi[1]], `a-`[psi[1]])+K*Physics:-`.`(`a+`[psi[2]], `a-`[psi[2]])/(C+D)

(22)

simplify(z2, size)

3*Physics:-`.`(`a+`[psi[1]], `a-`[psi[1]])+R*Physics:-`.`(`a+`[psi[1]], `a-`[psi[1]])+K*Physics:-`.`(`a+`[psi[2]], `a-`[psi[2]])/(C+D)

(23)

Note three things.

1) This result is not collected as we expect, but it is not wrong as you say.

2) The problem is not related to Physics.

 

For example, remove all Physics objects

subs(am[1] = am1, ap[1] = ap1, am[2] = am2, ap[2] = ap2, `.` = :-`*`, %)

3*(ap1*am1)+R*(ap1*am1)+K*(ap2*am2)/(C+D)

(24)

Try again and you see the same situation, the expression is not reduced in size

simplify(%, size)

3*ap1*am1+R*ap1*am1+K*ap2*am2/(C+D)

(25)

3) You can achieve what you want with the expression z2 directly using collect, as in

collect(z2, `.`)

(R+3)*Physics:-`.`(`a+`[psi[1]], `a-`[psi[1]])+K*Physics:-`.`(`a+`[psi[2]], `a-`[psi[2]])/(C+D)

(26)

____________________________________

 

In your third and last reply, "The coefficients of operator expression", you show that Coefficients is not working as expected. Indeed it takes coefficients in noncommutative products expressed using Physics:-`*`, not Physics:-`.`. This is explained in the help page.


On the other hand I agree with you, Coefficients should work with both Physics:-`*` and Physics:-`.` in equal footing. This is implemented now (to appear tomorrow), so that we get

Coefficients(z2, ap[1])

K*Physics:-`.`(`a+`[psi[2]], `a-`[psi[2]])/(C+D), (R+3)*`a-`[psi[1]]

(27)

Coefficients(z2, ap[1], 0)

K*Physics:-`.`(`a+`[psi[2]], `a-`[psi[2]])/(C+D)

(28)

It works also with a product as second argument, as in:

Coefficients(z2, [ap[1].am[1]], 1)

R+3

(29)

Coefficients(z2, [ap[1].am[1]], 0)

K*Physics:-`.`(`a+`[psi[2]], `a-`[psi[2]])/(C+D)

(30)

Coefficients(z2, [ap[1].am[1]])

K*Physics:-`.`(`a+`[psi[2]], `a-`[psi[2]])/(C+D), R+3

(31)

The key here was adapting Library:-Degree to work withPhysics:-`*` and Physics:-`.` in equal footing, and with that we got Coefficients working as expected:

Library:-Degree(z2, ap[2], minmax = both)

0, 1

(32)

Library:-Degree(z2, [ap[2], am[2]], minmax = both)

0, 2

(33)

Library:-Degree(z2, [ap[2].am[2]], minmax = both)

0, 1

(34)

In summary for your three replies: 1) AntiCommutator got improved (that resolved your issue with SortProducts), 2) use collect (simplify/size works very well but in this example requires a tweak and I am short of time to fix that); 3) Library:-Degree, and through it also Coefficients now work with Physics:-`.` the same way they do with Physics:-`*`.

 

All these changes and some others as a new Library:-Add and Library:-TensorCoefficients are already included in the update being prepared for tomorrow. Regarding Library:-Collect, yes it is part of the plan, but other things are a bit higher in the priority list (e.g. Hausorff's formula, relevant in basic QM). Anyway, at some point soon enough Collect will be in place together with a new Factor handling noncommutative products.

NULL



Download MaplePrimesSortProducts.mw

 

Edgardo S. Cheb-Terrab
Physics, Maplesoft

First 51 52 53 54 55 56 57 Page 53 of 58