Numerical Method’s Program codes
Numerical Method’s Program codes
Horner’s Method Numerical Method
# include<conio.h>
# include<math.h>
# define E 0.0001
struct comp
{
float x,y;
};
struct comp mult( struct comp a, struct comp b)
{
struct comp c;
c.x = (a.x*b.x- a.y*b.y); // x for real & y for imaginary
c.y = (a.y*b.x+ a.x*b.y);
return c;
}
struct comp add ( struct comp a, struct comp b)
{
struct comp c;
c.x = (a.x+b.x);
c.y = (a.y+b.y);
return c;
}
struct comp div ( struct comp a, struct comp b)
{
struct comp c;
c.x = (a.x*b.x+a.y*b.y)/(b.x*b.x+b.y*b.y) ;
c.y = (a.y*b.x-a.x*b.y)/(b.x*b.x+b.y*b.y) ;
return c;
}
struct comp sub ( struct comp a, struct comp b)
{
struct comp c;
c.x = (a.x-b.x);
c.y = (a.y-b.y);
return c;
}
void main ()
{
int n,i;
float s;
struct comp a[10], b[10], c[10],tem,x0,x1;
clrscr ();
printf("Enter the degree of the polynomial \n");
scanf("%d",&n);
printf("Enter the n+1 cofficients\n");
for (i=0;i<=n;i++)
scanf ("%f %f",&a[i].x, &a[i].y);
printf ("Enter the initial case x0\n");
scanf ("%f %f", &x0.x,&x0.y);
do
{
printf("\n%d",n);
b[0]=a[0];
c[0]=a[0];
do
{
for(i=1; i<=n; i++)
{
b[i]= add ( mult( x0, b[i-1] ), a[i]);
printf("\n%d",n);
}
for(i=1; i<n; i++)
c[i]= add ( mult( x0, c[i-1] ), b[i]);
x1 = sub ( x0, div( b[n], c[n-1] ));
tem.x=x1.x-x0.x;
tem.y=x1.y-x0.y;
s= sqrt ( tem.x*tem.x + tem.y*tem.y );
x0=x1;
}while (s>E);
printf ("The root is %f+%fi\n",x0.x,x0.y);
n--;
for(i=0;i<=n;i++)
b[i]=a[i];
printf("%d",n);
}while(n!=0);
getch ();
}
Lagrange Interpolation Numerical Method
#include<stdio.h>#include<conio.h>
#include<math.h>
void main()
{
clrscr();
float x[10],f[10];
int n;
float sum=0;
float finalprod=0;
float prod,m;
printf("how many datas do you have?\n");
scanf("%d",&n);
for(int i=0; i<n; i++)
{
printf("x[%d]:",i);
scanf("%f",&x[i]);
printf("f[%d]:",i);
scanf("%f",&f[i]);
}
printf("Input the estimating point(x):");
scanf("%f",&m);
for(int l=0; l<n; l++)
{
prod=1;
for(int k=0; k<n; k++)
{
if(l!=k)
prod*=(m-x[k])/(x[l]-x[k]);
}
finalprod=prod*f[l];
sum+=finalprod;
}
printf("the value of f(x) is : %f",sum);
getch();
}
Divided difference Method in Numerical Method
#include<stdio.h>#include<conio.h>
main()
{
clrscr();
int i,j,n;
float xp,fp,sum,pi,x[10],a[10],f[10][10];
printf("\n Input number of data points \n");
scanf("%d",&n);
printf("\n Input values of x ---f(x) \n");
for(i=0;i<n;i++)
{
scanf("%f",&x[i]);
scanf("%f",&f[i][0]);
}
/*construct difference table*/
for(j=1;j<n;j++)
for(i=0;i<=n-j;i++)
f[i][j]=(f[i+1][j-1]-f[i][j-1])/(x[i+j]-x[i]);
/* set the coefficients of interpolation polynomial*/
for(j=0;j<n;j++)
a[j]=f[0][j];
/*compute interpolation value*/
printf("\n input xp where interpolation is required \n");
scanf("%f",&xp);
sum=a[0];
for(i=1;i<n;i++)
{
pi=1.0;
for(j=0;j<i;j++)
pi=pi*(xp-x[j]);
sum=sum+a[i]*pi;
}
fp=sum;
/*write results*/
printf("\n");
printf(" newton interpolation\n");
printf("interpolated value\n");
printf("at x=%f is %f \n",xp,fp);
getch();
return(0);
}
Newtons polynomial for equally spaced data in Numerical Method
#include <stdio.h>#include <conio.h>
#include <math.h>
long fact(int);
float facta(float,int);
main()
{
float y,h,s,sum,x[10],f0[10],f[10][10];
int n,i,j,k;
clrscr();
x[0]=0.0;f[0][0]=0.0;
x[1]=0.2;f[1][0]=0.203;
x[2]=0.4;f[2][0]=0.423;
x[3]=0.6;f[3][0]=0.684;
x[4]=0.8;f[4][0]=1.030;
x[5]=1.0;f[5][0]=1.557;
x[6]=1.2;f[6][0]=2.572;
n=6;
for(i=1;i<=n;i++)
{
for(j=0;j<=n-i;j++)
f[j][i]=f[j+1][i-1]-f[j][i-1];
}
printf("Table for equally spaced data \n");
for(i=0;i<=n;i++)
{
printf("%.3f",x[i]);
for(j=0;j<=n-i;j++)
printf("\t%.3f",f[i][j]);
printf("\n");
}
for(k=1;k<=n;k++)
f0[k]=f[0][k];
y=0.5;
h=(x[1]-x[0]);
s=(y-x[0])/h;
sum=f[0][0];
for(i=1;i<=n;i++)
sum=sum+facta(s,i)*f0[i]/fact(i);
printf("\n at x=0.5 the interpolated value=%f",sum);
getch();
return 0;
}
long fact(int n)
{
int i;
long fact;
fact=1;
for(i=1;i<=n;i++)
fact=fact*i;
return(fact);
}
float facta(float s,int i)
{
float prod=1;
int k;
for(k=0;k<=i-1;k++)
prod=prod*(s-k);
return(prod);
}
Linear Regression (Least square)in Numerical Method
#include<conio.h>#include<math.h>
#include<stdio.h>
void main()
{
int n,i;
float x[10],y[10],a,b,sumx=0,sumy=0,sumx2=0,sumxy=0;
clrscr();
printf("enter the no of data");
scanf("%d",&n);
printf("enter the values of x");
for( i=1;i<=n;i++)
scanf("%f",&x[i]);
printf("enter the values of y");
for( i=1;i<=n;i++)
scanf("%f",&y[i]);
for( i=1;i<=n;i++)
{
sumx=sumx+x[i];
sumy=sumy+y[i];
sumx2=sumx2+pow(x[i],2);
sumxy=sumxy+x[i]*y[i];
}
b=((n*sumxy)-(sumx*sumy))/((n*sumx2)-(pow(sumx,2)));
a=(sumy-b*sumx)/n;
printf("y=%fx+%f\n",b,a);
printf("for checking the function\n");
for(i=1;i<=n;i++)
{
y[i]=1.6+1.2*x[i];
printf("%f\n",y[i]);
}
getch();
}
Transcedental (for y=a*x^b)in Numerical Method
#include<stdio.h>#include<conio.h>
#include<math.h>
void main()
{
int i,n;
float A, x[10],y[10],a,b,sumlnx=0,sumlny=0,sumlnxlny=0,sumlnx2=0;
clrscr();
printf("enter the no of data\n");
scanf("%d",&n);
printf("enter the values of x\n");
for(i=1;i<=n;i++)
scanf("%f",&x[i]);
printf("enter the value of y\n");
for(i=1;i<=n;i++)
scanf("%f",&y[i]);
for(i=1;i<=n;i++)
{
sumlnx=sumlnx+log(x[i]);
sumlny=sumlny+log(y[i]);
sumlnxlny=sumlnxlny+log(x[i])*log(y[i]);
sumlnx2=sumlnx2+pow(log(x[i]),2);
}
b=(n*sumlnxlny-(sumlnx*sumlny))/(n*sumlnx2-pow(sumlnx,2));
A=(sumlny-b*sumlnx)/n;
a=exp(A);
printf("required equation is:\n");
printf("y=%fx^%f\n",a,b);
printf("\nfor checking");
for(i=1;i<=n;i++)
{
y[i]=a*pow(x[i],b);
printf("%f\n",y[i]);
}
getch();
}
Trapezoidal Method in Numerical Method
#include<conio.h>#include<stdio.h>
#include<math.h>
#define f(x) sqrt(sin(x))
void main()
{
float a,b,it,h,f1,f2;
clrscr();
printf("enter two sampling points:\n");
scanf("%f%f",&a,&b);
f1=f(a);
f2=f(b);
h=b-a;
it=h*(f1+f2)/2;
printf("the value is: %f",it);
getch();
}
Composite Trapezoidal
#include<conio.h>#include<stdio.h>
#include<math.h>
#define f(x) sqrt(sin(x))
void main()
{
int n,i;
float a,b,fa,fb,it,h,x[10],sumf1=0;
clrscr();
printf("enter two sampling points: a & b\n");
scanf("%f%f",&a,&b);
printf("enter the no of segment:\n");
scanf("%d",&n);
h=(b-a)/n;
fa=f(a);
fb=f(b);
for(i=1;i<=n-1;i++)
{
x[i]=a+i*h;
sumf1=sumf1+f(x[i]);
}
it=h*(fa+fb+2*sumf1)/2;
printf("the value is: %f\n",it);
getch();
}
//Simpson’s 1/3 rd Method
#include<stdio.h>
#include<conio.h>
#include<math.h>
#define f(x) exp(x)
void main()
{
int i;
float h,fa,fb,f1,it,a,b,x1;
clrscr();
printf("enter the two sampling points:\n");
scanf("%f%f",&a,&b);
h=(b-a)/2;
x1=(a+b)/2;
fa=f(a);
fb=f(b);
f1=f(x1);
it=h*(fa+fb+4*f1)/3;
printf("value is:%f",it);
getch();
}
Composite Simpson’s 1/3 rd Method in Numerical Method
#include<stdio.h>#include<conio.h>
#define f(x) 1/(1+x)
void main()
{
clrscr();
float a,l,b,h,t,x1,s,f1=0,f2=0;
int n,i,j;
printf("enter the limit of integration:\n");
scanf("%f%f",&a,&b);
printf("\nenter the no of segment:");
scanf("%d",&n);
h=(b-a)/n;
if(n%2==0)
{
t=0;
l=b;
j=n;
}
else
{
x1=a+(n-1)*h;
t=h/2*(f(x1)+f(b));
l=x1;
j=n-1;
}
for(i=1;i<j;i++)
{
if(i%2==0)
{
f1=f1+f(a+i*h);
}
else
{
f2=f2+f(a+i*h);
}
}
s=(h/3*(f(a)+f(l)+4*f2+2*f1))+t;
printf("com sim:%f",s);
getch();
}
Composite Simpson’s 3/8 th Method in Numerical Method
#include<stdio.h>#include<conio.h>
#define f(x) 1/(1+x)
void main()
{
clrscr();
float a,s1,l,b,h,t,x1,x2,fx=0,s,f1=0,f2=0;
int n,i,j,choice;
printf("enter the limit of integration:\n");
scanf("%f%f",&a,&b);
printf("\nenter the no of segment:");
scanf("%d",&n);
h=(b-a)/n;
if(n%3==0)
{
t=0;
s=0;
l=b;
j=n;
}
else if(n%3==2)
{
t=0;
x1=a+(n-2)*h;
x2=a+(n-1)*h;
s=h/3*(f(x1)+f(b)+4*f(x2));
l=x1;
j=n-2;
}
else
{
x1=a+(n-1)*h;
s=0;
t=h/2*(f(x1)+f(b));
l=x1;
j=n-1;
}
for(i=1;i<j;i++)
{
if(i%3==0)
{
f1=f1+f(a+i*h);
}
else
{
f2=f2+f(a+i*h);
}
}
s1=3*h/8*(f(a)+f(l)+3*f2+2*f1)+t+s;
printf("com sim:%f",s1);
getch();
}
All Method of Integration (Combined)in Numerical Method
#include<stdio.h>#include<conio.h>
#define f(x) 1/(1+x)
void main()
{
clrscr();
float a,s1,l,b,h,t,m,x1,x2,fx=0,s,f1=0,f2=0;
int n,i,j,choice;
printf("enter the limit of integration:\n");
scanf("%f%f",&a,&b);
printf("\nenter the no of segment:");
scanf("%d",&n);
printf("enter your own choice from 1 to 6:");
scanf("%d",&choice);
m=f(a)+f(b);
switch(choice)
{
case 1:
h=b-a;
t=h*m/2;
printf("using trapezoidal rule:%f",t);
break;
case 2:
h=(b-a)/n;
for(i=1;i<=n-1;i++)
{
fx+=f(a+i*h);
}
t=h*(m+(2*fx))/2;
printf("comtrap:%f",t);
break;
case 3:
h=(b-a)/2;
t=h*(m+4*f(a+h))/3;
printf("simpson's 1/3 rule :%f",t);
break;
case 4:
h=(b-a)/n;
if(n%2==0)
{
t=0;
l=b;
j=n;
}
else
{
x1=a+(n-1)*h;
t=h/2*(f(x1)+f(b));
l=x1;
j=n-1;
}
for(i=1;i<j;i++)
{
if(i%2==0)
{
f1=f1+f(a+i*h);
}
else
{
f2=f2+f(a+i*h);
}
}
s=(h/3*(f(a)+f(l)+4*f2+2*f1))+t;
printf("com sim:%f",s);
break;
case 5:
h=(b-a)/3;
t=3*h/8*(m+3*(f(a+h)+f(a+2*h)));
printf("simpson's 3/8 :%f",t);
break;
case 6:
h=(b-a)/n;
if(n%3==0)
{
t=0;
s=0;
l=b;
j=n;
}
else if(n%3==2)
{
t=0;
x1=a+(n-2)*h;
x2=a+(n-1)*h;
s=h/3*(f(x1)+f(b)+4*f(x2));
l=x1;
j=n-2;
}
else
{
x1=a+(n-1)*h;
s=0;
t=h/2*(f(x1)+f(b));
l=x1;
j=n-1;
}
for(i=1;i<j;i++)
{
if(i%3==0)
{
f1=f1+f(a+i*h);
}
else
{
f2=f2+f(a+i*h);
}
}
s1=3*h/8*(f(a)+f(l)+3*f2+2*f1)+t+s;
printf("com sim:%f",s1);
break;
}
getch();
}
Gauss Elimination Method in Numerical Method
#include<stdio.h>#include<conio.h>
void main()
{
clrscr();
int n,i,j,k;
float a[10][10],x[10],sum,u;
printf("enter the order of matrix:");
scanf("%d",&n);
printf("\ncoefficient of matrix:");
for(i=1;i<=n;i++)
{
for(j=1;j<=n+1;j++)
{
scanf("%f",&a[i][j]);
}
}
printf("the given matrix :\n");
for(i=1;i<=n;i++)
{
for(j=1;j<=n+1;j++)
{
printf("%f",a[i][j]);
printf("\t");
}
printf("\n");
}
for(k=1;k<=n-1;k++)
{
for(i=k+1;i<=n;i++)
{
u=a[i][k]/a[k][k];
for(j=k;j<=n+1;j++)
{
a[i][j]=a[i][j]-u*a[k][j];
}
}
}
printf("upper tri matrix:\n");
for(i=1;i<=n;i++)
{
for(j=1;j<=n+1;j++)
{
printf("%f",a[i][j]);
printf("\t");
}
printf("\n");
}
printf("\nthe values of variables are:");
x[n]=a[n][n+1]/a[n][n];
printf("x[%d]=%f\n",n,x[n]);
for(i=n-1;i>=1;i--)
{
sum=0;
for(j=i+1;j<=n;j++)
{
sum=sum+a[i][j]*x[j];
}
x[i]=(a[i][n+1]-sum)/a[i][i];
printf("x[%d]=%f\n",i,x[i]);
}
getch();
}
Gauss Jordan in Numerical Method
#include<stdio.h>#include<math.h>
#include<conio.h>
void main()
{ clrscr();
int n,i,j,k;
float pivot,factor,a[10][10],b[10],x[10];
printf("Enter number of equation:");
scanf("%d",&n);
for(i=0;i<n;i++)
{
printf("coeficients of aug matrix:");
for(j=0;j<=n;j++)
{
scanf("%f",&a[i][j]);
}
}
for(i=0;i<n;i++)
{
pivot=a[i][i];
for(k=0;k<=n;k++)
{
a[i][k]=a[i][k]/pivot;
}
for(j=0;j<n;j++)
{
if(j!=i)
{
factor=a[j][i]/a[i][i];
for(k=0;k<=n;k++)
{
a[j][k]=a[j][k]-factor*a[i][k];}
}
}
}
printf("roots:");
for(i=0;i<n;i++)
{
printf("%f\n",a[i][n]);
}
getch();
}
Inverse of the matrix in Numerical Method.
#include<stdio.h>#include<conio.h>
#include<math.h>
void main()
{
int i,j,k,n;
float a[10][10],l,u,norm;
clrscr();
printf("Enter the order of matrix:\n");
scanf("%d",&n);
printf("Enter the coeff. a[i][j]:\n");
for(i=1;i<=n;i++)
{
for(j=1;j<=n;j++)
{
scanf("%f",&a[i][j]);
}
for(j=n+1;j<=2*n;j++)
{
if(j==n+i)
{
a[i][j]=1;
}
else
{
a[i][j]=0;
}
}
}
//making upper triangular matrix//
for(k=1;k<=n-1;k++)
{
for(i=k+1;i<=n;i++)
{
u=a[i][k]/a[k][k];
for(j=k;j<=2*n;j++)
{
a[i][j]=a[i][j]-(u*a[k][j]);
}
}
}
making diagonal matrix in Numerical Method
for(k=n;k>=2;k--){
for(i=k-1;i>=1;i--)
{
l=a[i][k]/a[k][k];
for(j=k;j<=2*n;j++)
{
a[i][j]=a[i][j]-(l*a[k][j]);
}
}
}
for(i=1;i<=n;i++) //normalization begins.
{
norm=a[i][i];
for(j=i;j<=2*n;j++)
{
a[i][j]=a[i][j]/norm;
}
}
printf("the inverse of matrix is\n");
for(i=1;i<=n;i++)
{
printf("\n");
for(j=n+1;j<=2*n;j++)
{
printf("%0.3f\t",a[i][j]);
}
}
getch();
}
Power Method in Numerical Method
#include<stdio.h>#include<conio.h>
#include<math.h>
void main()
{
float x[10][10],y[10],k[10],max1=0,max,diff;
int i,j,n;
clrscr();
printf("enter the order of Matrix:");
scanf("%d",&n);
for(i=0;i<n;i++)
{
for(j=0;j<n;j++)
{
printf("\nEnter coefficient matrix X:");
scanf("%f",&x[i][j]);
}
}
for(i=0;i<n;i++)
{
printf("Enter Matrix K:");
scanf("%f",&k[i]);
}
do
{
for (i=0;i<n;i++)
{
y[i]=0;
for(j=0;j<n;j++)
{
y[i]+=x[i][j]*k[j];
}
}
max=y[0];
for(i=0;i<n;i++)
{
if(fabs(y[i])>fabs(max));
max=y[i];
}
diff=max-max1;
for(i=0;i<n;i++)
{
y[i]=(1/max)*y[i];
k[i]=y[i];
}
max1=max;
}while(diff!=0);
printf("eigen value = %f", max);
printf("Eigen vector:");
for(i=0;i<n;i++)
printf("\n%f",k[i]);
getch();
}
RK 4th Order in Numerical Method
#include<stdio.h>#include<conio.h>
#include<math.h>
#define f(a,b) a*a+b*b
void main()
{
float x,x0,y,y0,h,x1,m,m1,m2,m3,m4;
int n,i;
clrscr();
printf("input initial values of x & y");
scanf("%f%f",&x,&y);
printf("enter step size");
scanf("%f",&h);
printf("enter the value of x at which y is to be calculated");
scanf("%f",&x1);
n=(x1-x)/h;
for(i=0; i<n; i++)
{
y0=y;
x0=x;
m1=f(x,y);
x=x0+h/2;
y=y0+(h*m1)/2;
m2=f(x,y);
x=x0+h/2;
y=y0+(h*m2)/2;
m3=f(x,y);
y=y0+h*m3;
x=x0+h;
m4=f(x,y);
m=(m1+2*m2+2*m3+m4)/6;
y=y0+h*m;
}
printf("the requored value is %f", y);
getch();
}
RK 4 for 2nd order Differential Equation in Numerical Method
#include<stdio.h>#include<conio.h>
#include<math.h>
float ODE(float);
float ODE1(float,float,float);
void main()
{
float x0,y10,y20,y1,xf,h,m11,m21,m31,m41,m12,m22,m32,m42,y11;
int i,n;
clrscr();
printf("Initial Guesses:\t");
scanf("%f %f %f",&x0,&y10,&y20);
printf("Enter Interval Size:\t");
scanf("%f",&h);
printf("Enter Final Point:\t");
scanf("%f",&xf);
n=(xf-x0)/h;
for (i=0;i<n;i++)
{
m11=ODE(y20);
m12=ODE1(x0,y10,y20);
m21=ODE(y20+m12*h/2);
m22=ODE1(x0+h/2,y10+m11*h/2,y20+m12*h/2);
m31=ODE(y20+m22*h/2);
m32=ODE1(x0+h/2,y10+m21*h/2,y20+m22*h/2);
m41=ODE(y20+m32*h);
m42=ODE1(x0+h,y10+m31*h,y20+m32*h);
y1=y10+((m11+2*m21+2*m31+m41)*h)/6;
y11=y20+((m12+2*m22+2*m32+m42)*h)/6;
y10=y1;
y20=y11;
x0=x0+h;
printf("\n\n%f, %f",y10,y20);
}
getch();
}
float ODE(float x)
{
float fxn=x;
return fxn;
}
float ODE1(float x, float y,float z)
{
float fxn=pow(x,2)*z+2*x*y;
return fxn;
}
Heun’s method in Numerical Method
#include<stdio.h>#include<conio.h>
#include<math.h>
#define f(x,y) 2*y/x
void main()
{
float x0,y1,y0,xf,h,m1,m2,m;
int i,n;
clrscr();
printf("Initial Guesses:\t");
scanf("%f %f",&x0,&y0);
printf("Enter Interval Size:\t");
scanf("%f",&h);
printf("Enter Final Point:\t");
scanf("%f",&xf);
n=(xf-x0)/h;
for (i=0;i<n;i++)
{
m1=f(x0,y0);
y1=y0+h*m1;
x0=x0+h;
m2=f(x0,y1);
m=(m1+m2)/2;
y0=y0+m*h;
}
printf("\n\ny(%f)= %f",x0,y0);
getch();
}
Heun’s Method for 2nd order differential Equation in Numerical Method
#include<stdio.h>#include<conio.h>
#include<math.h>
#define f1(y2) y2
#define f2(x,y1,y2) -y1-2*y2
void main()
{
float x0,y10,y20,y1,xf,h,m11,m12,m21,m22,m1,m2,y11;
int i,n;
clrscr();
printf("Initial Guesses:\t");
scanf("%f %f %f",&x0,&y10,&y20);
printf("Enter Interval Size:\t");
scanf("%f",&h);
printf("Enter Final Point:\t");
scanf("%f",&xf);
n=(xf-x0)/h;
for (i=0;i<n;i++)
{
m11=f1(y20);
m12=f2(x0,y10,y20);
m21=f1(y20+m12*h);
m22=f2((x0+h),(y10+m11*h),(y20+m12*h));
m1=(m11+m21)/2;
m2=(m12+m22)/2;
y1=y10+m1*h;
y11=y20+m2*h;
y10=y1;
y20=y11;
x0=x0+h;
printf("\n\n%f, %f",y10,y20);
}
getch();
}
Cubic Spline in Numerical Method
#include <stdio.h>#include <conio.h>
#include <math.h>
main()
{
float A[10],co[10][10],a[10],b[10],c[10],d[10],p,s[10],h[10],f[10],x[10],y;
int n,i,j,k,l;
clrscr();
x[1]=-0.6; f[1]=-0.5646;
x[2]=-0.3; f[2]=-0.2955;
x[3]=0.0; f[3]=0.0;
x[4]=0.3; f[4]=0.2955;
x[5]=0.6; f[5]=0.5646;
x[6]=0.9; f[6]=0.7833;
x[7]=1.2; f[7]=0.932;
n=7;
//find and print h[i]
for(i=1;i<=n-1;i++)
{
h[i]=x[i+1]-x[i];
printf("h[%d]=%f\n",i,h[i]);
}
getch();
c[1]=0;c[n]=0;
//calculate a[]
for(i=2;i<=n-1;i++)
A[i-1]=3*((f[i+1]-f[i])/h[i])-3*((f[i]-f[i-1])/h[i-1]);
//calculate co[][]
for(i=0;i<=n;i++)
{
for(j=0;j<=n;j++)
co[i][j]=0;
}
for(i=2;i<=n-1;i++)
{
co[i][i]=h[i-1];
co[i][i+1]=2*(h[i-1]+h[i]);
co[i][i+2]=h[i];
}
for(i=2;i<=n-1;i++)
co[i][n+1]=A[i-1];
//print co[][]
for(i=2;i<=n-1;i++)
{
for(j=3;j<=n+1;j++)
printf("%.3f\t",co[i][j]);
printf("\n");
} getch();
//find c[] by gauss elimination
for(k=2;k<=n-2;k++)
{
for(i=k+1;i<=n-1;i++)
{
if(fabs(co[i][k+1])<fabs(co[i+1][k+1]))
{
for(l=i;l<=n-1;l++)
{
for(j=k+2;j<=n+1;j++)
{
y=co[l][j];
co[l][j]=co[l+1][j];
co[l+1][j]=y;
}
}
}
}
for(i=k+1;i<=n-1;i++)
{
for(j=k+2;j<=n+1;j++)
{
co[i][j]=co[i][j]-(co[i][k]/co[k][k])*co[k][j];
}
}
}
c[n-1]=co[n-1][n]/co[n-1][n-1];
for(k=n-2;k>=2;k--)
{
y=0;
for(j=k+2;j<=n;j++)
{
y+=co[k][j]*c[j];
}
c[k]=(1/co[k][k])*(co[k][n]-y);
}
//print c[] and calculate a[],b[],d[]
for(i=1;i<=n-1;i++)
{
printf("\nc[%d]=%f",i,c[i]);
a[i]=f[i];
b[i]=((f[i+1]-f[i])/h[i])-h[i]*(2*c[i]+c[i+1])/3;
d[i]=(c[i+1]-c[i])/(3*h[i]);
}
k=1;
do{
a:printf("\n enter the value of x:");
scanf("%f",&p);
if((p<x[1])||(p>x[n])){
printf("\nnot in range");
goto a;
}
for(j=1;j<=n;j++){
if((p>=x[j])&&(p<x[j+1]))
i=j;
}
y=p-x[i];
s[i]=a[i]+b[i]*y+c[i]*pow(y,2)+d[i]*pow(y,3);
printf("the interpolated value is:%f",s[i]);
k++;
}
while(k<=3);
getch();
return 0;
}
Shooting Method in Numerical Method
#include<stdio.h>#include<conio.h>
#include<math.h>
#define f1(x,y,z) z
#define f2(x,y,z) 3*x+4*y
void main()
{
float s11,s12,s21,s22,x0,y0,z0,y[5],z[5],xn,b[5],yn,m[5],x,y1,h;
int i,j,k,n;
clrscr();
printf("enter the range of x:");
scanf("%f%f",&x0,&xn);
printf("enter the interval width h:");
scanf("%f",&h);
printf("enter the initial condition:\n");
printf("x0=");
scanf("%f",&x0);
x=x0;
printf("y0=");
scanf("%f",&y0);
y1=y0;
printf("yn=");
scanf("%f",&yn);
n=(xn-x0)/h;
for(j=0;j<=2;j++)
{
printf("iteration %d:",j+1);
if(j<2)
{
printf("enter the initial guess for z0=");
scanf("%f",&m[j]);
}
z0=m[j];
for(i=0;i<n;i++)
{ //if the given differential eqation is:(d2y/dx2)=3x+4y
s11=f1(x0,y0,z0); //let dy/dx=z => f1
s21=f2(x0,y0,z0); //then (d2y/dx2)=3x+4y =>f2
s12=f1((x0+h),(y0+h*s11),(z0+h*s21));
s22=f2((x0+h),(y0+h*s11),(z0+h*s21));
y[i]=y0+h*(s11+s12)/2;
z[i]=z0+h*(s21+s22)/2;
x0=x0+h;
y0=y[i];
z0=z[i];
}
b[j]=y0;//store the value of y at x=1 at the array b[]
x0=x; //initialize x0 back to the initial value
y0=y1; //initialize y0 back to the initial value
if(j==1) //find the better approximate for the initial value of z
{
m[2]=m[1]-((b[1]-yn)*(m[1]-m[0])/(b[1]-b[0]));
}
for(i=0;i<n;i++)
printf("y at x=%f is %f\n",(x+h*(i+1)),y[i]);
}
getch();
}
Programming of numerical method
c programming of numerical method
Post a Comment