NSM Practical
NSM Practical
𝑥
𝑥2 𝑥3 𝑥𝑛
𝑒 = 1+𝑥+ + +. . . . . . . . . . . . . . . . +
2 3 𝑛
Evaluate e0.5 and compare with the value 1.648721. Use six terms to
evaluate each series and compute true and approximate relative errors
(upto three significant figures) as terms are added.
Scilab Code :
clc;
n=3;//number of significant figutes
es=0.5*(10^(2-n));//percent,specified error criterion
x=0.5;
f(1)=1; //first estimate f=e^x=1
ft=1.648721; //true value of e^0.5=f
et(1)=(ft-f(1))*100/ft;
ea(1)=100;
i=2;
while ea(i-1)>=es
f(i)=f(i-1)+(x^(i-1)) //(factorial(i-1));
et(i)=(ft-f(i))*100/ft;
ea(i)=(f(i)-f(i-1))*100/f(i);
i=i+1;
end
for j=1:i-1
disp(ea(j),"Approximate estimate of error(%)=",et(j),"True % relative error
=",f(j),"Result=",j,"term number=")
disp(".........................................")
end
Output:
term number=
1.
Result=
1.
39.346924
100.
.........................................
term number=
2.
Result=
1.5
9.0203861
33.333333
.........................................
term number=
3.
Result=
1.75
-6.1428829
14.285714
.........................................
term number=
4.
Result=
1.875
-13.724517
6.6666667
.........................................
term number=
5.
Result=
1.9375
-17.515335
3.2258065
.........................................
term number=
6.
Result=
1.96875
-19.410743
1.5873016
.........................................
term number=
7.
Result=
1.984375
-20.358448
0.7874016
.........................................
term number=
8.
Result=
1.9921875
-20.8323
0.3921569
.........................................
term number=
9.
Result=
1.9960938
-21.069226
0.1956947
.........................................
term number=
10.
Result=
1.9980469
-21.187689
0.0977517
.........................................
term number=
11.
Result=
1.9990234
-21.24692
0.048852
.........................................
Practical No 1 b.: Program to calculate the roots of quadratic equation using the formula.
Problem Statement: Write a Scilab code to evaluate roots of quadratic equation x2-5x+6=0.
Scilab Code :
clc;
response=1;
while response==1
a=input("input the value of a:")
b=input("input the value of b:")
c=input("input the value of c:")
if a==0 then
if b~=0 then
r1=-c/b;
disp(r1,"the root:")
else disp("trivial solution")
end
else
discr=b^2-4*a*c;
if discr>=0 then
r1=(-b+sqrt(discr))/(2*a);
r2=(-b-sqrt(discr))/(2*a);
disp(r2,"and",r1,"the roots are:")
else
r1=-b/(2*a)
r2=r1;
i1=sqrt(abs(discr))/(2*a);
i2=-i1;
disp(r2+i2*sqrt(-1),r1+i1*sqrt(-1),"the roots are:")
end
end
response=input ("Do you want to continue (press 1 for yes and 2 for no)?")
if response==2 then exit;
end
end
Output:
input the value of a:1
a =
1.
input the value of b:-5
b =
-5.
input the value of c:6
c =
6.
3.
and
2.
Do you want to continue (press 1 for yes and 2 for no)?
Practical No 3: Solution of algebraic and transcendental equation by bisection method.
a.
Problem Statement: Write a Scilab code to find the real root of the equation x3-x-1=0 using bisection
method correct to four places of decimal.
Scilab Code :
clc;
deff('y=f(x)','y=x^3-x-1');
c=1;
m=(x1+x2)/2;
if f(m)*f(x1)>0
x1=m;
else
x2=m;
end
end
Output:
successive approximations
x1 x2 m f(m)
Problem Statement: Write a Scilab code to find the real root of the equation 2x3-x-5=0 using false
position method correct to five places of decimal.
Scilab Code :
clc;
deff('y=f(x)','y=x^3-2*x-5');
for i=1:25
x1=b*f(a)/(f(a)-f(b))+a*f(b)/(f(b)-f(a))
if (f(a)*f(x1))>0
b=x1;
else
a=x1;
end
if abs(f(x1))<d
break
end
printf('\t%f %f %f %f %f\n',a,b,f(a),f(b),x1);
end
successive iterations
a b f(a) f(b) x1
Problem Statement: Write a Scilab code to find the real root of the equation x3-2x-5=0 using Secant
method correct to six places of decimal.
Scilab Code :
clc;
deff('y=f(x)','y=x^3-2*x-5');
x1=2,x2=3//initial values
n=1;
c=0;
printf('successive iterations\nx1\tx2\tx3\tf(x3)\n')
while n==1
x3=(x1*f(x2)-x2*f(x1))/(f(x2)-f(x1));
printf('\t%f\t%f\t%f\t%f\n',x1,x2,x3,f(x3));
if f(x3)*f(x1)>0 then
x2=x3;
else
x1=x3;
end
if abs(f(x3))<0.000001 then
break;
end
c=c+1;
end
successive iterations
x1 x2 x3 f(x3)
x 1 3 5 7
f(x) 24 120 336 720
Scilab Code:
clc;
x=[1 3 5 7];
y=[24 120 336 720];
h=2
c=1;
for i=1:3
d1(c)=y(i+1)-y(i);
c=c+1;
end
c=1;
for i=1:2
d2(c)=d1(i+1)-d1(i);
c=c+1;
end
c=1;
for i=1:1
d3(c)=d2(i+1)-d2(i);
c=c+1;
end
d=[d1(1) d2(1) d3(1)];
x0=8;
pp=1;
y_x=y(1);
p=(x0-1)/2;
for i=1:3
pp=1;
for j=1:i
pp=pp*(p-(j-1))
end
y_x=y_x+(pp*d(i))/factorial(i);
end
printf('value of function at %f is:%f',x0,y_x);
Output:
Problem statement: Write a Scilab code to find Sin (38°) using Newton’s backward difference
interpolation formula for the following data.
x(in degrees) 15 20 25 30 35 40
Scilab code:
clc;
x=[15 20 25 30 35 40];
h=5
c=1;
for i=1:5
d1(c)=y(i+1)-y(i);
c=c+1;
end
c=1;
for i=1:4
d2(c)=d1(i+1)-d1(i);
c=c+1
end
c=1;
for i=1:3
d3(c)=d2(i+1)-d2(i);
c=c+1;
end
c=1;
for i=1:2
d4(c)=d3(i+1)-d3(i);
c=c+1;
end
c=1;
for i=1:1
d5(c)=d4(i+1)-d4(i);
c=c+1 ;
end
x0=38;
pp=1;
y_x=y(6);
p=(x0-x(6))/h;
for i=1:5
pp=1;
for j=1:i
pp=pp*(p-(j-1))
end
y_x=y_x+(pp*d(i))/factorial(i);
end
Output:
If y1 =4, y3 = 12, y4 =19 and yx =7, find x using lagrange’s interpolation formula.
clc;
y=[4 12 19];
x=[1 3 4];
y_x=7;
Y_X=0;
poly(0,'y');
for i=1:3
p=x(i);
for j=1:3
if i~=j then
p=p*((y_x-y(j))/(y(i)-y(j)))
end
end
Y_X=Y_X+p;
end
disp(Y_X,'Y_X=');
Output:
Y_X=
1.8571429
Practical No. 4
a.
Aim : Program for solving linear system of equations using Gauss Jordan method
2x+y+z=10
3x+2y+3z=18
x+4y+9z=16
code
clc;
A=[2,1,1,10;3,2,3,18;1,4,9,16];
for i=1:3
j=i
while (A(i,i)==0&j<=3)
for k=1:4
B(1,k)=A(j+1,k)
A(j+1,k)=A(i,k)
A(i,k)=B(1,k)
end
disp(A);
j=j+1;
end
disp(A);
for k=4:-1:i
A(i,k)=A(i,k)/A(i,i)
end
disp(A)
for k=1:3
if (k~=i) then
l=A(k,i)/A(i,i)
for m=i:4
A(k,m)=A(k,m)-l*A(i,m)
end
end
end
disp(A)
end
for i=1:3
printf('\nx(%i)=%g\n',i,A(i,4))
end
Output:
x(1)=7
x(2)=-9
x(3)=5
b. Program for solving linear system of equations using Gauss Seidel method.
2x-y-0.2z=7.85
-x+2y-z=-19.3
-y+2z=71.4
code:
clc;
function [X]=gaussseidel(A, n, N, X, b)
L=A;
U=A;
D=A;
for i=1:1:n
for j=1:1:n
D(i,j)=0;
end
D(i,j)=0;
end
end
end
end
for k=1:1:N
X=(D+L)^-1*(-U*X+b);
disp(X)
end
endfunction
A=[2 -1 0; -1 2 -1; 0 -1 2]
b=[7;1;1]
N=3;
n=3;
X=[0;0;0]
gaussseidel(A,n,N,X,b)
Output:
A =
2. - 1. 0.
- 1. 2. - 1.
0. - 1. 2.
b =
7.
1.
1.
X =
0.
0.
0.
1.75
0.6875
0.421875
1.046875
0.2734375
0.1074219
1.2949219
0.4638672
0.3122559
ans =
1.2949219
0.4638672
0.3122559
Practical No. 5
Programming to obtain derivatives numerically
Problem statement : Write and execute Scilab code for the following:
𝑑𝑦 𝑑2 𝑦
From the following table of values of x and y obtain 𝑑𝑥 and 𝑑𝑥 2 for x=1.2
Code:
clc;
x=[1.0 1.2 1.4 1.6 1.8 2.0 2.2];
y=[2.7182 3.3201 4.552 4.9530 6.0496 7.3891 9.0250];
c=1;
for i=1:6
d1(c)=y(i+1)-y(i);
c=c+1;
end
c=1;
for i=1:5
d2(c)=d1(i+1)-d1(i);
c=c+1;
end
c=1;
for i=1:4
d3(c)=d2(i+1)-d2(i);
c=c+1;
end
c=1;
for i=1:3
d4(c)=d3(i+1)-d3(i);
c=c+1;
end
c=1;
for i=1:2
d5(c)=d4(i+1)-d4(i);
c=c+1;
end
c=1;
for i=1:1
d6(c)=d5(i+1)-d5(i);
c=c+1;
end
x0=1.2
h=0.2;
f1=((d1(2)-d2(2)/2+d3(2)/3-d4(2)/4+d5(2)/5)/h);
printf('the first derivative of function at 1.2 is %f\n',f1)
f2=(d2(2)-d3(2)+(11*d4(2))/12-(5*d5(2))/6)/h^2;
printf('the second derivative of function at 1.2 is %f\n',f2)
OUTPUT:
the first derivative of function at 1.2 is 15.740317
the second derivative of function at 1.2 is -156.070833
Practical No. 6
a. Program for numerical integration using Trapezoidal rule.
Problem Statement :Write and execute scilab code for the following:
1
Evaluate ∫0 𝑥 4 dx, with n=6. Using Trapezoidal rule.
Code:
clc;
deff('y=f(x)','y=x^4')
a=input("Enter Lower Limit: ")
b=input("Enter Upper Limit: ")
n=input("Enter number of sum intervals: ")
h=(b-a)/n
add1=0
add2=0
for i=0:n
x=a+i*h
y=f(x)
disp([x y])
if(i==0)|(i==n)
add1=add1+y
else
add2=add2+y
end
end
I=(h/2)*(add1+2*add2)
disp(I,"Integration by Trapezoidal Rule is:")
Output:
Enter Lower Limit: 0
0. 0.
0.1666667 0.0007716
0.3333333 0.0123457
0.5 0.0625
0.6666667 0.1975309
0.8333333 0.4822531
1. 1.
0.2092335
b. Program for numerical integration using Simpson’s 1/3rd rule.
Problem Statement: Write and execute Scilab code for the following:
1
Evaluate ∫0 x 4 dx, with n=6, Simpson’s 1/3rd rule.
Code:
deff('y=f(x)','y=x^4')
a=input("Enter Lower Limit: ")
b=input("Enter Upper Limit: ")
n=input("Enter number of sum intervals: ")
h=(b-a)/n
add1=0
add2=0
add3=0
for i=0:n
x=a+i*h
y=f(x)
disp([x y])
if (i==0)|(i==n) then
add1=add1+y
else if (modulo(i,2)==0) then
add2=add2+y
else
add3=add3+y
end
end
end
I=(h/3)*(add1+2*add2+4*add3)
disp(I,"Integration by Simpsons (1/3)rd Rule is:")
Output:
0. 0.
0.1666667 0.0007716
0.3333333 0.0123457
0.5 0.0625
0.6666667 0.1975309
0.8333333 0.4822531
1. 1.
0.2001029
c. Program for numerical integration using Simpson’s 3/8th rule.
Problem statement: Write and execute Scilab code for the following.
1
Evaluate ∫0 𝑥 4 dx, with n=6. Using Simpson’s 3/8th rule.
Code:
deff('y=f(x)','y=x^4')
a=input("Enter Lower Limit: ")
b=input("Enter Upper Limit: ")
n=input("Enter number of sum intervals: ")
h=(b-a)/n
add1=0
add2=0
add3=0
for i=0:n
x=a+i*h
y=f(x)
disp([x y])
if (i==0)|(i==n) then
add1=add1+y
else if (modulo(i,3)==0) then
add2=add2+y
else
add3=add3+y
end
end
end
I=((3*h)/8)*(add1+2*add2+3*add3)
disp(I,"Integration by Simpsons (3/8)th Rule is:")
Output:
0. 0.
0.1666667 0.0007716
0.3333333 0.0123457
0.5 0.0625
0.6666667 0.1975309
0.8333333 0.4822531
1. 1.
0.2002315
Practical no. 7
Program statement: Write and execute Scilab code for the following:
Solve the following differential equation to find y(0.2) using Euler’s method.
𝑑𝑦
𝑑𝑥
= y-xy+x, y(0) = 1.
Code:
1.432
1.68384
1.9385472
2.1760891
b. Program to solve differential equation using modified Euler’s method
Program Statement: Write and execute Scilab code for the following:
Determine the value of y when x=0.1 using Modified Euler’s method given that y(0) =1 and
2
y’=𝑥 + 𝑦 take h=0.05.
Code:
Output
2.3096946
c. Program to solve differential equation using Runge- kutta 2nd order and 4th order
methods
Problem Statements: Write and execute Scilab Code for the following:
Determine the value of y when x=0.1 using Runge - Kutta 2nd order and 4th order methods
given that y(0) = 2 and y’ = y-x take h=0.1.
Code:
function [y]=f(a, b)
y=b-a;
endfunction
x0=0;
y0=2;
h=0.1;
for n=1:4
k1=h*f(x0,y0);
k2=h*f(x0+h,y0+k1);
y0=y0+(k1+k2)/2;
x0=x0+h;
printf('values of x0=%g\t and y0=%g\n',x0,y0);
end
function [y]=f(a, b)
y=b-a;
endfunction
x0=0;
y0=1;
h=0.25;
for n=1:4
k1=h*f(x0,y0);
k2=h*f(x0+h/2,y0+k1/2);
k3=h*f(x0+h/2,y0+k2/2);
k4=h*f(x0+h,y0+k3);
y0=y0+(k1+2*k2+2*k3+k4)/6;
x0=x0+h;
printf('values of x0=%g\t and y0=%g\n',x0,y0);
end
Output:
values of x0=0.1 and y0=2.205
values of x0=0.2 and y0=2.42103
values of x0=0.3 and y0=2.64923
values of x0=0.4 and y0=2.8909
values of x0=0.25 and y0=1.25
values of x0=0.5 and y0=1.5
values of x0=0.75 and y0=1.75
values of x0=1 and y0=2
Practical No. 8
Problem statement: Write and execute Scilab code for the following:
x 1 2 3 4 5 6 7
Code:
x=[1,2,3,4,5,6,7];
y=[0.5,2.5,2,4,3.5,6,5.5];
n=7;
s=0;
xsq=0;
xsum=0;
ysum=0;
for i =1.7
s=s+(det(x(1,i)))*(det(y(1,i)));
xsq=xsq+(det(x(1,i))^2);
xsum=xsum+det(x(1,i));
ysum=ysum+det(y(1,i));
end
disp(s,"sum of product of x and y=")
disp(xsq,"sum of square of x=")
disp(xsum,"sum of all the x=")
disp(ysum,"sum of all the y=")
a=xsum/n;
b=xsum/n;
a1=(n*s-xsum*ysum)/(n*xsq-xsum^2);
a0=b-a*a1;
disp(a1,"a1=")
disp(a0,"a0=")
disp("The equation of the line obtained is y=a0+a1*x")
Output:
sum of product of x and y=
0.5
sum of square of x=
1.
sum of all the x=
1.
sum of all the y=
0.5
a1=
0.5
a0=
0.0714286
The equation of the line obtained is y=a0+a1*x