Dear Scilabers
Just FYI, I had a scientific article published this weekend. You can see
it here:
http://www.aes.org/e-lib/browse.cfm?elib=19566
I'm not allowed to freely copy the article, but I hope nobody minds that I
show the Scilab code included in the article. Here it is:
Qts = 0.5;alpha = 0; // alpha = Vas/Vb
// alpha = 0, open baffleh = 0; // bass reflex system tuning ratio
// h = 0 for closed box systemsfunction y=resp(s)y =
s^2/(s^2+s/Qts+1+alpha);endfunction // Eq (15)
N0 = 12; // resolutionsteps = 8; // steps per unit time,t_max = 5;
// time (dimensionless)t_step = 1/steps;t = t_step:t_step:t_max;
mu_c = max([1 h]);t_c = %pi*N0/(12*mu_c);nt = length(t);
mprintf(" t - Contour");mprintf(" method - Exact\n");for i=1:nt do
ti = t(i);
// Selection defined in
// Eqs (31) and (32)
if ti < t_c then
mu = %pi*N0/(12*ti);
N = N0;
d = 3/N;
else
mu = mu_c;
N = ceil(N0*ti/t_c);
d = 3/N;
end
// Perform sum in Eq (24)
xsk = 0.0;
for k=-N:N do
u = k*d;
s = mu*(%i*u+1)^2;
dsdu = 2.0*mu*(%i-u);
xsk = xsk+imag(exp(s*ti)* ..
resp(s)/s*dsdu);
end
xsk = xsk*d/(2*%pi);
xs(i) = xsk; // Collect x-values
exact(i) = exp(-ti)*(1-ti);
// Collect exact result, compare
mprintf("%5.4e ",ti);
mprintf("%+12.11e ",xsk);
mprintf("%+12.11e\n",exact(i));end
g_tr = scf(); // Graph Time Responseplot([0 t],[1;
xs],"-xr");xgrid(color("grey70"));// Plot contour method numerical
errorerr = abs(xs - exact);g_e = scf(); // Graph Errorplot([0 t],[0;
err],"-r");xgrid(color("grey70"));
Best regards,
Claus
_______________________________________________
users mailing list
users@lists.scilab.org
http://lists.scilab.org/mailman/listinfo/users
