//此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);
/*----------------------------问题代码处*---------------------------*/
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);
}
Smith-waterman算法 openmp+mpi实现,布布扣,bubuko.com
原文地址：http://blog.csdn.net/gg_wang/article/details/28102373
