我这里有一段程序关于LU分解的，贴出来，大家给点意见，感激~

[caijinqing]

vv= new double[rows];//vv stores the implicit scaling of each row
indx= new int[rows];//indx is an output vector which records the row permutation effected by the partial pivoting
int rmax=0;
double TINY=1.0E-305;

for(int r=0;r<rows;r++)// find out the max value per row
{
double Max=0.0;

for(int c=0;c<cols;c++)
{
double temp=fabs(Get(r,c));
if(temp>Max)

Max=temp;

}
if(Max<TINY)// if the max value less than the tiny, it's a singular matrix
{
cout<<"Singular matrix"<<endl;
exit(0);
}
vv[r]=1.0/Max;// row r's scale factor
}

for (int c=0;c<cols;c++)// for r<c: find out the elements of U, and put them in the posation of the former matrix
{
if(c>0)
{
for(int r=0;r<c;r++)
{
double sum=Get(r,c);
if(r>0)
{
for(int k=0;k<r;k++)sum=sum-Get(r,k)*Get(k,c);
Set(r,c,sum);
}
}
}

double Max=0.0;
for(int r=c;r<rows;r++)//for r=c r<rows: find out the L,and put them in the posation of the former matrix
{

double sum=Get(r,c);
if(c>0)
{
for(int k=0;k<c;k++) sum=sum- Get(r,k)*Get(k,c);
Set(r,c,sum);
}
double dum=vv[r]*fabs(sum);
if(dum>= Max)
{
Max=dum;
rmax=r;
}


}
if(Max<TINY)
{
cout<<"Singular matrix"<<endl;
exit(0);
}
if(c!=rmax)// change rows, put the max value in the right persation
{
for (int k=0;k<cols;k++)
{
double a,b;

a=Get(rmax,k);
b=Get(c,k);

Set(rmax,k,b);
Set(c,k,a);
}
double dum= vv[c];
vv[c]=vv[rmax];
vv[rmax]=dum;

/*vv[rmax]=vv[c];*/ //also change the scale factor
}

indx[c]=rmax;
if( c!=cols-1)//divide by the pivot element
{
double a=Get(c,c);
if(a==0) a=TINY;
double dum=1.0/Get(c,c);
for(int r=c+1; r<rows;r++)

Set(r,c,Get(r,c)*dum);

}
}

[会飞的老鱼]

我这有一个C语言的

/*{ ************************************************************** }*/
/*{ 下面函数使用按列选主元作出矩阵A的三角分解(ＬＲ分解) }*/
/*{ ************************************************************** }*/

int LR_Method01(int n,double A[][MAX_SIZE])
{ int i,j,k,s,row;
double sum,max,t;
for(k=1;k<n;k++)
{ max=0; row=k;
for(i=k;i<=n;i++) /* 选列主元 */
{ sum=A[i][k];
for(s=1;s<k;s++) sum=sum-A[i][s]*A[s][k];
if (max<fabs(sum))
{ max=fabs(sum); row=i; }
}
if (row!=k) /*若列主元所在的行row不等于k则交换矩阵Ａ的第k行和第row行*/
{ for(s=1;s<=n;s++) /* 以及常数项b(k)和b(row) */
{ t=A[k][s]; A[k][s]=A[row][s]; A[row][s]=t; }
}
/* 计算r(k,k) ,k=1,2,…,n-1 */
for(j=1;j<k;j++) A[k][k]=A[k][k]-A[k][j]*A[j][k];
for(i=k+1;i<=n;i++)
{ for(s=1;s<k;s++)
{ A[i][k]=A[i][k]-A[i][s]*A[s][k]; /* 计算l(i,k) ,i=k+1,k+2,…,n */
A[k][i]=A[k][i]-A[k][s]*A[s][i]; /* 计算r(k,i) ,i=k+1,k+2,…,n */
}
A[i][k]=A[i][k]/A[k][k];
}
}
for(s=1;s<n;s++) A[n][n]=A[n][n]-A[n][s]*A[s][n]; /* 计算r(n,n) */
return(1);
} /* 函数结束 */

/*{ 参数说明： }
/*{ n 为系数矩阵Ａ的行数(或列数)，即Ａ的大小，二维数组Ａ用来 }*/
/*{ 存放系数矩阵Ａ的各元素，MATRIX为存放系数矩阵Ａ的二维数组 }*/
/*{ 的类型，其定义如下： }*/
/*{ #define MAX_SIZE 50 }*/
/*{ double A[MAX_SIZE][MAX_SIZE]; }*/
/*{ 注：上面过程利用系数矩阵Ａ的上三角元素位置来存放矩阵Ｒ的上三角}*/
/*{ 部分的元素，其主对角线以下的零元不存储，利用Ａ的主对角线以}*/
/*{ 下元素的位置来存放矩阵Ｌ的下三角部分的元素，而其主对角线上}*/
/*{ 的1和上三角部分的零元不存储。 }*/

[gbb21]

[Quote=引用楼主 caijinqing 的回复:]
vv= new double[rows];//vv stores the implicit scaling of each row
    indx= new int[rows];//indx is an output vector which records the row permutation effected by the partial pivoting
int rmax=0;
double TINY=1.0E-305;

for(int r=0;r <rows;r++)// find out the max value per row
{
double Max=0.0;

for(int c=0;c <cols;c++)
{
double temp=fabs(Get(r,c));
if(temp>Max)

Max=temp;

}
if(Max <TINY)// if the max value less than the tiny, it's a singular matrix
{
cout < <"Singular matrix" < <endl;
exit(0);
}
vv[r]=1.0/Max;// row r's scale factor
}

for (int c=0;c <cols;c++)// for r <c: find out the elements of U, and put them in the posation of the former matrix
{
if(c>0)
{
for(int r=0;r <c;r++)
{
double sum=Get(r,c);
if(r>0)
{
for(int k=0;k <r;k++)sum=sum-Get(r,k)*Get(k,c);
Set(r,c,sum);
}
}
}

double Max=0.0;
for(int r=c;r <rows;r++)//for r=c r <rows: find out the L,and put them in the posation of the former matrix
{

double sum=Get(r,c);
if(c>0)
{
for(int k=0;k <c;k++) sum=sum- Get(r,k)*Get(k,c);
Set(r,c,sum);
}
double dum=vv[r]*fabs(sum);
if(dum>= Max)
{
Max=dum;
rmax=r;
}


}
if(Max <TINY)
{
cout < <"Singular matrix" < <endl;
exit(0);
}
if(c!=rmax)// change rows, put the max value in the right persation
{
for (int k=0;k <cols;k++)
{
double a,b;

a=Get(rmax,k);
b=Get(c,k);

Set(rmax,k,b);
Set(c,k,a);
}
double dum= vv[c];
vv[c]=vv[rmax];
vv[rmax]=dum;

/*vv[rmax]=vv[c];*/    //also change  the scale factor
}

indx[c]=rmax;
if( c!=cols-1)//divide by the pivot element
{
double a=Get(c,c);
if(a==0) a=TINY;
double dum=1.0/Get(c,c);
for(int r=c+1; r <rows;r++)

  Set(r,c,Get(r,c)*dum);

}
}
[/Quote]
what kind of comments do you want?

[guoyuqing1988]

good

[c370141743]

好

[谷穗儿机构]

哈哈，看不明白

[new_universe]

good

[popxiha]

for (int c=0;c <cols;c++)// for r <

[glibfox]

学习了

[livem]

啊

[caijinqing]

抱歉，只贴了程序出来，没有提问题，我是希望大家能够给点意见，把这段程序并行。现在虽然并行出来了但是效率很低时间节省不到一半，下面这段程序应该可以并行 我原来使用的single dircetive, 请问大家这段程序如何并行呢
for (int c=0;c <cols;c++)// for r <c: find out the elements of U, and put them in the posation of the former matrix
{
if(c>0)
{
for(int r=0;r <c;r++)
{
double sum=Get(r,c);
if(r>0)
{
for(int k=0;k <r;k++)sum=sum-Get(r,k)*Get(k,c);
Set(r,c,sum);
}
}
}

[AAA20090987]

mark............