gnu-octave / symbolic

A Symbolic Package for Octave using SymPy
https://octave.sourceforge.io/symbolic/
GNU General Public License v3.0
152 stars 36 forks source link

subs() removing equations #835

Closed DanielMartensson closed 7 years ago

DanielMartensson commented 7 years ago

Hi!

I want to report a bug. When I doing subs and want to replace som symbolic variables with other symbolic variables, the subs command remove equations.

Look at the diffence between the variable LA1 and LA2 and you will se that LA1 contains full equation but LA2 have lost a lot of variables. I just want to replace d(theta1(t))/dt with Q1. Just to make it easier to use the jacobian later below.

I think that OctSymPy cannot handle LA1 because when I printing out LA1 so stops everything if I press (f) for forward. It's says that it needs to collect data.

"Waiting for data... [interupt to abort]"

% Ladda paket
pkg load symbolic

% Variabler
syms m1 m2 m3 L1 L2 L3 theta1(t) theta2(t) theta3(t) g

% Positioner
x1 = L1*sin(theta1(t));
y1 = L1*cos(theta1(t));
x2 = L1*sin(theta1(t)) + L2*sin(theta2(t));
y2 = L1*cos(theta1(t)) + L2*cos(theta2(t));

% Hastigheter
dx1 = simplify(diff(x1));
dy1 = simplify(diff(y1));
dx2 = simplify(diff(x2));
dy2 = simplify(diff(y2));

 % Resultathastigheter:
 dyx1 = simplify(dx1^2 + dy1^2);
 dyx2 = simplify(dx2^2 + dy2^2);

 % Tröghetsmoment
 J1 = simplify(m1*(x1^2 + y1^2) + m2*(x2^2 + y2^2));
 J2 = simplify(m2*((L2*sin(theta2(t)))^2 + (L2*cos(theta2(t)))^2));
 J3 = simplify(m1*y1^2 + m2*y2^2 + 1/2*m3*(L3/2)^2);

 % Potentionella energier
 PE1 = simplify(m1*g*x1);
 PE2 = simplify(m2*g*(x1+x2));

 % Langrange ekvationen
 LA1 = simplify((1/2*m1*dyx1 + J1*diff(theta1(t))^2 + 1/2*m2*dyx2 + J2*diff(theta2(t))^2 + J3*diff(theta3(t))^2) - (PE1 + PE2));
 syms Q1 Q2 Q3 Q4 Q5 Q6
 LA2 = subs(LA1, {theta1(t) theta2(t) theta3(t) diff(theta1(t)) diff(theta2(t)) diff(theta3(t))}, {Q1 Q2 Q3 Q4 Q5 Q6})
 % Lagrangian till Newtonian
 diff(jacobian( LA2, [Q4 Q5 Q6])) - jacobian(LA2, [Q1 Q2 Q3])

OS: Windows 10 GNU Octave 4.2.1 Symbolic version 2.6.0

cbm755 commented 7 years ago
  1. Is this really the simplest possible example? https://en.wikipedia.org/wiki/Minimal_Working_Example

  2. I guess this is related to diff(Q1, t) being formally zero. Does it help if you first sub in for the derivatives, then sub in the for non-derivatives?

I think that OctSymPy cannot handle LA1 because when I printing out LA1 so stops everything if I press (f) for forward. It's says that it needs to collect data.

This seems like an unrelated issue (?) I haven't been able to reproduce it yet...

cbm755 commented 7 years ago

I guess this is a MWE of what you're talking about. Next time, please do this step yourself.

>> syms theta(t) Q1 Q2
>> A = 2*theta(t) + 3*diff(theta)
warning: worked around octave bug #42735
warning: called from
    plus at line 42 column 5
A(t) = (symfun)

             d       
  2⋅θ(t) + 3⋅──(θ(t))
             dt      

>> subs(A, theta, Q1)    # ok
ans = (sym) 2⋅Q₁
>> subs(A, diff(theta), Q2)     # ok
ans = (sym) 3⋅Q₂ + 2⋅θ(t)
>> subs(A, {theta diff(theta)}, {Q1 Q2})    # oh no
ans = (sym) 2⋅Q₁
DanielMartensson commented 7 years ago

It does not helping if i first subs in the derivatives and then subs in the non-derivatives. The problem is when I sub in Q2 and Q3, not Q1, Q3, Q4, Q5, Q6.

cbm755 commented 7 years ago

MWE please.

cbm755 commented 7 years ago

To elaborate: in my MWE, you can workaround the problem by doing:

>> temp = subs(A, diff(theta), Q2);
temp = (sym) 3⋅Q₂ + 2⋅θ(t)
>> subs(temp, theta, Q1)
ans = (sym) 2⋅Q₁ + 3⋅Q₂

Why doesn't that trick work for your problem?

cbm755 commented 7 years ago

I'm beginning to think this isn't a bug.

It seems to me that subs(A, {theta diff(theta)}, {Q1 Q2}) and subs(A, {diff(theta) theta}, {Q2 Q1}) really should give different answers: the order in which you do those substitution matters!

(for awhile I thought this shouldn't be the case because we use simultaneous=True but I've reread the SymPy docs:

If the keyword simultaneous is True, the subexpressions will not be evaluated until all the substitutions have been made.

(notably, it does not promise the result is input-order independent). [edit, accidently had this in the quote]

DanielMartensson commented 7 years ago

Did the equation LA2 become very small for you, compered to LA1?

DanielMartensson commented 7 years ago

@cbm755

Hi! Here you got a print screen of my 1x2 matrix N, which are very large. I did first sub every variable in the system to NQ1. Didn't work well. Then I sub first the derivatives, then sub the position variables. Works OK. But this is not how it should be done. This is a deep bug.

https://image.ibb.co/euSXSF/Untitled2.png

I cannot trust the subs() command anymore!

Edit: If I do one variable for each subs, then I can trust it.

cbm755 commented 7 years ago

But this is not how it should be done. This is a deep bug.

As I tried to explain above, I do not think its a bug. Result depends on the order the substitutions are made...

Then I sub first the derivatives, then sub the position variables. Works OK.

Good, that is what I expected. I believe that is the correct thing to do here.

Perhaps we should try to come up with a nice example for "help subs" as a warning for others.

DanielMartensson commented 7 years ago

Perhaps we should try to come up with a nice example for "help subs" as a warning for others.

Everytime a user uses subs, the user got a notation from the terminal prompt: "Notice that it's not good to subs all variable at the sime time. Write "help subs" for more information to close this message".

That would be great?

cbm755 commented 7 years ago

no that would be terrible. More like: "Caution, substitutions happen in the order you specify them. If a latter substitution depends on a previous one, it might be better to use two steps in stead".

cbm755 commented 7 years ago

I mean in "help subs" we add: "Caution, substitutions happen in the order you specify them. If a latter substitution depends on a previous one in an ambiguous way, it might be better to call subs twice. A cautionary example:"

>> syms v(x)
>> A = x + v(x)
A = (sym) x + v(x)
>> subs(A, {x, v(x)}, {6, 42})
ans = (sym) v(6) + 6
>> subs(A, {v(x), x}, {42, 6})
ans = (sym) 48

"Here it would be safer to specify exactly what you mean:"

>> temp = subs(A, v(x), 42);
>> subs(temp, x, 6)
ans = (sym) 48
cbm755 commented 7 years ago

@DanielMartensson would you proofread #839 for me? Give that PR at "+1" if you feel those comments would help other people avoid this.

Also: "so stops everything if I press (f) for forward. It's says that it needs to collect data." is this something you can reproduce?

DanielMartensson commented 7 years ago

Where do I "+1" on comments?

Also: "so stops everything if I press (f) for forward. It's says that it needs to collect data." is this something you can reproduce?

Yes! Try this:

A = [1 3 4; 2 4 5; 0 3 2] syms t expm(A*t)

You will get so large equation that the computer cannot understanding it. Right know it takes a lot of time to repoduce this "everything if I press (f) for forward. It's says that it needs to collect data" due to the expm() command.

cbm755 commented 7 years ago

Just go to #839, read over the changes, then leave a comment. "+1" is a common comment on github for approving something to be merged.