# Kapitel 2_i: Kinematik eines ebenen Viergelenkgetriebes
# © 2015 Friedrich U. Mathiak
# mathiak@mechanik-info.de
# 
> restart:with(LinearAlgebra):with(plottools):with(plots):
# Beispiel 2-16: (Ebenes Vielgelenkgetriebe)
# Wir notieren die vier Vektoren der Getriebeglieder. Die Vektororientierungen werden mathematisch positiv entgegengesetzt dem Uhrzeigersinn eingeführt. Dasselbe gilt für die Drehwinkel. 
> R1:=l1*Vector([cos(phi1),sin(phi1)]);
> R2:=l2*Vector([cos(phi2),sin(phi2)]);
> R3:=l3*Vector([cos(phi3),sin(phi3)]);
> R4:=l4*Vector([cos(phi4),sin(phi4)]);
# Damit lautet die Schleifengleichung 
> R:=R1+R2+R3+R4;
> phi1:=0;
> R;
# Es ist folgendes Gleichungssystem zu lösen: 
> gl1:=R[1];gl2:=R[2];
> loe:=solve({gl1,gl2},{phi3,phi4},AllSolutions=true);
# In der Lösung treten die unbestimmen  Ausdrücke `~`[_Z1]*und*`~`[_Z2];  auf, die wir mit dem Befehl indets ansehen und denen wir Werte  zuweisen können. 
> indets(loe);
# Mit dem Befehl about (_Z~) erhalten wir Informationen über die Werte von _Z~.
> about(_Z1);about(_Z2);
# Wir weisen den unbestimmten Ausdrücken den Wert 0 zu.
> loe:=subs(op(1,indets(loe))=0,op(2,indets(loe))=1,loe): assign(%);
# Mit den Werten des Beispiels sind die Längen des Viergelenkgetriebes:
> l1:=150;l2:=60;l3:=160;l4:=100;
# Wir prüfen die Umlauffähigkeit des Getriebes
> ls:=sort([l1,l2,l3,l4],`>`); S1:=ls[1]+ls[4]; S2:=ls[2]+ls[3];
> if S1 > S2 then
>   print(`System ist nicht umlauffähig`);   else
>   print(`System ist umlauffähig`);
> end if;
# Maple liefert uns für phi3 und phi4  jeweils zwei Lösungen:
> phi3:=allvalues(phi3):
> phi4:=allvalues(phi4):
# Wir stellen die Lösungen grafisch dar. 
> eps:=10.^(-3):
> plot([phi3,phi4],phi2=0..2*Pi-eps,title= "Drehwinkel eines Viergelenkgetriebes\n",axes=boxed,gridlines=true,legend = ["phi3_1","phi3_2","phi4_1","phi4_2"], titlefont = ["ARIAL", 14], labels = ["Antriebswinkel phi2\n", "Winkel [rad]"], labeldirections = ["horizontal","vertical"],labelfont = ["HELVETICA", 10],color=[black,blue,red,green], axesfont = ["HELVETICA", "ROMAN", 8],legendstyle = [font = ["HELVETICA", 9], location = bottom]);
# Die Lösungskurven besitzen an der Stelle phi2 = À eine Unstetigkeit. Dort ändern sich beim Übergang die Vorzeichen der Winkel phi3 und phi4.
> phi3_1:=phi3[1]: phi3_2:=phi3[2]:
> phi4_1:=phi4[1]: phi4_2:=phi4[2]:
# Wir bilden aus den obigen Lösungen  mittels piecewise stetige und stetig differenzierbare Funktionen für die Winkel phi3 und phi4, die zu der im Buch TM III skizzierten Einbaulage führen.
> phi3:=piecewise(phi2 < Pi,phi3_2,phi2 > Pi,phi3_1):
> phi4:=piecewise(phi2 < Pi,phi4_2,phi2 > Pi,phi4_1):
> plot([phi3,phi4],phi2=0..2*Pi,scaling=constrained,title= "Drehwinkel eines Viergelenkgetriebes\n",axes=boxed,gridlines=true,legend = ["phi3","phi4"], titlefont = ["ARIAL", 14], labels = ["Antriebswinkel phi2\n","Winkel [rad]"], labeldirections = ["horizontal","vertical"],labelfont = ["HELVETICA", 10],color=[black,blue], axesfont = ["HELVETICA", "ROMAN", 8],legendstyle = [font = ["HELVETICA", 9], location=bottom]);
# Wir prüfen, ob mit den so zusammengesetzten Winkeln die Schleifengleichungen erfüllt sind:
> gl1w:=evalf(gl1):gl2w:=evalf(gl2):
> plot([gl1w,gl2w],phi2=0..2*Pi-eps);
# Animation der Viergelenkkette.
> ra:=   [R1[1],R1[2]]:
> rb:=ra+[R2[1],R2[2]]:
> rc:=rb+[R3[1],R3[2]]:
> linie:={}:
> linie:=linie union {[ra,rb]} union {[rb,rc]} union {[rc,[0,0]]}:
> plots[animate](plot,[linie,color = "NavyBlue",thickness=2,scaling=constrained],phi2=0..2*Pi,axes=boxed,gridlines=true,frames=180);
# Bei schwingenden Getriebegliedern (hier Glied 4)  wird die Bewegung durch die Umkehrlagen (Totlagen) begrenzt.  Zu deren Berechnung bilden wir die 1. Ableitung des Abtriebswinkels phi4 und bestimmen die beiden Nullstellen:
> phi4str:=diff(phi4,phi2):
> plot(phi4str,phi2=0..2*Pi-eps,axes=boxed,gridlines=true);
> loe1:=fsolve(phi4str,phi2=3/4*Pi..Pi);evalf(convert(loe1,degrees));
> loe2:=fsolve(phi4str,phi2=3/2*Pi..2*Pi);evalf(convert(loe2,degrees));
> phi4a:=evalf(subs(phi2=loe1,phi4)); evalf(convert(%,degrees));
> phi4i:=evalf(subs(phi2=loe2,phi4)); evalf(convert(%,degrees));
# Ermittlung der Koppelkurve eines Punktes K auf der Koppel (n: Anzahl der Unterteilungen des Koppelgliedes)
> n:=3: dlam:=1/n: lambda:=-dlam:
> for i to n+1 do
>   lambda:=lambda+dlam:
>   rK:=rb+lambda*(rc-rb):
>   pl[i]:=plot([rK[1],rK[2],phi2=0..2*Pi-eps],color="NavyBlue",filled = [color = "Blue", transparency = .5]):
> end do:
> plots[display](seq(pl[i],i=1..n+1),title ="\nKoppelkurven",titlefont = ["ARIAL", 14],axes=boxed,gridlines=true,scaling=constrained,labels = ["x1", "x2"], labeldirections = ["horizontal", "vertical"], labelfont = ["HELVETICA", 10]);
# Wir bilden die 1. Ableitung der Drehwinkel Æ3 und Æ4  nach Æ2 und erhalten
> phi3_p:=simplify(diff(phi3,phi2)):
> phi4_p:=simplify(diff(phi4,phi2)):
> plot([phi3_p,phi4_p],phi2=0..2*Pi,gridlines=true,axes=boxed,legend = ["phi3p","phi4p"], titlefont = ["ARIAL", 14],color=[black,blue], axesfont = ["HELVETICA", "ROMAN", 8],legendstyle = [font = ["HELVETICA", 9], location=bottom]);
# Wir bilden die 2. Ableitung der Drehwinkel Æ3 und Æ4  nach Æ2 und erhalten
> phi3_pp:=simplify(diff(phi3_p,phi2)):
> phi4_pp:=simplify(diff(phi4_p,phi2)):
> plot([phi3_pp,phi4_pp],phi2=0..2*Pi,gridlines=true,axes=boxed,legend = ["phi3pp","phi4pp"], titlefont = ["ARIAL", 14],color=[black,blue], axesfont = ["HELVETICA", "ROMAN", 8],legendstyle = [font = ["HELVETICA", 9], location=bottom]);
# Wir bilden die 3. Ableitung der Drehwinkel nach der Zeit t , um die Änderung der Beschleunigung (Ruck) aufzuzeigen, die für die Belastung des Getriebes von Bedeutung ist.
> phi3_ppp:=simplify(diff(phi3_pp,phi2)):
> phi4_ppp:=simplify(diff(phi4_pp,phi2)):
> plot([phi3_ppp,phi4_ppp],phi2=0..2*Pi-eps,gridlines=true,axes=boxed,legend = ["phi3ppp","phi4ppp"], titlefont = ["ARIAL", 14], color=[black,blue], axesfont = ["HELVETICA","ROMAN",8],legendstyle=[font=["HELVETICA", 9], location=bottom]);
# Wie wir der obigen Grafik entnehmen können, ergeben sich die größten Beschleunigungsänderungen in der Umgebung von phi2 = À.
# Wir berechnen die Geschwindigkeiten der Punkte  B und C und wählen dazu diff(varphi2, t) = omega0 and omega0 = 1.; 
> vB:=omega0*Vector([-R2[2],R2[1]]):VBN:=VectorNorm(vB,Euclidean):
> vC:=phi4_p*Vector([R4[2],-R4[1]]):VCN:=VectorNorm(vC,Euclidean):
> omega3:=Transpose(vC+vB).R3/l3^2:
> omega0:=1.;
> plot([vB[1],vB[2]],phi2=0..2*Pi,gridlines=true,axes=boxed,title ="\nGeschwindigkeitkomponenten des Punktes B\n",titlefont = ["ARIAL", 14],legend = ["vB_1","vB_2"], titlefont = ["ARIAL", 14], color=[black,blue], axesfont = ["HELVETICA","ROMAN",8],legendstyle=[font=["HELVETICA", 9], location=bottom]);
> plot([vC[1],vC[2]],phi2=0..2*Pi,gridlines=true,axes=boxed,title ="\nGeschwindigkeitkomponenten des Punktes C\n",titlefont = ["ARIAL", 14],legend = ["vC_1","vC_2"], titlefont = ["ARIAL", 14], color=[black,blue], axesfont = ["HELVETICA","ROMAN",8],legendstyle=[font=["HELVETICA", 9], location=bottom]);
# Wir geben noch die Beträge der Geschwindigkeite für die Punkte B und C aus
> plot([VCN,VBN],phi2=0..2*Pi,gridlines=true,axes=boxed,title ="\nGeschwindigkeitkomponenten des Punktes C\n",titlefont = ["ARIAL", 14],legend = ["vC","vB"], titlefont = ["ARIAL", 14],color=[black,blue], axesfont = ["HELVETICA","ROMAN",8],legendstyle=[font=["HELVETICA", 9], location=bottom]);
# Die betragsmäßig größte Geschwindigkeit des Punktes C ergibt sich an der Stelle
> phi2max:=fsolve(diff(VCN,phi2),phi2=3..4);
> maxVCN:=evalf(subs(phi2=phi2max,VCN));
# Wir berechnen die Winkelgeschwindigkeit der Koppel
> plot(omega3,phi2=0...2*Pi,gridlines=true,axes=boxed,title ="\nWinkelgeschwindigkeit der Koppel\n",titlefont = ["ARIAL", 14]);
> 
;
# 
