Practical No. 1: Iterative Calculation
a. Program for iterative calculation
Aim: Write a code in Scilab to find the sum of first n natural number.
Scilab Code:
clc
i=1;
sum=0;
n=input("Enter the no.of iterations: ");
printf("\n Iteration Sum");
while(i<=n)
sum=sum+i;
printf("\n %d %d",i,sum);
i=i+1;
disp("Result:"+string(sum));
end
Output:
Enter the no.of iterations: 5
Iteration Sum
1 1
2 3
3 6
4 10
5 15
b. Program to calculate the roots of a quadratic equation using the formula.
Aim: Write a code in Scilab to find the root of a quadratic equation \( x^2+4x+3=0 \)
Scilab Code:
clc
disp("Enter the value of");
a=input("a: ");
b=input("b: ");
c=input("c: ");
x1=(-b+sqrt(b^2-4*a*c))/2*a;
x2=(-b-sqrt(b^2-4*a*c))/2*a;
printf("x1= %.2f, x2= %.2f",x1,x2);
Output:
"Enter the value of"
a: 1
b: 4
c: 3
x1= -1.00, x2= -3.00
c. Program to evaluate \( 𝑒^𝑥 \) using infinite series.
Aim: Write a code in Scilab to evaluate the value of \( 𝑒^𝑥 \)
using infinite series
Scilab Code:
clc
x=input("x = ");
term=1;
sum=0;
n=0;
tolerance=1e-6;
while(abs(term)>tolerance)
sum=sum+term;
n=n+1;
term=term*x/n;
end
printf("e^x = %f",sum)
Output:
x = 1
e^x = 2.718282
Practical No. 2: Solution of algebraic and transcendental equations
a. Program to solve algebraic and transcendental equation by bisection method.
Aim: Write a code in Scilab to find the roots of \( x^3-x-4 \) using bisection method
Scilab Code:
clc;
clear;
function f=f(x)
f=x^3-x-4;
endfunction
a=input("a = ");
b=input("b = ");
n=input("n = ");
tolerance=1e-6;
p=0;
printf("Iteration a b c f(c)")
for i=1:n
c=(a+b)/2;
printf("\n %d %.4f %.4f %.4f %.4f\n",i,a,b,c,f(c));
if(abs(f(c))<tolerance|abs(c-p)<0.0001)
break;
end
if(f(c)*f(a)<0)then
b=c;
else
a=c;
end
p=c;
end
Output:
a = 0
b = 1
n = 5
Iteration a b c f(c)
1 0.0000 1.0000 0.5000 -4.3750
2 0.5000 1.0000 0.7500 -4.3281
3 0.7500 1.0000 0.8750 -4.2051
4 0.8750 1.0000 0.9375 -4.1135
5 0.9375 1.0000 0.9688 -4.0596
b. Program to solve algebraic and transcendental equation by false position method.
Aim: Write a code in Scilab to find the roots of \( x^3-9x+1 \) using false position method
Scilab Code:
clc;
clear;
function f=f(x)
f=x^3-9*x+1;
endfunction
a=input("a = ");
b=input("b = ");
n=input("n = ");
tolerance=1e-6;
p=0;
printf("Iteration a b x f(x)")
for i=1:n
x=(a*f(b)-b*f(a))/(f(b)-f(a));
printf("\n %d %.4f %.4f %.4f %.4f\n",i,a,b,x,f(x));
if(abs(f(x))<tolerance|abs(x-p)<0.0001)
break;
end
if(f(x)*f(a)<0)then
b=x;
else
a=x;
end
p=x;
end
Output:
a = 0
b = 1
n = 5
Iteration a b x f(x)
1 0.0000 1.0000 0.1250 -0.1230
2 0.0000 0.1250 0.1113 -0.0004
c. Program to solve algebraic and transcendental equation by Secant method.
Aim: Write a code in Scilab to find the roots of \( cos(x)+2sin(x)+x^2 \) using Secant method
Scilab Code:
clc;
clear;
function f=f(x)
f=cos(x)+2*sin(x)+x^2;
endfunction
a=input("a = ");
b=input("b = ");
n=input("n = ");
tolerance=1e-6;
p=0;
printf("Iteration a b x f(x)")
for i=1:n
x=b-((b-a)/(f(b)-f(a)))*f(b);
printf("\n %d %.4f %.4f %.4f %.4f\n",i,a,b,x,f(x));
if(abs(f(x))<tolerance|abs(x-p)<0.0001)
break;
end
if(f(x)*f(a)<0)then
b=x;
else
a=x;
end
p=x;
end
Output:
a = -1
b = 0
n = 20
Iteration a b x f(x)
1 -1.0000 0.0000 -0.8752 -0.1285
2 -0.8752 0.0000 -0.7755 -0.0847
3 -0.7755 0.0000 -0.7150 -0.0449
4 -0.7150 0.0000 -0.6842 -0.0211
5 -0.6842 0.0000 -0.6701 -0.0093
6 -0.6701 0.0000 -0.6639 -0.0040
7 -0.6639 0.0000 -0.6612 -0.0017
8 -0.6612 0.0000 -0.6601 -0.0007
9 -0.6601 0.0000 -0.6596 -0.0003
10 -0.6596 0.0000 -0.6594 -0.0001
11 -0.6594 0.0000 -0.6593 -0.0001
d. Program to solve algebraic and transcendental equation by Newton Raphson method.
Aim: Write a code in Scilab to find the roots of \( 3x-cos(x)-1 \) using Newton Raphson method
method.
Scilab Code:
clc;
clear;
function f=f(x)
f=3*x-cos(x)-1;
endfunction
a=input("a = ");
b=input("b = ");
n=input("n = ");
tolerance=1e-6;
p=0;
printf("Iteration xn x(n+1) f(x)")
for i=1:n
x=b-(f(b)/numderivative(f,b));
printf("\n %d %.4f %.4f %.4f\n",i,b,x,f(x));
if(abs(f(x))<tolerance|abs(x-p)<0.0001)
break;
end
b=x;
p=x;
end
printf("\n");
printf("The approximate value of root is: %.4f",x);
Output:
a = 0
b = 1
n = 5
Iteration xn x(n+1) f(x)
1 1.0000 0.6200 0.0462
2 0.6200 0.6071 0.0001
3 0.6071 0.6071 0.0000
The approximate value of root is: 0.6071
Practical No. 3: Interpolation
a. Program for Newton’s forward interpolation.
Aim: Write a code in Scilab to find y(0.5) using Newton’s forward interpolation.
| x | 0 | 1 | 2 | 3 |
| y | 1 | 2 | 1 | 10 |
Scilab Code:
clc;clear;
x=[0, 1, 2, 3];
y=[1, 2, 1, 10];
n=length(x);
h=x(2)-x(1);
T=x';
for i=0:n-1
T(1:n-i,i+2)=diff(y',i);
end
printf("Forward Difference Table\n");
printf("\n x y dy d2y d3y ");
disp(T);
printf("\n")
xi=input("Enter the value of xi: ");
p=(xi-x(1))/h; yi=y(1); q=p;
for j=1:n-1
yi=yi+(q*T(1,j+2))/factorial(j);
q=q*(p-j);
end
disp("Interpolated value of yi is: ",yi);
Output:
Forward Difference Table
x y dy d2y d3y
0. 1. 1. -2. 12.
1. 2. -1. 10. 0.
2. 1. 9. 0. 0.
3. 10. 0. 0. 0.
Enter the value of xi: 0.5
"Interpolated value of yi is: "
2.5
b. Program for Newton’s backward interpolation.
Aim: Write a code in Scilab to find y(7.5) using Newton’s backward interpolation.
| x | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
| y | 1 | 8 | 27 | 64 | 125 | 216 | 343 | 512 |
Scilab Code:
clc;clear;
x=[1,2,3,4,5,6,7,8];
y=[1,8,27,64,125,216,343,512];
n=length(x);
h=x(2)-x(1);
T=x';
for i=0:n-1
T(i+1:n,i+2)=diff(y',i);
end
printf("\n Forward Difference Table\n");
printf("\n x y dy d2y d3y d4y d5y d6y d7y");
disp(T);
printf("\n")
xi=input("Enter the value of xi: ");
//p=max(find(T(:,1)<xi));
p=(xi-x(n))/h; yi=y(n); q=p;
for j=1:n-1
yi=yi+(q*T(n,j+2))/factorial(j);
q=q*(p+j);
end
disp("Interpolated value of yi is: ",yi);
Output:
Forward Difference Table
x y dy d2y d3y d4y d5y d6y d7y
1. 1. 0. 0. 0. 0. 0. 0. 0.
2. 8. 7. 0. 0. 0. 0. 0. 0.
3. 27. 19. 12. 0. 0. 0. 0. 0.
4. 64. 37. 18. 6. 0. 0. 0. 0.
5. 125. 61. 24. 6. 0. 0. 0. 0.
6. 216. 91. 30. 6. 0. 0. 0. 0.
7. 343. 127. 36. 6. 0. 0. 0. 0.
8. 512. 169. 42. 6. 0. 0. 0. 0.
Enter the value of xi: 7.5
"Interpolated value of yi is: "
421.875
C. Program for Lagrange’s interpolation.
Aim: Write a code in Scilab to find y(3) using Lagrange’s interpolation formula.
| x | 1 | 2 | 5 | 6 |
| y | 10 | 15 | 22 | 35 |
Scilab Code:
clc;clear;
x=[1,2,5,6]
y=[10,15,22,35]
n=length(x);
L=0;
xi=input("Enter the value of xi: ")
for i=1:n
N=1;D=1;
for j=1:n
if (i==j)
continue;
else
N=N*(xi-x(j));
D=D*(x(i)-x(j));
end
end
L=L+(N/D)*y(i);
end
disp(L);
Output:
Enter the value of xi: 3
16.
Practical No. 4: Solving linear system of equations by iterative methods
a. Program for solving linear system of equations using Gauss Jordan method.
Aim: Write a code in Scilab to find the solution of following system of linear equations using Gauss Jordan method.
\( x+y+z=90 \)
\( 2x+3y+6z=370 \)
\( 3x-8y-4z=-340 \)
Scilab Code:
clc;clear;
a=[1,1,1;2,3,6;3,-8,-4]; b=[90;370;-340];
n=length(b);
a=[a b]
disp(a);
for j=1:n
for i=j+1:n
mf=a(i,j)/a(j,j);
a(i,:)=a(i,:)-mf*a(j,:);
end
end
//Triangularization
disp("The triangularization");
disp(a);
x=zeros(1,n);
disp("Solution values are: ")
//Backward Substitution
for i=n:-1:1
ax=sum(a(i,1:n).*x);
x(i)=(a(i,n+1)-ax)/a(i,i);
disp("x("+string(i)+") = "+string(x(i)));
end
Output:
1. 1. 1. 90.
2. 3. 6. 370.
3. -8. -4. -340.
"The triangularization"
1. 1. 1. 90.
0. 1. 4. 190.
0. 0. 37. 1480.
"Solution values are: "
"x(3) = 40"
"x(2) = 30"
"x(1) = 20"
b. Program for solving linear system of equations using Gauss Seidel method.
Aim: Write a code in Scilab to find the solution of following system of linear equations Gauss Seidel method.
\( 6x+y+z=105 \)
\( 4x+8y+3z=155 \)
\( 5x+4y-10z=-65 \)
Scilab Code:
clc;clear;
a=[6,1,1;4,8,3;5,4,-10]; b=[105;155;65];
l=length(b);
n=input("Enter the no. of iterations: ");
x=zeros(1,l);
for j=1:n
disp("After "+string(j)+" iterations, Solutions are: ");
xp=x;
for i=1:l
Sum(i)=sum(a(i,:).*x);
x(i)=(b(i)-Sum(i)+a(i,i)*x(i))/a(i,i);
disp("x("+string(i)+") =",x(i));
end
var=find(abs(x-xp)>0.0001);
if var==[]
break;
end
end
Output:
Enter the no. of iterations: 20
"After 1 iterations, Solutions are: "
"x(1) ="
17.5
"x(2) ="
10.625
"x(3) ="
6.5
"After 2 iterations, Solutions are: "
"x(1) ="
14.645833
"x(2) ="
9.6145833
"x(3) ="
4.66875
"After 3 iterations, Solutions are: "
"x(1) ="
15.119444
"x(2) ="
10.064497
"x(3) ="
5.0855208
"After 4 iterations, Solutions are: "
"x(1) ="
14.974997
"x(2) ="
9.9804311
"x(3) ="
4.9796710
"After 5 iterations, Solutions are: "
"x(1) ="
15.006650
"x(2) ="
10.004299
"x(3) ="
5.0050442
"After 6 iterations, Solutions are: "
"x(1) ="
14.998443
"x(2) ="
9.9988870
"x(3) ="
4.9987762
"After 7 iterations, Solutions are: "
"x(1) ="
15.000389
"x(2) ="
10.000264
"x(3) ="
5.0003004
"After 8 iterations, Solutions are: "
"x(1) ="
14.999906
"x(2) ="
9.9999344
"x(3) ="
4.9999267
"After 9 iterations, Solutions are: "
"x(1) ="
15.000023
"x(2) ="
10.000016
"x(3) ="
5.0000179
"After 10 iterations, Solutions are: "
"x(1) ="
14.999994
"x(2) ="
9.9999961
"x(3) ="
4.9999956
Practical No. 5: Numerical Differentiation
a. Program to obtain derivatives numerically.
Aim: Write a code in Scilab to find the f'(1.1) for the following data.
| x | 1 | 1.1 | 1.2 | 1.3 | 1.4 | 1.5 | 1.6 |
| y | 7.989 | 8.403 | 8.781 | 9.129 | 9.451 | 9.750 | 10.031 |
Scilab Code:
clc;clear
x=[1,1.1,1.2,1.3,1.4,1.5,1.6];
y=[7.989,8.403,8.781,9.129,9.451,9.750,10.031];
n=length(x);
h=x(2)-x(1);
T=x';
for i=0:n-1
T(1:n-i,i+2)=diff(y',i);
end
printf("Forward Difference Table\n");
printf("\n x y dy d2y d3y ");
disp(T);
printf("\n")
xi=input("Enter the value of xi: ");
p=max(find(T(:,1)<=xi));
yi=0;
for j=1:n-1
yi=yi+(-1)^(j+1)*(T(p,j+2)/(h*j));
end
disp("Interpolated value of y is: ",yi);
Output:
Forward Difference Table
x y dy d2y d3y d4y d5y d6y
1. 7.989 0.414 -0.036 0.006 -0.002 0.001 0.002
1.1 8.403 0.378 -0.03 0.004 -0.001 0.003 0.
1.2 8.781 0.348 -0.026 0.003 0.002 0. 0.
1.3 9.129 0.322 -0.023 0.005 0. 0. 0.
1.4 9.451 0.299 -0.018 0. 0. 0. 0.
1.5 9.75 0.281 0. 0. 0. 0. 0.
1.6 10.031 0. 0. 0. 0. 0. 0.
Enter the value of xi: 1.1
"Interpolated value of y is: "
3.9518333
Practical No. 6: Numerical Integration
a. Program for numerical integration using Trapezoidal rule.
Aim: Write a code in Scilab to find numerical integration by using trapezoidal rule for the following data:
| x | 0 | 1 | 2 | 3 | 4 | 5 | 6 |
| y | 1.000 | 0.500 | 0.333 | 0.250 | 0.200 | 0.1667 | 0.1429 |
Scilab Code:
clc;clear;
fx=input("Enter the function: ","s")
function y=f(x)
y=evstr(fx);
endfunction
a=input("Enter the value of lower limit a = ")
b=input("Enter the value of upper limit b = ")
n=input("Enter the number of iterations n = ")
h=(b-a)/n;
x=a:h:b
y=f(x);
I=y(1)+y(n+1);
for i=2:1:n
I=I+2*y(i);
end
I=(h/2)*I;
disp("The value of the integral using Trapezoidal rule is: ",I);
Output:
Enter the function: x^(-1)
Enter the value of lower limit a = 1
Enter the value of upper limit b = 2
Enter the number of iterations n = 4
"The value of the integral using Trapezoidal rule is: "
0.6970238
b. Program for numerical integration using Simpson’s 1/3rd rule.
Aim: Write a code in Scilab to find numerical integration by using Simpson’s 1/3rd rule for the following data:
| x | 0 | 1 | 2 | 3 | 4 | 5 | 6 |
| y | 1.00 | 0.500 | 0.333 | 0.250 | 0.200 | 0.1667 | 0.1429 |
Scilab Code:
clc;clear;
fx=input("Enter the function: ","s")
function y=f(x)
y=evstr(fx);
endfunction
a=input("Enter the value of lower limit a = ")
b=input("Enter the value of upper limit b = ")
n=input("Enter the number of iterations n = ")
h=(b-a)/n;
x=a:h:b
y=f(x);
I=y(1)+y(n+1);
for i=2:2:n
I=I+4*y(i);
end
for i=3:2:n-1
I=I+2*y(i);
end
I=(h/3)*I;
disp("The value of the integration by Simpsons 1/3rd rule is: ",I);
Output:
Enter the function: x^(-1)
Enter the value of lower limit a = 1
Enter the value of upper limit b = 2
Enter the number of iterations n = 4
"The value of the integration by Simpsons 1/3rd rule is: "
0.6932540
c. Program for numerical integration using Simpson’s 3/8th rule.
Aim: Write a code in Scilab to find numerical integration by using Simpson’s 3/8th rule for the following data:
| x | 0 | 1 | 2 | 3 | 4 | 5 | 6 |
| y | 1.00 | 0.500 | 0.333 | 0.250 | 0.200 | 0.1667 | 0.1429 |
Scilab Code:
clc;clear;
fx=input("Enter the function: ","s")
function y=f(x)
y=evstr(fx);
endfunction
a=input("Enter the value of lower limit a = ")
b=input("Enter the value of upper limit b = ")
n=input("Enter the number of iterations n = ")
h=(b-a)/n;
x=a:h:b
y=f(x);
I=y(1)+y(n+1);
y(1,n+1:n+2)=0;
for i=2:3:n
I=I+3*y(i)+3*y(i+1)+2*y(i+2);
end
I=(3*h/8)*I;
disp("The value of the integration by Simpsons 1/3rd rule is: ",I);
Output:
Enter the function: x^(-1)
Enter the value of lower limit a = 1
Enter the value of upper limit b = 2
Enter the number of iterations n = 4
"The value of the integration by Simpsons 1/3rd rule is: "
0.6602679
Practical No. 7: Solution of differential equations
a. Program to solve differential equations using Euler’s method.
Aim: Write a scilab code to find y'(0.3) using Euler’s method where \( \frac{dy}{dx}=1-y, \ y(0)=0 \).
Scilab Code:
//solve the ODE (1st order) using Eulers Method
clc;clear;
fxy=input("Enter the function: ","s")
function z=f(x,y)
z=evstr(fxy);
endfunction
x0=input("Enter the value of x0 = ");
y0=input("Enter the value of y0 = ");
xm=input("Enter the value of xm = ");
h=input("Take h = ");
x=x0:h:xm;
n=length(x);
y(1)=y0;
disp("x y");
disp([x0 y0]);
for i=1:n-1
y(i+1)=y(i)+h*f(x(i),y(i));
disp([x(i+1) y(i+1)]);
end
Output:
Enter the function: 1-y
Enter the value of x0 = 0
Enter the value of y0 = 0
Enter the value of xm = 0.3
Take h = 0.1
"x y"
0. 0.
0.1 0.1
0.2 0.19
0.3 0.271
