codice:
#include <stdio.h>
#include <malloc.h>
main ()
{
double **coeff, *x, *y, *xval, *yval, *lambda;
int n, nn, i, j, sizef, flag , txt;
FILE *nodi,*punti;
/* Scegli come leggere la dimensione e i nodi*/
txt=0;
printf("Vuoi leggere la dimensione e i nodi da file(seleziona 1) o da tastiera?\n");
scanf("%d",&txt);
if (txt==1)
{
/* Lettura matrice da file */
//Associo l'identificativo al file
nodi=fopen("nodi.txt","r");
fscanf(nodi,"%d",&n);
//Allocazione dei vettori dei nodi
x=(double *)calloc(n,sizeof(double));
y=(double *)calloc(n,sizeof(double));
//Allocazione del vettore dei lambda
lambda=(double *)calloc(n,sizeof(double));
//Lettura dei punti da interpolare
for(i=0;i<n;i++)
{
fscanf(nodi,"%lf",&x[i]);
fscanf(nodi,"%lf",&y[i]);
}
}
else
{
/* Lettura numero di nodi */
printf("Inserire il numero di nodi\n");
scanf("%d",&n);
//Allocazione dei vettori dei nodi
x=(double *)calloc(n,sizeof(double));
y=(double *)calloc(n,sizeof(double));
//Allocazione del vettore dei lambda
lambda=(double *)calloc(n,sizeof(double));
//Lettura dei nodi da tastiera
for(i=0;i<n;i++)
{
printf("Inserisci l'ascissa del nodo %d-smo \n",i);
scanf("%lf" , &x[i]);
printf("Inserisci l'ordinata del nodo %d-smo \n",i);
scanf("%lf" , &y[i]);
}
}
//Allocazione della matrice dei coefficienti
coeff=(double **)calloc(n-1,sizeof(double *));
sizef=sizeof(double);
for(i=0;i<n;i++)
{
coeff[i]=(double *)calloc(4,sizef);
}
i=0;
do
{
j=i+1;
do
{
if(x[i]==x[j])
flag=0;
j=j+1;
}
while((j<n)&&(flag));
i=i+1;
}
while((i<n-1)&&(flag));
if(flag!=0)
{
//ordina i nodi con la funzione creata
ordina(n,x,y);
flag=coefficienti(n, x, y, coeff,lambda);//calcola i coeff della spline
for(i=0;i<n-1;i++)
{
for(j=0;j<4;j++)
{
printf("La matrice dei coefficienti della spline ha come elemento nelle posizione %d %d coeff=%f\n",i,j,coeff[i][j]);
}
}
/* Scegli come leggere i punti in cui valutare la spline*/
txt=0;
printf("Vuoi leggere i punti in cui valutare la spline da file(seleziona 1) o da tastiera?\n");
scanf("%d",&txt);
/* Lettura numero di pt in cui valutare la spline */
printf("Inserire il numero di punti in cui valutare la spline\n");
scanf("%d",&nn);
if (txt==1)
{
//apre il file dei punti in cui valutare la spline
punti=fopen("punti.txt","r");
xval=(double *)calloc(nn,sizeof(double));
yval=(double *)calloc(nn,sizeof(double));
for(i=0;i<nn;i++)
{
fscanf(punti,"%lf",&xval[i]);
}
}
else
{
/* Lettura numero di pt in cui valutare la spline */
printf("Inserire il numero di punti in cui valutare la spline\n");
scanf("%d",&nn);
//Allocazione dei vettori contenenti ascisa ed ordinata dei punti di valutazione
xval=(double *)calloc(nn,sizeof(double));
yval=(double *)calloc(nn,sizeof(double));
//Lettura dei punti di valutazione da tastiera
for(i=0;i<nn;i++)
{
printf("Inserisci l'ascissa del %d-smo punto di valutazione\n",i);
scanf("%lf" , &xval[i]);
}
}
flag=val(n,x,y[0],y[n-1],coeff,lambda[n-1],nn,xval,yval,lambda,y);//valuta la spline nei punti dati
}
else
{
printf("Due nodi sono uguali\n");
}
}
//calcola i coefficienti della spline
int coefficienti(int n, double *x, double *y,double **coeff,double *lambda)
{
int i,flag;
double bs[n-1],bi[n-1], diffdiv[n-1], dp[n],h[n-1];
/* calcolo delle ampiezze degli intervalli */
for (i=0;i<n;i++)
{
h[i]=x[i+1]-x[i];
}
// costruzione delle diagonali superiore ed inferiore della matrice B
bs[0]=1;
for(i=1;i<n-1;i++)
{
bs[i]=h[i-1]/(x[i+1]-x[i-1]); //=delta[i]
bi[i-1]=1-bs[i];
}
bi[n-2]=1;
/* calcolo delle differenze divise del primo ordine */
for(i=0;i<n;i++)
{
diffdiv[i]=(y[i+1]-y[i])/h[i];
}
flag=sistema(n, bs, bi, diffdiv, lambda);
//Calcolo dei coefficienti
for(i=0;i<n-1;i++)
{
coeff[i][0]=y[i];
coeff[i][1]=lambda[i];
coeff[i][2]=(diffdiv[i]-lambda[i])/h[i];
coeff[i][3]=(-1)*(coeff[i][2]/h[i])-((diffdiv[i]-lambda[i+1])/h[i]*h[i]);
}
return 1;
}
//risolve il sistema tridiagonle
int sistema(int n, double *d, double *e, double *diffdiv, double *lambda)
{
int flag,i;
double y[n],m[n-1],dp[n];
i=1;
dp[0]=2;
do
{
m[i-1]=e[i-1]/dp[i-1];
dp[i]=2-m[i-1]*d[i-1];
flag=dp[i];
i++;
}
while((i<n)&&(flag!=0));
//risolvo il sistema 3Bdeltaf
y[0]=3*diffdiv[0];
for(i=1;i<n-1;i++)
{
y[i]=3*(e[i-1]*diffdiv[i-1]+d[i]*diffdiv[i])-(y[i-1]*m[i-1]);
}
y[n-1]=3*diffdiv[n-2]-y[n-2]*m[n-2];
//Calcolo i lambda
lambda[n-1]=y[n-1]/dp[n-1];
for(i=2;i<=n;i++)
{
lambda[n-i]=y[n-i]-d[n-i]*lambda[n-i+1];
lambda[n-i]=lambda[n-i]/dp[n-i];
}
return flag;
}
//ricerca l'intervallo a cui appartiene il punto
int ricint(int n, double *x, double xval)
{
int i,j,m,flag;
;
if(xval<x[0])
{
flag=0;
}
else if(xval>x[n-1])
{
flag=n;
}
else
{
i=0;
j=n-1;
do
{
m=(i+j)/2;//pt medio dell'intervallo di definizione
if(xval<x[m])
{
j=m;
}
else
{
i=m;
}
}
while((j-i)!=1);
flag=j;
}// quì termina l'else
return flag;
}
//funzione che valuta la spline in un vettore di punti
int val(int n,double *x,double y1,double yn,double **coeff,double l,int nn,double *xval,double *yval, double *lambda, double *y)
{
int i,j, intapp;
FILE *f_mat;
for(i=0;i<nn;i++)
{
intapp=ricint(n,x,xval[i]);
if(xval[i]<x[0])
{
yval[i]=y[0]+lambda[0]*(xval[i]-x[0]);
}
else if(xval[i]>x[n-1])
{
yval[i]=y[n-1]+lambda[n-1]*(xval[i]-x[n-1]);
}
else
{
yval[i]=coeff[intapp-1][3]*(xval[i]-x[intapp])+coeff[intapp-1][2];
for(j=0;j<2;j++)
{
yval[i]=yval[i]*(xval[i]-x[intapp-1])+coeff[intapp-1][1-j];
}
printf("La spline nel punto di valutazione %d-smo, cioè in corrispondenza di xval=%f vale yval=%f\n",i,xval[i],yval[i]);
}
}
//salvo sul file valori.txt
f_mat=fopen("valori.txt","w");
for(i=0;i<nn;i++)
fprintf( f_mat, "%lf\n",yval[i]);
return 1;
}
//funzione che uso per ordinare i nodi di interpolazione e le corrispondenti ordinate
int ordina(int n, double *x, double *y)
{
int i,stop,flag;
double app;
flag=1;
stop=n-1;
do
{
flag=0;
for(i=0;i<stop;i++)
{
if(x[i]>x[i+1])
{
/* scambio dei nodi consecutivi */
app=x[i];
x[i]=x[i+1];
x[i+1]=app;
/* scambio delle ordinate corrispondenti ai nodi scambiati */
app=y[i];
y[i]=y[i+1];
y[i+1]=app;
flag=1;
}
}
stop=stop-1;
}
while(flag!=0);
return flag;
}