多项式与连分式函数的计算

最新推荐文章于 2022-10-03 15:46:11 发布

转载最新推荐文章于 2022-10-03 15:46:11 发布 · 689 阅读

1 ·

CC 4.0 BY-SA版权

原文链接：https://my.oschina.net/RapidBird/blog/3424

文章标签：

#python #c/c++

2019独角兽企业重金招聘Python工程师标准>>>

#include "stdafx.h"
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include "ComplexNumber.h"
#include "Multinomial.h"

//复系数多项式求值
void ncply(double ar[], double ai[], int n, double x, double y, double *u, double *v)
{
int i;
double p,q,s,t;
s=ar[n-1];
t=ai[n-1];
for (i=n-2; i>=0; i--)
{
  ocmul(s,t,x,y,&p,&q);
  s=p+ar[i];
  t=q+ai[i];
}
*u=s;
*v=t;
return;
}

//一维多项式多组求值
void nplys(double a[], int n, double x[], int m, double p[])
{
int i,j,mm,nn,ll,t,s,kk,k;
double *b,y,z;
b=(double *)malloc(2*n*sizeof(double));
y=a[n-1];
for (i=0; i<=n-1; i++) b[i]=a[i]/y;
k=log(n-0.5)/log(2.0)+1;
nn=1;
for (i=0; i<=k-1; i++) nn=2*nn;
for (i=n; i<nn-1; i++) b[i]=0.0;
b[nn-1]=1.0;
t=nn;
s=1;
for (i=1; i<=k-1; i++)
{
  t=t/2;
  mm=-t;
  for (j=1; j<=s; j++)
  {
   mm=mm+t+t;
   b[mm-1]=b[mm-1]-1.0;
   for (kk=2; kk<=t; kk++)
    b[mm-kk]=b[mm-kk]-b[mm-1]*b[mm+t-kk];
  }
  s=s+s;
}
for (kk=1; kk<=m; kk++)
{
  for (i=0; i<=(nn-2)/2; i++)
   a[i]=x[kk-1]+b[2*i];
  mm=1;
  z=x[kk-1];
  for (i=1; i<=k-1; i++)
  {
   mm=mm+mm;
   ll=mm+mm;
   z=z*z;
   for (j=0; j<=nn-1; j=j+ll)
    a[j/2]=a[j/2]+a[(j+mm)/2]*(z+b[j+mm-1]);
  }
  z=z*z/x[kk-1];
  if (nn!=n) a[0]=a[0]-z;
  p[kk-1]=a[0]*y;
}
return;
}

//复系数多项式相乘
void ncmul(double pr[], double pi[], int m, double qr[], double qi[], int n, double sr[], double si[], int k)
{
int i,j;
double a,b,c,d,u,v;
for (i=0; i<=k-1; i++)
{
  sr[i]=0.0;
  si[i]=0.0;
}
for (i=0; i<=m-1; i++)
  for (j=0; j<=n-1; j++)
  {
   a=pr[i];
   b=pi[i];
   c=qr[j];
   d=qi[j];
   ocmul(a,b,c,d,&u,&v);
   sr[i+j]=sr[i+j]+u;
   si[i+j]=si[i+j]+v;
  }
return;
}

//多项式相乘
void npmul(double p[], int m, double q[], int n, double s[], int k)
{
int i,j;
for (i=0; i<=k-1; i++) s[i]=0.0;
for (i=0; i<=m-1; i++)
for (j=0; j<=n-1; j++)
s[i+j]=s[i+j]+p[i]*q[j];
return;
}

//多项式相除
void npdiv(double p[], int m, double q[], int n, double s[], int k, double r[], int l)
{
int i,j,mm,ll;
for (i=0; i<=k-1; i++) s[i]=0.0;
if (q[n-1]+1.0==1.0) return;
ll=m-1;
for (i=k; i>=1; i--)
{
  s[i-1]=p[ll]/q[n-1];
  mm=ll;
  for (j=1; j<=n-1; j++)
  {
   p[mm-1]=p[mm-1]-s[i-1]*q[n-j-1];
   mm=mm-1;
  }
  ll=ll-1;
}
for (i=0; i<=l-1; i++) r[i]=p[i];
return;
}

//复系数多项式相除
void ncdiv(double pr[], double pi[], int m, double qr[], double qi[], int n, double sr[], double si[], int k, \
     double rr[], double ri[], int l)
{
int i,j,mm,ll;
double a,b,c,d,u,v;
for (i=0; i<=k-1; i++)
{
  sr[i]=0.0;
  si[i]=0.0;
}
d=qr[n-1]*qr[n-1]+qi[n-1]*qi[n-1];
if (d+1.0==1.0) return;
ll=m-1;
for (i=k; i>=1; i--)
{
  a=pr[ll];
  b=pi[ll];
  c=qr[n-1];
  d=qi[n-1];
  ocdiv(a,b,c,d,&u,&v);
  sr[i-1]=u;
  si[i-1]=v;
  mm=ll;
  for (j=1; j<=n-1; j++)
  {
   a=sr[i-1];
   b=si[i-1];
   c=qr[n-j-1];
   d=qi[n-j-1];
   ocmul(a,b,c,d,&u,&v);
   pr[mm-1]=pr[mm-1]-u;
   pi[mm-1]=pi[mm-1]-v;
   mm=mm-1;
  }
  ll=ll-1;
}
for (i=0; i<=l-1; i++)
{
  rr[i]=pr[i];
  ri[i]=pi[i];
}
return;
}

//一维多项式求值
double nplyv(double a[], int n, double x)
{
int i;
double u;
u=a[n-1];
for (i=n-2; i>=0; i--)
u=u*x+a[i];
return(u);
}

//函数连分式的计算
double nfpqv(double x[], double b[], int n, double t)
{
int k;
double u;
u=b[n-1];
for (k=n-2; k>=0; k--)
{
  if (fabs(u)+1.0==1.0)
   u=1.0e+35*(t-x[k])/fabs(t-x[k]);
  else
   u=b[k]+(t-x[k])/u;
}
return(u);
}

//二维多项式求值
double nbply(double a[], int m, int n, double x, double y)
{
int i,j;
double u,s,xx;
u=0.0;
xx=1.0;
for (i=0; i<=m-1; i++)
{
  s=a[i*n+n-1]*xx;
  for (j=n-2; j>=0; j--)
   s=s*y+a[i*n+j]*xx;
  u=u+s;
  xx=xx*x;
}
return(u);

}

----根据《C语言常用算法程序集》整理

投票>

转载于:https://my.oschina.net/RapidBird/blog/3424