[Scilab-users] Scilab code published in scientific article

classic Classic list List threaded Threaded
1 message Options
Claus Futtrup Claus Futtrup
Reply | Threaded
Open this post in threaded view
|

[Scilab-users] Scilab code published in scientific article

Dear Scilabers

Just FYI, I had a scientific article published this weekend. You can see it  here:

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 baffle
h = 0; // bass reflex system tuning ratio
       // h = 0 for closed box systems
function y=resp(s)
y = s^2/(s^2+s/Qts+1+alpha);
endfunction // Eq (15)

N0    = 12; // resolution
steps = 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 Response
plot([0 t],[1; xs],"-xr");
xgrid(color("grey70"));
// Plot contour method numerical error
err = abs(xs - exact);
g_e = scf(); // Graph Error
plot([0 t],[0; err],"-r");
xgrid(color("grey70"));


Best regards,
Claus

_______________________________________________
users mailing list
[hidden email]
http://lists.scilab.org/mailman/listinfo/users