openmp+mpi混合编程实现Smith-Waterman算法

GG_wang 2014-06-02 02:05:11
此Smith-Waterman 算法分别用mpi与openmp实现是没问题的,但是两个混合编程的时候就会出各种问题,希望懂的能够给指条明路。。。万分感谢


#include <omp.h>
#include <time.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <mpi.h>
/*File pointers*/
FILE *ptr_file_1, *ptr_file_2,*ptr_file_3; //test file pointers

/*Definitions*/
#define TRUE 1
#define FALSE 0
#define Match 1
#define MissMatch -1
#define GapPenalty 2
#define MAXLEN 20

/*Global Variables*/
char inputC[5]; //user input character
int inputI;
int StrLen1,StrLen2;
int intcheck = TRUE;

char holder, ch;
int filelen1 = 0;
int filelen2 = 0;
int filelen3 = 0;
int i,j,k,l,m,n,lenA,lenB,compval,top,left,left_top,rows,columns,contrl;
char FASTA1[MAXLEN],FASTA2[MAXLEN];
char dash = '-';

char strA[MAXLEN]; //holds 1st string to be aligned in character array
char strB[MAXLEN]; //holds 2nd string to be aligned in character array
int HiScore; //holds value of highest scoring alignment(s).
int HiScorePos[2]; //holds the position of the HiScore
int SWArray[MAXLEN][MAXLEN]; //S-W Matrix

char MaxA[MAXLEN];
char MaxB[MAXLEN];
char OptA[MAXLEN];
char OptB[MAXLEN];

int MaxAcounter = 1; //MaxA counter
int MaxBcounter = 1; //MaxB counter
int cont = TRUE;
int check;


int rank,size;
MPI_Status status;
int buff[3];



/*ALIGNMENT FUNCTION*/
int Align(int PosA,int PosB) {

/*Function Variables*/
int relmax = -1; //hold highest value in sub columns and rows
int relmaxpos[2]; //holds position of relmax

if(SWArray[PosA-1][PosB-1] == 0) {
cont = FALSE;
}

while(cont == TRUE) { //until the diagonal of the current cell has a value of zero

/*Find relmax in sub columns and rows*/
for(i=PosA; i>0; --i) {

if(relmax < SWArray[i-1][PosB-1]) {

relmax = SWArray[i-1][PosB-1];
relmaxpos[0] = i-1;
relmaxpos[1] = PosB-1;
}
}

for(j=PosB; j>0; --j) {

if(relmax < SWArray[PosA-1][j-1]) {

relmax = SWArray[PosA-1][j-1];
relmaxpos[0]=PosA-1;
relmaxpos[1]=j-1;
}
}

/*Align strings to relmax*/
if((relmaxpos[0] == PosA-1) && (relmaxpos[1] == PosB-1)) { //if relmax position is diagonal from current position simply align and increment counters

MaxA[MaxAcounter] = strA[relmaxpos[0]-1];
++MaxAcounter;
MaxB[MaxBcounter] = strB[relmaxpos[1]-1];
++MaxBcounter;

}

else {

if((relmaxpos[1] == PosB-1) && (relmaxpos[0] != PosA-1)) { //maxB needs at least one '-'

for(i=PosA-1; i>relmaxpos[0]-1; --i) { //for all elements of strA between PosA and relmaxpos[0]

MaxA[MaxAcounter]= strA[i-1];
++MaxAcounter;
}
for(j=PosA-1; j>relmaxpos[0]; --j) { //set dashes to MaxB up to relmax

MaxB[MaxBcounter] = dash;
++MaxBcounter;
}

MaxB[MaxBcounter] = strB[relmaxpos[1]-1]; //at relmax set pertinent strB value to MaxB
++MaxBcounter;
}

if((relmaxpos[0] == PosA-1) && (relmaxpos[1] != PosB-1)) { //MaxA needs at least one '-'

for(j=PosB-1; j>relmaxpos[1]-1; --j) { //for all elements of strB between PosB and relmaxpos[1]

MaxB[MaxBcounter] = strB[j-1];
++MaxBcounter;
}
for(i=PosB-1; i>relmaxpos[1]; --i) { //set dashes to strA

MaxA[MaxAcounter] = dash;
++MaxAcounter;
}

MaxA[MaxAcounter] = strA[relmaxpos[0]-1];
++MaxAcounter;
}
}

//printf("(%i,%i)",relmaxpos[0],relmaxpos[1]);
Align(relmaxpos[0],relmaxpos[1]);
}

return(cont);
}


int SWMax(int num1,int num2,int num3)
{
int max=-1;
if(num1>num2)
{
max = num1;
}
else if(num2>num3)
{
max = num2;
}else
{
max =num3;
}
return max;
}

/*MAIN FUNCTIONS*/
int main(int argc,char ** argv)
{

MPI_Init(&argc,&argv);//MPI初始化
MPI_Comm_rank(MPI_COMM_WORLD,&rank);//得到当前进程标识
MPI_Comm_size(MPI_COMM_WORLD,&size);//总的进程个数

printf("size=%d\n",size);


/*open first file*/
ptr_file_1 = fopen("d:/mpi/str1.fa", "r");
/*check to see if it opened okay*/
if(ptr_file_1 == NULL) {
printf("Error opening 'str1.fa'\n");
system("PAUSE");
exit(8);
}

/*open second file*/
ptr_file_2 = fopen("d:/mpi/str2.fa", "r");
/*check to see if it opened okay*/
if(ptr_file_2 == NULL) {
printf("Error opening 'str2.fa'\n");
system("PAUSE");
exit(8);
}
/*measure file1 length*/
while(holder != EOF) {
holder = fgetc(ptr_file_1);
if(filelen1 < MAXLEN && holder !=EOF) {
strA[filelen1] = holder;
++filelen1;
}
}
strA[filelen1] = -1;
holder = '0';
/*measure file2 length*/
while(holder != EOF) {
holder = fgetc(ptr_file_2);
if(filelen2 < MAXLEN && holder !=EOF) {
strB[filelen2] = holder;
++filelen2;
}
}
strB[filelen2] = -1;
fclose(ptr_file_1);
fclose(ptr_file_2);

lenA = strlen(strA)-1;
lenB = strlen(strB)-1;





time_t t1 = time(NULL);
printf("t1=%d\n",t1);


omp_set_num_threads(2);

/*----------------------------问题代码处*---------------------------*/
for(contrl=0;contrl<lenA+lenB+1;contrl++)
{
if(contrl<=lenA)
{
rows = contrl;
j = 0;
}
else
{
rows = lenA;
j = contrl-lenA;
}

if(rank==0)
{

//#pragma omp parallel for
for(columns = j;columns <= lenB;columns++)
{
//printf("%d\n",columns);
if(rows == 0||columns == 0)
{
SWArray[rows][columns]=0;
}
else
{
left = SWArray[rows][columns-1] - GapPenalty;
top = SWArray[rows-1][columns] - GapPenalty;
if(strA[rows-1] == strB[columns-1])
{
left_top = (SWArray[rows-1][columns-1] + Match);
}else
{
left_top = SWArray[rows-1][columns-1] + MissMatch;
}
compval = SWMax(left,top,left_top);
if(compval<0) compval = 0;
SWArray[rows][columns] = compval;
buff[0]=compval;
buff[1]=rows;
buff[2]=columns;
MPI_Send(buff,4,MPI_INT,rank+1,99,MPI_COMM_WORLD);
MPI_Recv(buff,4,MPI_INT,rank+1,99,MPI_COMM_WORLD,&status);
SWArray[buff[1]][buff[2]]=buff[0];
compval = 0;
}
}
}
else
{
//#pragma omp parallel for
for(columns = j;columns <= lenB;columns++)
{
//printf("%d\n",columns);
if(rows == 0||columns == 0)
{
SWArray[rows][columns]=0;
}
else
{
left = SWArray[rows][columns-1] - GapPenalty;
top = SWArray[rows-1][columns] - GapPenalty;
if(strA[rows-1] == strB[columns-1])
{
left_top = (SWArray[rows-1][columns-1] + Match);
}else
{
left_top = SWArray[rows-1][columns-1] + MissMatch;
}
compval = SWMax(left,top,left_top);
if(compval<0) compval = 0;
SWArray[rows][columns] = compval;
buff[0]=compval;
buff[1]=rows;
buff[2]=columns;
MPI_Send(buff,4,MPI_INT,rank-1,99,MPI_COMM_WORLD);
MPI_Recv(buff,4,MPI_INT,rank-1,99,MPI_COMM_WORLD,&status);
SWArray[buff[1]][buff[2]]=buff[0];
compval = 0;
}
}

}
MPI_Barrier(MPI_COMM_WORLD);
}


/*PRINT S-W Table*/
ptr_file_3 = fopen("str3.fa", "w+");
if(ptr_file_3 == NULL) {
printf("Error opening 'str3.fa'\n");
system("PAUSE");
exit(8);
}
fprintf(ptr_file_3," 0");
for(i = 0; i <= lenB; ++i) {
fprintf(ptr_file_3," %c",strB[i]);
}
fprintf(ptr_file_3,"\n");

for(i = 0; i <= lenA; ++i) {
if(i < 1) {
fprintf(ptr_file_3,"0");
}

if(i > 0) {
fprintf(ptr_file_3,"%c",strA[i-1]);
}

for(j = 0; j <= lenB; ++j) {
fprintf(ptr_file_3,"%3i",SWArray[i][j]);
}
fprintf(ptr_file_3,"\n");
}
fclose(ptr_file_3);
/*MAKE ALIGNMENTT*/
for(i=0; i<=lenA; ++i) { //find highest score in matrix: this is the starting point of an optimal local alignment

for(j=0; j<=lenB; ++j) {

if(SWArray[i][j] > HiScore) {

HiScore = SWArray[i][j];
HiScorePos[0] = i;
HiScorePos[1] = j;
}
}
}

//send Position to alignment function
MaxA[0] = strA[HiScorePos[0]-1];
MaxB[0] = strB[HiScorePos[1]-1];

check = Align(HiScorePos[0],HiScorePos[1]);

//in the end reverse Max A and B
k=0;
for(i = strlen(MaxA)-1; i > -1; --i) {
OptA[k] = MaxA[i];
++k;
}

k=0;
for(j=strlen(MaxB)-1; j > -1; --j) {
OptB[k] = MaxB[j];
++k;
}

printf("%s\n%s \n",OptA,OptB);


time_t t2 = time(NULL);
printf("t2=%d\n",t2);
printf("Time-consuming is:%ds\n",(t2-t1));



MPI_Finalize();

return(0);

}
...全文
455 回复 打赏 收藏 转发到动态 举报
写回复
用AI写文章
回复
切换为时间正序
请发表友善的回复…
发表回复

3,881

社区成员

发帖
与我相关
我的任务
社区描述
C/C++ 其它技术问题
社区管理员
  • 其它技术问题社区
加入社区
  • 近7日
  • 近30日
  • 至今
社区公告
暂无公告

试试用AI创作助手写篇文章吧