//This program numerically solves the state equations from the system in Lab 3
//
//by Randy Baker  9-16-03


//open file to write to

u=file('open', 'C:\temp\scilabintegration.txt', 'unknown')


// User input of system component values
m1=input("\n\nInput Value for Mass 1\n\n");
m2=input("\n\nInput Value for Mass 2\n\n");
ks1=input("\n\nInput Spring Constant 1\n\n");
ks2=input("\n\nInput Spring Constant 2\n\n");
ks3=input("\n\nInput Spring Constant 3\n\n");
kd1=input("\n\nInput Value For Damper 1\n\n");
kd2=input("\n\nInput Value For Damper 2\n\n");
kd3=input("\n\nInput Value For Damper 3\n\n");
F=input("\n\nInput Value For Applied Force\n\n");

// User input of initial conditions

x1=input("\n\nInput Initial Value of X1\n\n");
v1=input("\n\nInput Initial Value of V1\n\n");
x2=input("\n\nInput Initial Value of X2\n\n");
v2=input("\n\nInput Initial Value of V2\n\n");


// define the state matrix function
// the values returned are [x1, v1, x2, v2]
X=[x1, v1, x2, v2];	

// Define the function derivative to evaluate the derivatives of the state variables
// The values returned will be the derivatives
function derivative=f(X)
	derivative = [ X($, 2), (-(kd1+kd2)*X($,2)-(ks1+ks2)*X($,1)+kd2*X($,4)+ks2*X($,3))/m1, X($,4), (-(kd2+kd3)*X($,4)-(ks2+ks3)*X($,3)+kd2*X($,2)+ks2*X($,1)+F)/m2 ]; // use state equations to calculate derivatives
endfunction


// Set the step size for the integration
h=0.001;


//
// Loop for printing the results of integration
//
fprintf(u,' t(s)           x1            v1            x2            v2\n\n')
for t=0:h:10,



fprintf(u,'%f     %f     %f      %f     %f', t, X($,1),	X($,2), X($,3),	X($,4));
X = [X ; X($,:) + h*f(X)];

end