並行求解三對角矩陣 CR方法

此函數求解三對角矩陣，傳入數組是從1開始，不是從0開始

此方法將求解分爲兩個部分，理論上能夠減小一半的計算時間

編譯器要支持OpenMP

int pTM2(double a[], double b[], double c[], double r[],  double x[],int n){
	//原矩陣形式 要求 對角佔優 c(1)=0 b(n)=0
	//b1 c1
	//a2 b2 c2
	//   a3 b3 c3
	//    ···
	//         an bn

	double *f, *g,*h;
	f=new double[n+1];
	h=new double[n+1];
	g=new double[n+1];
	int p=n/2;
#pragma omp parallel sections
{
#pragma omp section
{
	g[1]=c[1]/b[1];
	h[1]=r[1]/b[1];
	for(int i=2;i<=p;i++){
		g[i]=c[i]/(b[i]-a[i]*g[i-1]);
		h[i]=(r[i]-a[i]*h[i-1])/(b[i]-a[i]*g[i-1]);
	}
}
#pragma omp section
{
	f[n]=a[n]/b[n];
	h[n]=r[n]/b[n];
	for(int j=n-1;j>=(p+1);j--){
		f[j]=a[j]/(b[j]-c[j]*f[j+1]);
		h[j]=(r[j]-c[j]*h[j+1])/(b[j]-c[j]*f[j+1]);
	}
}
}
#pragma omp parallel sections
{
#pragma omp section
{
	x[p]=(h[p]-g[p]*h[p+1])/(1-g[p]*f[p+1]);
	for(int i=p-1;i>=1;i--){
		x[i]=h[i]-g[i]*x[i+1];
	}
}
#pragma omp section
{
	x[p+1]=(h[p+1]-f[p+1]*h[p])/(1-g[p]*f[p+1]);
	for(int j=p+2;j<=n;j++){
		x[j]=h[j]-f[j]*x[j-1];
	}
}
}
	delete[] f;
	delete[] g;
	delete[] h;
	return 0;
}