Zoals ik al mailde:simulink werkt niet met symbolische tekens (tenminste, niet bij mijn weten).
Als je parametrisch wilt werken: f*ck simulink en gebruik gewoon m-script.
Ik heb vroeger, toen ik nog op de TU zat, ooit nog eens zoiets met een massa veertje gedaan op een laptop regel-examen. Hieronder de code van dat examen....
Code: Select all
close all
clear all
clc
% Het systeem is simpel:
% twee massa's (parallel) verbonden met een veer en een demper.
% beide massa's hebben wrijving met de grond
% de eerste massa wordt aangedreven middels een kracht.
% Het probleem van dit voorbeeld zou een positioneringsprobleem kunnen
% zijn, waarbij:
% de tweede massa moet gepositioneerd worden
% de acceleratie van de eerste/tweede massa binnen bepaalde grenzen moet
% blijven, maar bijv. ook de veer+demper kracht
% (en deze moeten dus als output moeten worden opgenomen)
%
% Tekenafspraak: indrukken van de veer/damper levert een POSITIEVE kracht
% definieer eerst alle te gebruiken symbolen
syms ddz1 dz1 z1 ddz2 dz2 z2 % alle coordinaten (+afgeleiden)
syms u1 % alle inputs
syms m1 m2 k d d1 d2 % alle paramters
% definieer tussenvergelijkingen (alle krachten)
F_veer = (z1 -z2 ) * k;
F_damper = (dz1-dz2) * d;
F1_wrijf = dz1 * d1;
F2_wrijf = dz2 * d2;
% stel de bewegingsvergelijkingen op
ddz1 = (u1 - F_veer - F_damper - F1_wrijf)/m1;
ddz2 = ( F_veer + F_damper - F2_wrijf)/m2;
% definieer een toestand en input/output vector voor state space
x = [dz1 ; z1 ; dz2 ; z2];
u = [u1];
% stel nu f(x,u) samen
fx = [ddz1 ; dz1 ; ddz2 ; dz2 ]; % volgens: x' = f(x,u)
hx = [ z1 ; z2 ; m1*ddz1 ; m2*ddz2 ; F_veer+F_damper]; % volgens: y = h(x,u)
% dan nu de A,B,C,D matrix berekenen
A = jacobian(fx,x);
B = jacobian(fx,u);
C = jacobian(hx,x);
D = jacobian(hx,u);
% definieer de waarden voor de parameters!
k = 10;
d = 1;
d1 = 5;
d2 = 10;
m1 = 10;
m2 = 5;
% vul deze waarden in. en converteer deze van symbolic > double
A = double(eval(A));
B = double(eval(B));
C = double(eval(C));
D = double(eval(D));
% maak er een state-space van
sys = minreal(ss(A,B,C,D));
% en plot de afzonderlijke bode-diagrammen (input,output)=(1,2) betekent:
% laat de positie van massa 2 zien onder invloed van input 1, de kracht
input_output = [
1 1 ; % welke input wil je zien? [ 1 .. length(u) ]
1 2 ; % welke output wil je zien? [ 1 .. length(hx) ]
];
for i=1:length(input_output)
bode(sys(input_output(2,i),input_output(1,i))); hold on;
end
% we kunnen nu een stap op de ingang zetten van 50 N (tussen .1 en 2 s.):
dt = 1e-3;
t = [0:5000]*dt;
ut = zeros(size(t)); ut( find( (.1<t) & (t<2) ) )=50;
yt = lsim(sys,ut,t);
figure
subplot(311); plot(t,[yt(:,1) yt(:,2) ut(:)./10]); legend('z_1','z_2','0.1*F(t)',4);
subplot(312); plot(t,yt(:,2)-yt(:,1)); legend('z_2-z_1',0);
subplot(325); plot(t,[yt(:,3) yt(:,4)]); legend('z"_1','z"_2',0);
subplot(326); plot(t,[yt(:,5)]); legend('F_{veer}+F_{demper}',0);
, tussendoor werk ik met Simulink en daar loop ik nu even mee vast. Dus bij deze een vraag voor alle regelnichten onder ons:
