C++算法

C++算法

1.用MEX接口

这种方法是最常见的方法，相对来说也比较简单，主要涉及修改C++中入口函数和传递参数格式。

首先在matlab中输入

mex -setup

用matlab 调用C++
运行结果

连分式主代码

//连分式求积法.cpp
#include
#include
using namespace std;
//计算函数连分式值
double funpqv(double x[], double b[], int n, double t)
{
int k;
double u;
u = b[n];
for (k = n - 1; k >= 0; k–)
{
if (fabs(u) + 1.0 == 1.0)
u = 1.0e1 + 35 * (t - x[k]) / fabs(t - x[k]);
else
u = b[k] + (t - x[k]) / u;
}
return(u);

}
//计算连分式新的一节b[j]

void funpaj(double x[], double y[], double b[], int j)
{
int k,flag = 0;
double u;
u = y[j];
for (k = 0; (k < j) && (flag == 0); k++)
{
if ((u - b[k]) + 1.0 == 1.0)flag = 1;
else
u = (x[j] - x[k]) / (u - b[k]);
}
if (flag == 1) u = 1e1 + 35;
b[j] = u;
return;

}
double pqinteg(double a0, double b0, double eps, double (f) (double))
{
int k, j, il, flag, n;
double h, * g, * b, h0, g0, d, s0, s1, x;
b = new double[10];
h = new double[10];
g = new double[10];
il = 0; n = 1; h0 = (b0 - a0) / n; flag = 0;
g0 = h0 * ((*f) (a0) + (*f)(b0)) / 2.0;
while ((il < 20 && (flag == 0)))
{
il = il + 1;
h[0] = h0; g[0] = g0; b[0] = g[0];
j = 1; s1 = g[0];
while (j <= 7)
{
d = 0.0;
for (k = 0; k <= n - 1; k++)
{
x = a0 + (k + 0.5) * h[j - 1];
d = d + (*f)(x);

		}
		g[j] = (g[j - 1] + h[j - 1] * d) / 2.0;
		h[j] = h[j - 1] / 2.0;
		n = 2 * n;
		funpaj(h, g, b, j);
		s0 = s1;
		s1 = funpqv(h, b, j, 0.0);
		if (fabs(s1 - s0) >= eps)j = j + 1;
		else
		{
			cout << "最后一次迭代连分式节数=" << j << endl;
			j = 10;
		}
	}

	h0 = h[j - 1]; g0 = g[j - 1];
	if (j == 10) flag = 1;
}
cout << "迭代次数=" << il << endl;
delete[] b; delete[] h; delete[] g;
return(s1);

}

#include"mex.h"

//例
#include
#include
#include “连分式求积分.cpp”
using namespace std;
int main()
{
double d, eps, pqintegf(double);
eps = 0.0000001;
d = pqinteg(0.0, 4.3, eps, pqintegf);
cout << “积分值s=” << d << endl;
return 0;

}
double pqintegf(double x)
{
return(exp(-x * x));
}

高斯求积分

调用

#include
#include
using namespace std;
void part(double a, double b, int m，int n, double fa[], double fb[], double s[2])
{
int mm, k, j;
double sa[4], sb[4], ca[4], cb[4], sma, smb, cma, cmb;
sma - sin(m * a); smb = sin(m * b);
cma = cos(m * a); cmb = cos(m * b);
sa[0] = sma; sa[1] = cma; sa[2] = -sma; sa[3] = -cma;
sb[0] = smb; sb[1] = cmb; sb[2] = -smb; sb[3] = -cmb;
ca[0] = cma; ca[1] = -sma; ca[2] = -cma; ca[3] = sma;
cb[0] = cmb; cb[1] = -smb; cb[2] = -cmb; cb[3] = smb;
s[0] = 0.0; s[1] = 0.0;
mm = 1;
for (k = 0; k <= n - 1; k＋＋)
{
j = k;
while (>= 4)j = j - 4;
mm = mm * m;
s[0] = s[0] + (fb[k] * sb[j] - fa[k] * sa[j]) / (1.0 * mm);
s[1] = s[1] + (fb[k] * cb[j] - fa[k] * ca[j]) / (1.0 * mm);
}
s[1]–s[1];
return.

}

#include"mex.h"

例

include

#include
#include"高振荡函数求积法.cpp"
using namespace std;
int main()
{
int n, m;
double a, b;
double s[2], fa[4] = { 0.0,1.0,0.0,-3.0 };
double fb[4] = { 6.2831852,1.0,-6.2831852,-3.0 };
a = 0.0; b = 6.2831852;
m = 30; n = 4;
part(a, b, m, n, fa, fb, s);
cout << “s(0)=” << s[0] << endl;
cout << “s(1)=” << s[1] << endl;
return 0;

}

结果