INT16 CFWTproc::OnGetCoef() /* DO NOT CALL THIS FUNCTION FROM C++ SCOPE. */ /* IT MAY INTERFERE WITH THE INTERPRETER SESSION */ { INT16 __nErr = O_K; INT32 __nErrCnt = 0; MIC_CHECK; __nErrCnt = CDlpObject_GetErrorCount(); data* idCoef = MIC_GET_I_EX(idCoef,data,1,1); if (CDlpObject_GetErrorCount()>__nErrCnt) return NOT_EXEC; __nErr = GetCoef(idCoef); return __nErr; }
int main( int argc, char ** argv) { int me, NTasks; #ifndef FORCE_NOMPI int required, provided; #endif INTS id, m; /* CSC Data */ INTS n, dof; INTL nnzeros, edgenbr; COEF val; /* Local Data */ INTS localnodenbr; INTS * nodelist; INTS root; INTS base; COEF * lrhs; COEF * globrhs; COEF * globrhs_recv; COEF * globx; COEF * globx2; COEF * globprod; /* Other data */ COEF * matElem; double prec, xmin, xmax, sum1, sum2; INTS i, j, k; INTS solver; INTS zero=0; INTS one=1; INTS nb_threads; root = -1; base = 1; #ifndef FORCE_NOMPI required=MPI_THREAD_MULTIPLE; MPI_Init_thread(&argc, &argv, required, &provided); MPI_Comm_size(MPI_COMM_WORLD, &NTasks); MPI_Comm_rank(MPI_COMM_WORLD, &me); n = dof = 0; #else NTasks = 1; me = 0; #endif if (argc >= 2) { n = atoi(argv[1]); dof = atoi(argv[2]); } else { if (me == 0) { fprintf(stderr, "Usage: %s <size> <DofNumber>\n", argv[0]); return 1; } } xmin = 0.0; xmax = 1.0; /* Starting MURGE*/ CALL_MURGE(MURGE_Initialize(2)); id = 0; /* Set Options */ prec = 1e-7; /* Call MURGE_Get_Solver(solver) */ solver = MURGE_SOLVER_PASTIX; if ( solver == MURGE_SOLVER_PASTIX ) { CALL_MURGE(MURGE_SetDefaultOptions(id, zero)); CALL_MURGE(MURGE_SetOptionINT(id, IPARM_VERBOSE, API_VERBOSE_NO)); CALL_MURGE(MURGE_SetOptionINT(id, IPARM_MATRIX_VERIFICATION, API_YES)); /* CSCd Required for product in verification */ CALL_MURGE(MURGE_SetOptionINT(id, IPARM_FREE_CSCUSER, API_CSC_PRESERVE)); nb_threads = 1; #ifdef _OPENMP #pragma omp parallel shared(nb_threads) { nb_threads = omp_get_num_threads(); } #endif /* _OPENMP */ if (me == 0) { fprintf(stdout, "Running on %ld threads and %d MPI Tasks\n", (long)nb_threads, NTasks); } CALL_MURGE(MURGE_SetOptionINT(id, IPARM_THREAD_NBR, nb_threads)); } else if (solver == MURGE_SOLVER_HIPS) { #ifdef HIPS if ( method == 1 ) { CALL_MURGE(MURGE_SetDefaultOptions(id, HIPS_ITERATIVE)); } else { CALL_MURGE(MURGE_SetDefaultOptions(id, HIPS_HYBRID)); CALL_MURGE(MURGE_SetOptionINT(id, HIPS_PARTITION_TYPE, zero)); CALL_MURGE(MURGE_SetOptionINT(id, HIPS_DOMSIZE, domsize)); } CALL_MURGE(MURGE_SetOptionINT(id, HIPS_SYMMETRIC, zero)); CALL_MURGE(MURGE_SetOptionINT(id, HIPS_LOCALLY, zero)); CALL_MURGE(MURGE_SetOptionINT(id, HIPS_ITMAX, itmax)); CALL_MURGE(MURGE_SetOptionINT(id, HIPS_KRYLOV_RESTART, restart)); CALL_MURGE(MURGE_SetOptionINT(id, HIPS_VERBOSE, verbose)); CALL_MURGE(MURGE_SetOptionINT(id, HIPS_DOMNBR, NTasks)); CALL_MURGE(MURGE_SetOptionINT(id, HIPS_CHECK_GRAPH, one)); CALL_MURGE(MURGE_SetOptionINT(id, HIPS_CHECK_MATRIX, one)); #endif } CALL_MURGE(MURGE_SetOptionINT(id, MURGE_IPARAM_DOF, dof)); CALL_MURGE(MURGE_SetOptionINT(id, MURGE_IPARAM_SYM, MURGE_BOOLEAN_FALSE)); CALL_MURGE(MURGE_SetOptionINT(id, MURGE_IPARAM_BASEVAL, base)); CALL_MURGE(MURGE_SetOptionREAL(id, MURGE_RPARAM_EPSILON_ERROR, prec)); /* Set the graph : all processor enter some edge of the graph that corresponds to non-zeros location in the matrix */ /**************************************** ** Enter the matrix non-zero pattern ** ** you can use any distribution ** ****************************************/ /* this processor enters the A(myfirstrow:mylastrow, *) part of the matrix non-zero pattern */ if (me == 0) { edgenbr = 3*n-4; CALL_MURGE(MURGE_GraphBegin(id, n, edgenbr)); /* Dirichlet boundary condition */ CALL_MURGE(MURGE_GraphEdge(id, one, one)); CALL_MURGE(MURGE_GraphEdge(id, n, n)); /* Interior */ for (i = 2; i < n; i++) { for (j = -1; j <= 1; j++) { CALL_MURGE(MURGE_GraphEdge(id, i, i+j)); /* if (j != 0) { MURGE_GraphEdge(id, j+i, i); } */ } } } else { edgenbr = 0; CALL_MURGE(MURGE_GraphBegin(id, n, edgenbr)); } CALL_MURGE(MURGE_GraphEnd(id)); /* Get Local nodes */ CALL_MURGE(MURGE_GetLocalNodeNbr(id, &localnodenbr)); nodelist = (INTS*)malloc(localnodenbr*sizeof(INTS)); CALL_MURGE(MURGE_GetLocalNodeList(id, nodelist)); /* compute the number of non-zeros; */ nnzeros = 0; for (m = 0; m < localnodenbr; m++) { i = nodelist[m]; if (i == 1 || i == n) { /* Boundaries */ nnzeros = nnzeros + 1; } else { /* Interior */ for (k = -1; k <= 1; k++) { nnzeros = nnzeros + 1; } } } /* We are using dof so a non zero is in fact a block of size dof**2 */ nnzeros = nnzeros * dof*dof; /* You can enter only coefficient (i,j) that are in A(nodelist, nodelist) on this processor */ /* We enter the lower and upper triangular part of the matrix so sym = 0 */ /* matElem is the identity matrix of size 'dof' stored by line */ CALL_MURGE(MURGE_AssemblyBegin(id, n, nnzeros, MURGE_ASSEMBLY_OVW, MURGE_ASSEMBLY_OVW, MURGE_ASSEMBLY_FOOL, MURGE_BOOLEAN_FALSE)); #ifdef _OPENMP #pragma omp parallel default(none) private(m, i, k, matElem) \ shared(dof, localnodenbr, n, xmin, xmax, nodelist, id, stderr) #endif /* _OPENMP */ { matElem = (COEF*)malloc(dof*dof*sizeof(COEF)); #ifdef _OPENMP #pragma omp for #endif /* _OPENMP */ for (m = 0; m < localnodenbr; m++) { i = nodelist[m]; if ( i == 1 || i == n ) { /* Boundaries */ GetCoef(matElem, i,i,xmin,xmax,n, dof); CALL_MURGE(MURGE_AssemblySetNodeValues(id, i, i, matElem)); } else { for (k = -1; k <= 1; k++) { GetCoef(matElem,i+k,i,xmin,xmax,n, dof); CALL_MURGE(MURGE_AssemblySetNodeValues(id, i, i+k, matElem)); } } } free(matElem); } CALL_MURGE(MURGE_AssemblyEnd(id)); /* We matElem the rhs */ lrhs = (COEF*)malloc(localnodenbr*dof*sizeof(COEF)); globrhs = (COEF*)malloc(n*dof*sizeof(COEF)); for (k = 0; k < n*dof; k++) globrhs[k] = 0.0; for (m = 0; m < localnodenbr; m++) { GetRhs(&val,nodelist[m],xmin,xmax,n); for (k = 0; k < dof; k++) globrhs[(nodelist[m]-1)*dof+k] = val; for (k = 0; k < dof; k++) lrhs[m*dof+k] = val; } globrhs_recv = (COEF*)malloc(n*dof*sizeof(COEF)); #ifndef FORCE_NOMPI MPI_Allreduce(globrhs, globrhs_recv, n*dof, MURGE_MPI_COEF, MPI_SUM, MPI_COMM_WORLD); #else memcpy(globrhs_recv, globrhs, n*dof*sizeof(COEF)); #endif free(globrhs); CALL_MURGE(MURGE_SetLocalRHS(id, lrhs, MURGE_ASSEMBLY_OVW, MURGE_ASSEMBLY_OVW)); /* Get the global solution without refinement */ globx = (COEF*)malloc(n*dof*sizeof(COEF)); fprintf(stdout, "> Getting solution without refinement <\n"); CALL_MURGE(MURGE_SetOptionINT(id, IPARM_MURGE_MAY_REFINE, API_YES)); CALL_MURGE(MURGE_SetOptionINT(id, IPARM_MURGE_REFINEMENT, API_NO)); CALL_MURGE(MURGE_GetGlobalSolution(id, globx, root)); /* Get the global solution with refinement and using MURGE_GetSolution() * Test chaining two get solution. */ fprintf(stdout, "> Getting solution with refinement <\n"); globx2 = (COEF*)malloc(n*dof*sizeof(COEF)); CALL_MURGE(MURGE_SetOptionINT(id, IPARM_MURGE_REFINEMENT, API_YES)); CALL_MURGE(MURGE_GetSolution(id, n, nodelist, globx2, MURGE_ASSEMBLY_RESPECT)); /* note that we cannot check solution in-between getsolution calls because * setting X would overwrite RHS... * Maybe this will be handled in next generation of murge interface... */ fprintf(stdout, "> Check solution without refinement <\n"); CALL_MURGE(MURGE_SetGlobalRHS(id, globx, -one, MURGE_ASSEMBLY_OVW)); globprod = (COEF*)malloc(n*dof*sizeof(COEF)); CALL_MURGE(MURGE_GetGlobalProduct(id, globprod, -one)); sum1 = 0; sum2 = 0; for (k = 0; k < n*dof; k++) { sum1 = sum1 + globprod[k]*globprod[k]; sum2 = sum2 + (globprod[k] - globrhs_recv[k])*(globprod[k]-globrhs_recv[k]); } fprintf(stdout, "||AX - B||/||AX|| : %.15g\n", sqrt(sum2/sum1)); fprintf(stdout, "> Check solution with refinement <\n"); CALL_MURGE(MURGE_SetGlobalRHS(id, globx2, -one, MURGE_ASSEMBLY_OVW)); CALL_MURGE(MURGE_GetGlobalProduct(id, globprod, -one)); sum1 = 0; sum2 = 0; for (k = 0; k < n*dof; k++) { sum1 = sum1 + globprod[k]*globprod[k]; sum2 = sum2 + (globprod[k] - globrhs_recv[k])*(globprod[k]-globrhs_recv[k]); } fprintf(stdout, "||AX - B||/||AX|| : %.15g\n", sqrt(sum2/sum1)); /* Store in a file */ if (me == 0) store(globx,xmin,xmax,n,dof); { int iter; INTS n_coefs = n/NTasks; INTS *coef_idx; COEF *coef_vals; if (me < n%NTasks) n_coefs++; fprintf(stdout, "Now using MURGE_SetRHS and MURGE_ASSEMBLY_FOOL\n"); coef_idx = malloc(n_coefs*sizeof(INTS)); coef_vals = malloc(n_coefs*dof*sizeof(COEF)); /* cyclic distribution of RHS */ for (iter = 0; iter < n_coefs; iter++) { coef_idx[iter] = me + iter*NTasks + 1; /* baseval == 1 */ for (k = 0; k < dof; k++) coef_vals[iter*dof+k] = globrhs_recv[(me + iter*NTasks)*dof + k]; } CALL_MURGE(MURGE_SetRHS(id, n_coefs, coef_idx, coef_vals, MURGE_ASSEMBLY_OVW, MURGE_ASSEMBLY_OVW, MURGE_ASSEMBLY_FOOL)); free(coef_vals); free(coef_idx); CALL_MURGE(MURGE_GetGlobalSolution(id, globx, root)); CALL_MURGE(MURGE_SetGlobalRHS(id, globx, -one, MURGE_ASSEMBLY_OVW)); CALL_MURGE(MURGE_GetGlobalProduct(id, globprod, -one)); sum1 = 0; sum2 = 0; for (k = 0; k < n*dof; k++) { sum1 = sum1 + globprod[k]*globprod[k]; sum2 = sum2 + (globprod[k] - globrhs_recv[k])*(globprod[k]-globrhs_recv[k]); } fprintf(stdout, "||AX - B||/||AX|| : %.15g\n", sqrt(sum2/sum1)); } /* I'm Free */ CALL_MURGE(MURGE_Clean(id)); CALL_MURGE(MURGE_Finalize()); #ifndef FORCE_NOMPI MPI_Finalize(); #endif free(nodelist); free(lrhs); free(globx); free(globx2); free(globprod); free(globrhs_recv); return 0; }
int main( int argc, char ** argv) { INTS ierr; int me, NTasks, required, provided; INTS id, m; /* CSC Data */ INTS n, dof; INTL nnzeros, edgenbr; COEF val; /* Local Data */ INTS localnodenbr; INTS * nodelist; INTS new_localnodenbr; INTS * new_nodelist; INTS root; INTS base; COEF * lrhs; COEF * globrhs; COEF * globrhs_recv; COEF * globx; COEF * globprod; /* Other data */ COEF * matElem; double prec, xmin, xmax, sum1, sum2; INTS i, j, k; INTS solver; INTS zero=0; INTS one=1; INTS nb_threads; INTS id_seq; root = -1; base = 1; required=MPI_THREAD_MULTIPLE; MPI_Init_thread(&argc, &argv, required, &provided); MPI_Comm_size(MPI_COMM_WORLD, &NTasks); MPI_Comm_rank(MPI_COMM_WORLD, &me); n = dof = 0; if (argc >= 2) { n = atoi(argv[1]); dof = atoi(argv[2]); } else { if (me == 0) { fprintf(stderr, "Usage: %s <size> <DofNumber>\n", argv[0]); return 1; } } xmin = 0.0; xmax = 1.0; /* Starting MURGE*/ ierr = MURGE_Initialize(one); if (ierr != MURGE_SUCCESS) { fprintf(stderr, "Error %ld in MURGE_Initialize\n", (long)ierr); return 1; } id = 0; /* Set Options */ prec = 1e-7; /* Call MURGE_Get_Solver(solver) */ solver = MURGE_SOLVER_PASTIX; if ( solver == MURGE_SOLVER_PASTIX ) { MURGE_SetDefaultOptions(id, zero); MURGE_SetOptionINT(id, IPARM_VERBOSE, API_VERBOSE_NO); MURGE_SetOptionINT(id, IPARM_MATRIX_VERIFICATION, API_YES); nb_threads = 1; #ifdef _OPENMP #pragma omp parallel shared(nb_threads) { nb_threads = omp_get_num_threads(); } #endif /* _OPENMP */ if (me == 0) { fprintf(stdout, "Running on %ld threads and %d MPI Tasks\n", (long)nb_threads, NTasks); } MURGE_SetOptionINT(id, IPARM_THREAD_NBR, nb_threads); } else if (solver == MURGE_SOLVER_HIPS) { #ifdef HIPS if ( method == 1 ) { MURGE_SetDefaultOptions(id, HIPS_ITERATIVE); } else { MURGE_SetDefaultOptions(id, HIPS_HYBRID); MURGE_SetOptionINT(id, HIPS_PARTITION_TYPE, zero); MURGE_SetOptionINT(id, HIPS_DOMSIZE, domsize); } MURGE_SetOptionINT(id, HIPS_SYMMETRIC, zero); MURGE_SetOptionINT(id, HIPS_LOCALLY, zero); MURGE_SetOptionINT(id, HIPS_ITMAX, itmax); MURGE_SetOptionINT(id, HIPS_KRYLOV_RESTART, restart); MURGE_SetOptionINT(id, HIPS_VERBOSE, verbose); MURGE_SetOptionINT(id, HIPS_DOMNBR, NTasks); MURGE_SetOptionINT(id, HIPS_CHECK_GRAPH, one); MURGE_SetOptionINT(id, HIPS_CHECK_MATRIX, one); #endif } MURGE_SetOptionINT(id, MURGE_IPARAM_DOF, dof); MURGE_SetOptionINT(id, MURGE_IPARAM_SYM, MURGE_BOOLEAN_FALSE); MURGE_SetOptionINT(id, MURGE_IPARAM_BASEVAL, base); MURGE_SetOptionREAL(id, MURGE_RPARAM_EPSILON_ERROR, prec); /* Set the graph : all processor enter some edge of the graph that corresponds to non-zeros location in the matrix */ /**************************************** ** Enter the matrix non-zero pattern ** ** you can use any distribution ** ****************************************/ /* this processor enters the A(myfirstrow:mylastrow, *) part of the matrix non-zero pattern */ if (me == 0) { edgenbr = 3*n-4; MURGE_GraphBegin(id, n, edgenbr); /* Dirichlet boundary condition */ MURGE_GraphEdge(id, one, one); MURGE_GraphEdge(id, n, n); /* Interior */ for (i = 2; i < n; i++) { for (j = -1; j <= 1; j++) { MURGE_GraphEdge(id, i, i+j); /* if (j != 0) { MURGE_GraphEdge(id, j+i, i); } */ } } } else { edgenbr = 0; MURGE_GraphBegin(id, n, edgenbr); } MURGE_GraphEnd(id); /* Get Local nodes */ MURGE_GetLocalNodeNbr(id, &localnodenbr); nodelist = (INTS*)malloc(localnodenbr*sizeof(INTS)); MURGE_GetLocalNodeList(id, nodelist); new_localnodenbr = localnodenbr; new_nodelist = nodelist; if (NTasks > 1) { /* move column 1 to next proc */ int found = 0; int owner, owner_rcv, index; for (m = 0; m < localnodenbr; m++) { i = nodelist[m]; if (i == 1) found = 1; } if (found == 1) owner = me; else owner = 0; MPI_Allreduce(&owner, &owner_rcv, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD); if (owner_rcv == me) { new_localnodenbr--; new_nodelist = (INTS*)malloc(new_localnodenbr*sizeof(INTS)); index = 0; for (m = 0; m < localnodenbr; m++) { i = nodelist[m]; if (i != 1) new_nodelist[index++] = nodelist[m]; } } if (owner_rcv == (me+1)%NTasks) { index = 0; new_localnodenbr++; new_nodelist = (INTS*)malloc(new_localnodenbr*sizeof(INTS)); new_nodelist[index++] = 1; for (m = 0; m < localnodenbr; m++) { new_nodelist[index++] = nodelist[m]; } } } /* Compute the number of non-zeros; */ nnzeros = 0; for (m = 0; m < new_localnodenbr; m++) { i = new_nodelist[m]; if (i == 1 || i == n) { /* Boundaries */ nnzeros = nnzeros + 1; } else { /* Interior */ for (k = -1; k <= 1; k++) { nnzeros = nnzeros + 1; } } } /* We are using dof so a non zero is in fact a block of size dof**2 */ edgenbr = nnzeros; nnzeros = nnzeros * dof*dof; /* You can enter only coefficient (i,j) that are in A(nodelist, nodelist) on this processor */ /* We enter the lower and upper triangular part of the matrix so sym = 0 */ /* matElem is the identity matrix of size 'dof' stored by line */ { INTS * ROWs = malloc(edgenbr*sizeof(INTS)); INTS * COLs = malloc(edgenbr*sizeof(INTS)); int index = 0; for (m = 0; m < new_localnodenbr; m++) { i = new_nodelist[m]; if ( i == 1 || i == n ) { /* Boundaries */ ROWs[index] = i; COLs[index] = i; index++; } else { for (k = -1; k <= 1; k++) { ROWs[index] = i+k; COLs[index] = i; index++; } } } ierr = MURGE_AssemblySetSequence(id, edgenbr, ROWs, COLs, MURGE_ASSEMBLY_OVW, MURGE_ASSEMBLY_OVW, MURGE_ASSEMBLY_FOOL, MURGE_BOOLEAN_TRUE, &id_seq); if (ierr != MURGE_SUCCESS) { fprintf(stderr, "Error %ld in MURGE_AssemblySetSequence\n", (long)ierr); exit(1); } free(ROWs); free(COLs); } { COEF * values = malloc(nnzeros*sizeof(COEF)); int index = 0; matElem = (COEF*)malloc(dof*dof*sizeof(COEF)); for (m = 0; m < new_localnodenbr; m++) { i = new_nodelist[m]; if ( i == 1 || i == n ) { /* Boundaries */ GetCoef(matElem, i,i,xmin,xmax,n, dof); memcpy(&(values[dof*dof*index]), matElem, dof*dof*sizeof(COEF)); index++; } else { for (k = -1; k <= 1; k++) { GetCoef(matElem,i+k,i,xmin,xmax,n, dof); memcpy(&(values[dof*dof*index]), matElem, dof*dof*sizeof(COEF)); index++; } } } free(matElem); ierr = MURGE_AssemblyUseSequence(id, id_seq, values); if (ierr != MURGE_SUCCESS) { fprintf(stderr, "Error %ld in MURGE_AssemblySetSequence\n", (long)ierr); exit(1); } free(values); } /* We build the rhs */ lrhs = (COEF*)malloc(localnodenbr*dof*sizeof(COEF)); globrhs = (COEF*)malloc(n*dof*sizeof(COEF)); for (k = 0; k < n*dof; k++) globrhs[k] = 0.0; for (m = 0; m < localnodenbr; m++) { GetRhs(&val,nodelist[m],xmin,xmax,n); for (k = 0; k < dof; k++) globrhs[(nodelist[m]-1)*dof+k] = val; for (k = 0; k < dof; k++) lrhs[m*dof+k] = val; } globrhs_recv = (COEF*)malloc(n*dof*sizeof(COEF)); MPI_Allreduce(globrhs, globrhs_recv, n*dof, MURGE_MPI_COEF, MPI_SUM, MPI_COMM_WORLD); free(globrhs); MURGE_SetLocalRHS(id, lrhs, MURGE_ASSEMBLY_OVW, MURGE_ASSEMBLY_OVW); /* Get the global solution */ globx = (COEF*)malloc(n*dof*sizeof(COEF)); MURGE_GetGlobalSolution(id, globx, root); MURGE_SetGlobalRHS(id, globx, -one, MURGE_ASSEMBLY_OVW); globprod = (COEF*)malloc(n*dof*sizeof(COEF)); MURGE_GetGlobalProduct(id, globprod, -one); sum1 = 0; sum2 = 0; for (k = 0; k < n*dof; k++) { sum1 += globprod[k]*globprod[k]; sum2 += (globprod[k] - globrhs_recv[k])*(globprod[k]-globrhs_recv[k]); } fprintf(stdout, "||AX - B||/||AX|| : %.15g\n", sqrt(sum2/sum1)); /* Store in a file */ if (me == 0) store(globx,xmin,xmax,n,dof); { int iter; INTS n_coefs = n/NTasks; INTS *coef_idx; COEF *coef_vals; if (me < n%NTasks) n_coefs++; fprintf(stdout, "Now using MURGE_SetRHS and MURGE_ASSEMBLY_FOOL\n"); coef_idx = malloc(n_coefs*sizeof(INTS)); coef_vals = malloc(n_coefs*dof*sizeof(COEF)); /* cyclic distribution of RHS */ for (iter = 0; iter < n_coefs; iter++) { coef_idx[iter] = me + iter*NTasks + 1; /* baseval == 1 */ for (k = 0; k < dof; k++) coef_vals[iter*dof+k] = globrhs_recv[(me + iter*NTasks)*dof + k]; } MURGE_SetRHS(id, n_coefs, coef_idx, coef_vals, MURGE_ASSEMBLY_OVW, MURGE_ASSEMBLY_OVW, MURGE_ASSEMBLY_FOOL); free(coef_vals); free(coef_idx); MURGE_GetGlobalSolution(id, globx, root); MURGE_SetGlobalRHS(id, globx, -one, MURGE_ASSEMBLY_OVW); MURGE_GetGlobalProduct(id, globprod, -one); sum1 = 0; sum2 = 0; for (k = 0; k < n*dof; k++) { sum1 += globprod[k]*globprod[k]; sum2 += (globprod[k] - globrhs_recv[k])*(globprod[k]-globrhs_recv[k]); } fprintf(stdout, "||AX - B||/||AX|| : %.15g\n", sqrt(sum2/sum1)); } /* I'm Free */ MURGE_Clean(id); MURGE_Finalize(); MPI_Finalize(); free(nodelist); free(lrhs); free(globx); free(globprod); free(globrhs_recv); return 0; }